Re: [R] length of a string
João Fadista wrote: > Dear all, > > I would like to know how can I compute the length of a string in a dataframe. > Example: > > SEQUENCE ID > TGCTCCCATCTCCACGGHR04FS00645 > ACTGAACTCCCATCTCCAAT HR0595847847 > > I would like to know how to compute the length of each SEQUENCE. ?nchar > Best regards, > João Fadista > > [[alternative HTML version deleted]] > > > > > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Executing a DOS proram from R
Walter Paczkowski wrote:
> Good morning,
>
> I wrote a function to generate tables using LaTex. The function creates a
> .tex file with all the layouts I want. I would now like to call pdflatex.exe
> from inside this same function. Does anyone know an R command to do this?
> And also, what is the form of the pdflatex command to compile a .tex file? I
> always use WinEdt so I never actually call pdflatex.exe directly. Is the
> call just:
>
> pdflatex filename.tex
>
> or is there something else I need to know to create a pdf file using LaTex?
>
> Thanks,
>
> Walt Paczkowski
The following works for me and puts the resulting PDF into the R
current working directory:
system("pdflatex c:/myfolder/filename.tex")
?system
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Fwd: data frame row manipulation
Calle wrote:
> Hello,
>
> struggling with the very basic needs... :( any help appreciated.
>
> #using the package doBY
> #who drinks how much beer per day and therefor cannot calculate rowise
> maxvals
> evaluation=data.frame(date=c(1,2,3,4,5,6,7,8,9),
> name=c("Michael","Steve","Bob",
> "Michael","Steve","Bob","Michael","Steve","Bob"), vol=c(3,5,4,2,4,5,7,6,7))
> evaluation #
>
> maxval=summaryBy(vol ~ name,data=evaluation,FUN = function(x) { c(ma=max(x))
> } )
> maxval # over all days per person
>
> #function
> getMaxVal=function(x) { maxval$vol.ma[maxval$name==x] }
> getMaxVal("Steve") # testing the function for one name is ok
>
> #we want to add a column, that shows the daily drinkingvolume in relation to
> the persons max-vol.
> evaluation[,"relDrink"]= evaluation$vol/getMaxVal(evaluation$name)
> #
> # this brings the error:
> #
> #Warning message:
> # Korrupter Data Frame: Spalten werden abgeschnitten oder mit NAs
> # aufgefüllt in: format.data.frame(x, digits = digits, na.encode = FALSE)
>
> errortest= evaluation$vol/getMaxVal(evaluation$name)
> errortest
> # this brings:
> # numeric(0)
>
>
> #target was the following:
> #show in each line the daily consumed beer per person and in the next column
>
> #the all time max consumed beer for this person´(or divided by daily vol):
> #
> # datename vol relDrink
> #11 Michael 37
> #22 Steve 56
> #33 Bob 47
> #44 Michael 27
> #55 Steve 47
> #66 Bob 57
> #77 Michael 77
> #88 Steve 66
> #99 Bob 77
>
> # who can help???
Does this do what you want?
evl <- data.frame(date=c(1,2,3,4,5,6,7,8,9),
name=c("M","S","B","M","S","B","M","S","B"),
vol=c(3,5,4,2,4,5,7,6,7))
evl <- merge(evl, aggregate(evl$vol,
list(name = evl$name),
max),
all=TRUE)
names(evl)[4] <- "MaxDrink"
evl$RelDrink <- with(evl, vol / MaxDrink)
evl
name date vol MaxDrink RelDrink
1B6 57 0.7142857
2B3 47 0.5714286
3B9 77 1.000
4M1 37 0.4285714
5M7 77 1.000
6M4 27 0.2857143
7S2 56 0.833
8S5 46 0.667
9S8 6 6 1.000
>
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Est of SE of coefficients from lm() function
Arun Kumar Saha wrote: > Dear all R users, > > Can anyone tell me how I can get estimate of SE of coefficients from, lm() > function? > > I tried following : > > x = 1:10 > lm(x[-1]~x[-10]-1)$coefficients > > Here I got the est. of coefficient, however I also want to get some > "automated" way to get estimate of SE. > y <- 2:10 > x <- 1:9 > summary(lm(y ~ x - 1)) Call: lm(formula = y ~ x - 1) Residuals: Min 1Q Median 3Q Max -0.4211 -0.1053 0.2105 0.5263 0.8421 Coefficients: Estimate Std. Error t value Pr(>|t|) x 1.157890.02883 40.17 1.62e-10 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 0.4867 on 8 degrees of freedom Multiple R-Squared: 0.9951, Adjusted R-squared: 0.9944 F-statistic: 1613 on 1 and 8 DF, p-value: 1.624e-10 > coef(summary(lm(y ~ x - 1))) Estimate Std. Error t value Pr(>|t|) x 1.157895 0.02882750 40.16632 1.624006e-10 > coef(summary(lm(y ~ x - 1)))[,2] [1] 0.02882750 > Regards, > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Selective transformation
Allan Kamau wrote: > I am looking for a way to transform select observations based on a value > based criteria. > Why - Am learning r and would like to perform regression analysis of given > variables of the babies dataset (part of UsingR) for example babies$wt1, the > data in the variables does contain values which should be interpreted as > unknown, some variables have 999 for unknown and some have 99 for the same, > since lm() expects not available data to be marked using NA. > I would like to use a solution that does not employ loops (I think it may not > be the ideal way) > > I am looking at using apply() and supply the name of my function responsible > for transformation, but am unable to know now to reference the element of the > vector/list being currently processed by apply() so I may do in place > substitution (if value is 99 or 999) of the value with NA. Does this do what you want? babies$wt1 <- with(babies, replace(wt1, wt1 == 999, NA)) ?replace > > Pinpoint customers who are looking for what you sell. > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Covariance of data with missing values.
Ted Harding wrote: > Hi Rolf! > > Have a look at the 'norm' package. > > This does just what you;re asking for (assuming multivariate > normal, and allowing CAR missingness -- i.e. probability of > missing may depend on observed values, but must not depend on > unobserved). > > Read the documentation for the various function *very* carefully! > Drop me a line if you want more info. > > Best wishes, > Ted. Also, consider mlest() in the mvnmle package. > On 15-Aug-07 21:16:32, Rolf Turner wrote: >> I have a data matrix X (n x k, say) each row of which constitutes >> an observation of a k-dimensional random variable which I am willing, >> if not happy, to assume to be Gaussian, with mean ``mu'' and >> covariance matrix ``Sigma''. Distinct rows of X may be assumed to >> correspond to independent realizations of this random variable. >> >> Most rows of X (all but 240 out of 6000+ rows) contain one or more >> missing values. If I am willing to assume that missing entries are >> missing completely at random (MCAR) then I can estimate the covariance >> matrix Sigma via maximum likelihood, by employing the EM algorithm. >> Or so I believe. >> >> Has this procedure been implemented in R in an accessible form? >> I've had a bit of a scrounge through the searching facilities, >> and have gone through the FAQ, and have found nothing that I could >> discern to be directly relevant. >> >> Thanks for any pointers that anyone may be able to give. >> >> cheers, >> >> Rolf Turner > > > E-Mail: (Ted Harding) <[EMAIL PROTECTED]> > Fax-to-email: +44 (0)870 094 0861 > Date: 15-Aug-07 Time: 23:17:48 > -- XFMail -- > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problem with "by": does not work with ttest (but with lme)
Daniel Stahl wrote:
> Hello,
>
> I would like to do a large number of e.g. 1000 paired ttest using the
> by-function. But instead of using only the data within the 1000 groups, R
> caclulates 1000 times the ttest for the full data set(The same happens with
> Wilcoxon test). However, the by-function works fine with the lme function.
> Did I just miss something or is it really not working? If not, is there any
> other possibility to avoid loops?
> Thanks
> Daniel
>
> Here is the R help example for "by"
> require(stats)
> attach(warpbreaks)
> by(warpbreaks, tension, function(x) lm(breaks ~ wool, data = x))
> *->works great
> by(warpbreaks,tension,function(x)t.test(breaks ~ wool,data=warpbreaks,paired
> = TRUE))
> *Same output for each level of tension:
>
> tension: L
>
> Paired t-test
>
> data: breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
>
>
>
> tension: M
>
> Paired t-test
>
> data: breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
>
>
>
> tension: H
>
> Paired t-test
>
> data: breaks by wool
> t = 1.9956, df = 26, p-value = 0.05656
> alternative hypothesis: true difference in means is not equal to 0
> 95 percent confidence interval:
> -0.1735803 11.7291358
> sample estimates:
> mean of the differences
> 5.78
Try something like this:
library(MASS)
df <- mvrnorm(30, mu=c(-1,1), Sigma = diag(2))
df <- as.data.frame(df)
df$GROUP <- rep(1:3, each=10)
df.uni <- reshape(df, varying = list(c("V1","V2")), v.names="Y",
direction="long")
by(df.uni, df.uni$GROUP, function(x)t.test(Y ~ time,
data = x, paired = TRUE))
df.uni$GROUP: 1
Paired t-test
data: Y by time
t = -4.3719, df = 9, p-value = 0.001792
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-3.249894 -1.033530
sample estimates:
mean of the differences
-2.141712
---
df.uni$GROUP: 2
Paired t-test
data: Y by time
t = -6.4125, df = 9, p-value = 0.0001234
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-3.277425 -1.568074
sample estimates:
mean of the differences
-2.422749
---
df.uni$GROUP: 3
Paired t-test
data: Y by time
t = -4.4918, df = 9, p-value = 0.001507
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-3.581428 -1.182313
sample estimates:
mean of the differences
-2.381871
> --
>
> ______
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Change in R**2 for block entry regression
David Kaplan wrote: > Hi all, > > I'm demonstrating a block entry regression using R for my regression > class. For each block, I get the R**2 and the associated F. I do this > with separate regressions adding the next block in and then get the > results by writing separate summary() statements for each regression. > Is there a more convenient way to do this and also to get the change in > R**2 and associated F test for the change? > > Thanks in advance. > > David I'm not sure this is the best approach, but you might start with the data frame returned by applying anova() to several models and extend that to include the squared multiple correlation and increments: mod.0 <- lm(breaks ~ 1, data = warpbreaks) mod.1 <- lm(breaks ~ 1 + wool, data = warpbreaks) mod.2 <- lm(breaks ~ 1 + wool + tension, data = warpbreaks) mod.3 <- lm(breaks ~ 1 + wool * tension, data = warpbreaks) BlockRegSum <- anova(mod.0, mod.1, mod.2, mod.3) BlockRegSum$R2 <- 1 - (BlockRegSum$RSS / BlockRegSum$RSS[1]) BlockRegSum$IncR2 <- c(NA, diff(BlockRegSum$R2)) BlockRegSum$R2[1] <- NA BlockRegSum Analysis of Variance Table Model 1: breaks ~ 1 Model 2: breaks ~ 1 + wool Model 3: breaks ~ 1 + wool + tension Model 4: breaks ~ 1 + wool * tension Res.DfRSS Df Sum of Sq FPr(>F)R2 IncR2 1 53 9232.8 2 52 8782.1 1 450.7 3.7653 0.1 0.0488114 0.0488114 3 50 6747.9 22034.3 8.4980 0.0006926 0.3 0.2 4 48 5745.1 21002.8 4.1891 0.0210442 0.4 0.1 > BlockRegSum$R2 [1] NA 0.04881141 0.26914067 0.37775086 > BlockRegSum$IncR2 [1] NA 0.04881141 0.22032926 0.10861019 > summary(mod.1)$r.squared [1] 0.04881141 > summary(mod.2)$r.squared [1] 0.2691407 > summary(mod.3)$r.squared [1] 0.3777509 -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] tapply grand mean
Lauri Nikkinen wrote:
> Thanks Chuck but I would fancy the output made by tapply because the
> idea is to make a barplot based on those values.
>
> -Lauri
sum1 <- summary(y ~ x + f, data = d, fun=mean,
method="cross", overall=TRUE)
df <- data.frame(x = sum1$x, f = sum1$f, y = sum1$S)
df
xf y
11 lev1 6.452326
22 lev1 7.403041
33 lev1 6.117648
44 lev1 7.831390
55 lev1 6.746213
6 ALL lev1 6.910124
71 lev2 15.861256
82 lev2 17.296270
93 lev2 17.976864
10 4 lev2 19.696998
11 5 lev2 21.101952
12 ALL lev2 18.386668
13 1 lev3 61.393455
14 2 lev3 68.208299
15 3 lev3 73.479837
16 4 lev3 80.323382
17 5 lev3 87.430087
18 ALL lev3 74.167012
19 1 ALL 27.902346
20 2 ALL 30.969203
21 3 ALL 32.524783
22 4 ALL 35.950590
23 5 ALL 38.426084
24 ALL ALL 33.154601
library(lattice)
barchart(y ~ x | f, data = df, layout=c(4,1,1))
OR
barchart(S ~ x | f, data = sum1, layout=c(4,1,1))
> 2007/8/8, Chuck Cleland <[EMAIL PROTECTED]
> <mailto:[EMAIL PROTECTED]>>:
>
> Lauri Nikkinen wrote:
> > Hi R-users,
> >
> > I have a data.frame like this (modificated from
> > https://stat.ethz.ch/pipermail/r-help/2007-August/138124.html).
> >
> > y1 <- rnorm(20) + 6.8
> > y2 <- rnorm(20) + (1:20* 1.7 + 1)
> > y3 <- rnorm(20) + (1:20*6.7 + 3.7)
> > y <- c(y1,y2,y3)
> > x <- rep(1:5,12)
> > f <- gl(3,20, labels=paste("lev", 1:3, sep=""))
> > d <- data.frame(x=x,y=y, f=f)
> >
> > and this is how I can calculate mean of these levels.
> >
> > tapply(d$y, list(d$x, d$f), mean)
> >
> > But how can I calculate the mean of d$x 1 and 2 and the grand mean
> of d$x 1,
> > 2, 3, 4, 5 (within d$f) into a table?
>
> You might like the tables produced by summary.formula() in the Hmisc
> package:
>
> library(Hmisc)
>
> summary(y ~ x + f, data = d, fun=mean, method="cross", overall=TRUE)
>
> UseMethod by x, f
>
> +-+
> |N|
> |y|
> +-+
> +---+-+-+-+-+
> | x | lev1 | lev2 | lev3 | ALL |
> +---+-+-+-+-+
> |1 | 4 | 4 | 4 |12 |
> | | 6.452326|15.861256|61.393455|27.902346|
> +---+-+-+-+-+
> |2 | 4 | 4 | 4 |12 |
> | | 7.403041|17.296270|68.208299|30.969203|
> +---+-+-+-+-+
> |3 | 4 | 4 | 4 |12 |
> | | 6.117648|17.976864|73.479837|32.524783|
> +---+-+-+-+-+
> |4 | 4 | 4 | 4 |12 |
> | | 7.831390|19.696998|80.323382|35.950590|
> +---+-+-+-+-+
> |5 | 4 | 4 | 4 |12 |
> | | 6.746213|21.101952|87.430087|38.426084|
> +---+-+-+-+-+
> |ALL|20 |20 |20 |60 |
> | | 6.910124|18.386668|74.167012|33.154601|
> +---+-+-+-+-+
>
> summary(y ~ I(x %in% c(1,2)) + f, data = d, fun=mean, method="cross",
> overall=TRUE)
>
> UseMethod by I(x %in% c(1, 2)), f
>
> +-+
> |N|
> |y|
> +-+
> +-+-+-+-+-+
> |I(x %in% c(1, 2))| lev1 | lev2 | lev3 | ALL |
> +-+-+-+-+-+
> | FALSE |12 |12 |12 |36 |
> | | 6.898417|19.591938|80.411102|35.633819|
> +-+-+-+-+-+
> | TRUE | 8 | 8 | 8 |24 |
> | | 6.927684|16.578763|64.800877|29.435774|
> +-+-+-+-+-+
> | ALL|20 |20 |20 |60 |
> | | 6.910124|18.386668|74.167012|33.154601|
> +-+-+-+-+-+
>
> > Regards,
> > Lauri
> >
> > [[alternative HTML version deleted]]
> >
> > __
> > R-help@stat.math.ethz.ch <mailto:R-help@stat.math.ethz.ch> mailing
> list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> http://www.R-project.org/posti
Re: [R] tapply grand mean
Lauri Nikkinen wrote:
> Hi R-users,
>
> I have a data.frame like this (modificated from
> https://stat.ethz.ch/pipermail/r-help/2007-August/138124.html).
>
> y1 <- rnorm(20) + 6.8
> y2 <- rnorm(20) + (1:20*1.7 + 1)
> y3 <- rnorm(20) + (1:20*6.7 + 3.7)
> y <- c(y1,y2,y3)
> x <- rep(1:5,12)
> f <- gl(3,20, labels=paste("lev", 1:3, sep=""))
> d <- data.frame(x=x,y=y, f=f)
>
> and this is how I can calculate mean of these levels.
>
> tapply(d$y, list(d$x, d$f), mean)
>
> But how can I calculate the mean of d$x 1 and 2 and the grand mean of d$x 1,
> 2, 3, 4, 5 (within d$f) into a table?
You might like the tables produced by summary.formula() in the Hmisc
package:
library(Hmisc)
summary(y ~ x + f, data = d, fun=mean, method="cross", overall=TRUE)
UseMethod by x, f
+-+
|N|
|y|
+-+
+---+-+-+-+-+
| x | lev1 | lev2 | lev3 | ALL |
+---+-+-+-+-+
|1 | 4 | 4 | 4 |12 |
| | 6.452326|15.861256|61.393455|27.902346|
+---+-+-+-+-+
|2 | 4 | 4 | 4 |12 |
| | 7.403041|17.296270|68.208299|30.969203|
+---+-+-+-+-+
|3 | 4 | 4 | 4 |12 |
| | 6.117648|17.976864|73.479837|32.524783|
+---+-+-+-+-+
|4 | 4 | 4 | 4 |12 |
| | 7.831390|19.696998|80.323382|35.950590|
+---+-+-+-+-+
|5 | 4 | 4 | 4 |12 |
| | 6.746213|21.101952|87.430087|38.426084|
+---+-+-+-+-+
|ALL|20 |20 |20 |60 |
| | 6.910124|18.386668|74.167012|33.154601|
+---+-+-+-+-+
summary(y ~ I(x %in% c(1,2)) + f, data = d, fun=mean, method="cross",
overall=TRUE)
UseMethod by I(x %in% c(1, 2)), f
+-+
|N|
|y|
+-+
+-+-+-+-+-+
|I(x %in% c(1, 2))| lev1 | lev2 | lev3 | ALL |
+-+-+-+-+-+
| FALSE |12 |12 |12 |36 |
| | 6.898417|19.591938|80.411102|35.633819|
+-+-+-+-+-+
| TRUE | 8 | 8 | 8 |24 |
| | 6.927684|16.578763|64.800877|29.435774|
+-+-+-+-+-+
| ALL|20 |20 |20 |60 |
| | 6.910124|18.386668|74.167012|33.154601|
+-+-+-+-+-+
> Regards,
> Lauri
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Number of days in each month
Lauri Nikkinen wrote: > Hi R-users, > > What is the best way to achieve a table which contains all days and months > between years 2007-2020? I would like to calculate number of days in each > month within those years (to data frame). > > Regards, > > Lauri How about this? X <- seq(ISOdate(2007,1,1), ISOdate(2020,12,31), by=60*60*24) mytab <- table(substring(X, 3, 4), substring(X, 6, 7)) mytab 01 02 03 04 05 06 07 08 09 10 11 12 07 31 28 31 30 31 30 31 31 30 31 30 31 08 31 29 31 30 31 30 31 31 30 31 30 31 09 31 28 31 30 31 30 31 31 30 31 30 31 10 31 28 31 30 31 30 31 31 30 31 30 31 11 31 28 31 30 31 30 31 31 30 31 30 31 12 31 29 31 30 31 30 31 31 30 31 30 31 13 31 28 31 30 31 30 31 31 30 31 30 31 14 31 28 31 30 31 30 31 31 30 31 30 31 15 31 28 31 30 31 30 31 31 30 31 30 31 16 31 29 31 30 31 30 31 31 30 31 30 31 17 31 28 31 30 31 30 31 31 30 31 30 31 18 31 28 31 30 31 30 31 31 30 31 30 31 19 31 28 31 30 31 30 31 31 30 31 30 31 20 31 29 31 30 31 30 31 31 30 31 30 31 head(as.data.frame(mytab)) Var1 Var2 Freq 1 07 01 31 2 08 01 31 3 09 01 31 4 10 01 31 5 11 01 31 6 12 01 31 > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Multivariable correlation
John Sorkin wrote:
> Given a square matrix of variables X, is there any way to compute a
> multivariable correlation among all the variables in the array? It is
> possible to calculate the correlation of all pairs of variables in the array,
> but I want to know the correlation of all the variables taken together, i.e.
> for the matrix
>
> X=x1 x2 x3
> x4 x5 x6
> x7 x8 x9
>
> I don't want the pair-wise correlations corr(x1,x2), corr(x1,x3), corr(x1,x4)
> . . . etc., but rather
> corr(x1,x2,x3,x4,x5,x6,x7,x8,x9).
> Thanks,
> John
John:
Have you considered coefficient alpha or an intraclass correlation?
Both of these are readily available in R packages.
hope this helps,
Chuck
> John Sorkin M.D., Ph.D.
> Chief, Biostatistics and Informatics
> Baltimore VA Medical Center GRECC,
> University of Maryland School of Medicine Claude D. Pepper OAIC,
> University of Maryland Clinical Nutrition Research Unit, and
> Baltimore VA Center Stroke of Excellence
>
> University of Maryland School of Medicine
> Division of Gerontology
> Baltimore VA Medical Center
> 10 North Greene Street
> GRECC (BT/18/GR)
> Baltimore, MD 21201-1524
>
> (Phone) 410-605-7119
> (Fax) 410-605-7913 (Please call phone number above prior to faxing)
> [EMAIL PROTECTED]
>
> Confidentiality Statement:
> This email message, including any attachments, is for the so...{{dropped}}
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] extract columns of a matrix/data frame
yuvika wrote:
> Hello all,
>
> I have a matrix whose column names look like
>
> a1 a2 b1 b2 b3 c1 c2
> 1 23713 2
> 4 67814 3
>
> Now, I can have any number of a's. not just two as shown above and same
> goes for b's and c's. I need to extract all the a's columns and put them in
> another matrix, extract all b's columns and put them in some matrix and same
> goes for "c". How can I identify such pattern and get subsets of this matrix
> depending on columns names?
>
> I will appreciate a quick reply.
> Thanks a lot.
mymat <- matrix(runif(60), ncol=6)
colnames(mymat) <- c("a1","a2","b1","b2","c1","c2")
mymat
a1 a2 b1 b2 c1 c2
[1,] 0.73623481 0.25204019 0.332436396 0.36629507 0.39517285 0.62491949
[2,] 0.48867382 0.20933245 0.511805497 0.03142542 0.82168732 0.20550784
[3,] 0.89198874 0.24477456 0.629644977 0.23442137 0.17828551 0.29640615
[4,] 0.99222414 0.49044514 0.571213786 0.91068115 0.09484414 0.78108139
[5,] 0.66615787 0.13183354 0.004350679 0.32443025 0.38742483 0.76044740
[6,] 0.06642704 0.96257552 0.189716240 0.83969989 0.53470898 0.28319039
[7,] 0.31172264 0.20201281 0.577353264 0.62082694 0.31649255 0.40977000
[8,] 0.52890283 0.46576510 0.107363256 0.72534897 0.12038182 0.06295499
[9,] 0.55292555 0.76459699 0.212533012 0.73275529 0.98008863 0.85302931
[10,] 0.84320369 0.09958472 0.158443155 0.92321443 0.50935938 0.08514859
mymat[,grep("^a", colnames(mymat))]
a1 a2
[1,] 0.73623481 0.25204019
[2,] 0.48867382 0.20933245
[3,] 0.89198874 0.24477456
[4,] 0.99222414 0.49044514
[5,] 0.66615787 0.13183354
[6,] 0.06642704 0.96257552
[7,] 0.31172264 0.20201281
[8,] 0.52890283 0.46576510
[9,] 0.55292555 0.76459699
[10,] 0.84320369 0.09958472
?grep
> -
>
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] remove NA rows and columns
Antje wrote: > Hello, > > I guess, it's a rather simple thing but I cannot find a short way to reduce a > matrix, removing all rows and columns having just NA elements. > > testmatrix <- matrix(nrow=6, ncol=4) > testmatrix[2:5,2:3] <- seq(2) > > > testmatrix > [,1] [,2] [,3] [,4] > [1,] NA NA NA NA > [2,] NA11 NA > [3,] NA22 NA > [4,] NA11 NA > [5,] NA22 NA > [6,] NA NA NA NA > > the new matrix should look like this (by the way, I don't "know" which rows > and > columns are the one to be deleted... > > > testmatrix > [,1] [,2] > [1,] 11 > [2,] 22 > [3,] 11 > [4,] 22 testmatrix[!apply(is.na(testmatrix), 1, all), !apply(is.na(testmatrix), 2, all)] [,1] [,2] [1,]11 [2,]22 [3,]11 [4,]22 > Ciao, > Antje > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > > -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Extracting random parameters from summary lme and lmer
676 0.014 0.013 0.660 > typeSngl -0.421 0.080 0.007 0.033 -0.027 -0.001 0.001 > > Standardized Within-Group Residuals: > Min Q1 Med Q3 Max > -3.59074329 -0.63776965 0.03829878 0.67303837 3.33952680 > > Number of Observations: 4059 > Number of Groups: 65 > > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > > -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Plotting a smooth curve from predict
Wilson, Andrew wrote: > Probably a very simple query: > > When I try to plot a curve from a fitted polynomial, it comes out rather > jagged, not smooth like fitted curves in other stats software. Is there > a way of getting a smooth curve in R? > > What I'm doing at the moment (for the sake of example) is: > >> x <- c(1,2,3,4,5,6,7,8,9,10) > >> y <- c(10,9,8,7,6,6.5,7,8,9,10) > >> b <- data.frame(cbind(x,y)) > >> w <- gls(y ~ I(x)+I(x^2),correlation=corARMA(p=1),method="ML",data=b) > >> plot(predict(w),type="l") Make predictions for more than 10 values of x: x <- c(1,2,3,4,5,6,7,8,9,10) y <- c(10,9,8,7,6,6.5,7,8,9,10) b <- data.frame(cbind(x,y)) library(nlme) w <- gls(y ~ I(x)+I(x^2), correlation=corARMA(p=1), method="ML", data=b) plot(seq(1,10,len=100), predict(w, data.frame(x = seq(1,10, len=100))), xlab="x", ylab="Predicted y", type="l") > Many thanks, > > Andrew Wilson > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Q: extracting data from lm
D. R. Evans wrote: > Warning: I am a complete newbie to R. I have read ISwR, but I am still > finding myself completely stuck on some simple concepts. > > I have tried everything I can think of to solve this one, and finally > decided that enough was enough and I need a pointer to a solution. > > I have the following summary from lm(): > > > >> summary(lm(nu1~nu4)) > > Call: > lm(formula = nu1 ~ nu4) > > Residuals: > Min 1Q Median 3Q Max > -1572.62 -150.38 -21.70 168.57 2187.84 > > Coefficients: > Estimate Std. Error t value Pr(>|t|) > (Intercept) 29.88739 43.68881 0.6840.494 > nu4 1.000360.01025 97.599 <2e-16 *** > --- > Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 > > Residual standard error: 470.9 on 298 degrees of freedom > Multiple R-Squared: 0.9697, Adjusted R-squared: 0.9696 > F-statistic: 9526 on 1 and 298 DF, p-value: < 2.2e-16 > > > > But I want to access some of these numbers programmatically. I finally > figured out that to get the estimate of the nu4 coefficient I need to do: > > > >> lm(nu1~nu4)$coefficients[2] > nu4 > 1.000363 > > > > which to me as a long-time C++ programmer is close to black magic (I've > been programming since 1972; I have to say that R is unlike anything I've > ever seen, and it's far from trivial to get my head around some of it -- > for example, how I could have known a priori that the above is the way to > get the nu4 coefficient is beyond me). Anyway, having figured out how to > get the estimate of the coefficient, I not-unnaturally wanted also to find > a way to access the std. error of the estimate (the value 0.01025 in the > summary). But I am completely mystified as to how to do it :-( > > Any help gratefully (VERY gratefully) received, and I apologise if this is > a really, really stupid question and that the answer lies somewhere in some > documentation that I've obviously not properly taken on board. coef(summary(lm(nu1 ~ nu2)))[,2] Also, try the following which is often useful: str(summary(lm(nu1 ~ nu2))) > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Subscript out of bounds when using datadist() from Design library
Cody Hamilton wrote: > I am running R version 2.4.1 on Windows XP. I have a question regarding the > datadist() function from the Design library. I have a data.frame (call it > my.data) with 4 columns. When I submit the code > > datadist(data=my.data) > > I get the following error message: > > Error in X[[1]] : subscript out of bounds > > I suspect there may be something wrong with my data.frame (I'm certain there > is nothing wrong with datadist()), but I don't know what. Has anyone > experienced the mistake I seem to be making? If I follow the help page for datadist(), I think you want the following: datadist(my.data) Note the following in the description of the data argument: "Unless the first argument is a fit object, data must be an integer." A data frame is not a fit object, so I think that was the reason it did not work for you. > Regards, > Cody Hamilton > Edwards Lifesciences > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Operator > and
Yann Mauon wrote: > Hello, > > I'am new to the R world and have a lot of question but the first is : How to > deal with <> opertor in table objects? (Or how to deal with <> in > general...) I explain my problem. > > I read a file with the read.table expression. I then obtain a matrix. I read > the first line for example with the commande data[,1]. Then I would like to > select only the element in this line that are greater than 2. Is there an > elegant way to achieve that ? > > Thanks by advance... Note that read.table() returns a data frame, not a matrix. To subset, try this: subset(data, data[,1] > 2) OR subset(data, data[,1] > 2, select=1) Of course, it is always nice if each column in your data frame has a meaningful name, so it can be referred to by name rather than number. ?subset ?names -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] multinomial logit estimation
Walter Paczkowski wrote:
> Good morning,
>
> I'd like to estimate a simple multinomial logit model in R (not a McFadden
> conditional logit). For instance, I'd like to estimate the probability of
> someone having one of eight titles in a company with the independent
> variables being the company characteristics. A binary logit is well
> documented. What about the multinomial?
>
> Thanks,
>
> Walt Paczkowski
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
Have you considered the results of
RSiteSearch("multinomial", restrict="function")
which point to a number of relevant packages and functions?
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to combine presence only data sets to one presence/absence table
Patrick Zimmermann wrote:
> Problem: I have a Set of samples each with a list of observed species
> (presence only).
> Data is stored in a excel spreadsheet and the columns (spl) have
> different numbers of observations (spcs).
> Now I want to organize the data in a species by sample matrix with
> presence/absence style in R.
>
> data style (in excel):
>
> spl_A spl_B spl_C
> spcs1 spcs1 spcs2
> spcs2 spcs3 spcs3
> spcs4 spcs5
> spcs5
>
> desired style:
>
> spl_A spl_B spl_C
> spcs1 1 1 0
> spcs2 1 0 1
> spcs3 0 1 1
> .
> .
> .
>
> How and in which form do I import the data to R?
> (read.table() seems not to be appropriate, as data is not organized as a
> table)
>
> How can I create the species by sample matrix?
I'm not going to tackle how to read in the Excel data, but assuming
you had several vectors of species names gather together in a list, you
could construct a presence/absence data frame or matrix as follows:
spl_A <- c("spcs1","spcs2","spcs4","spcs5")
spl_B <- c("spcs1","spcs3")
spl_C <- c("spcs2","spcs3","spcs5")
mylist <- list(spl_A = spl_A, spl_B = spl_B, spl_C = spl_C)
mymat <- sapply(mylist,
function(x){as.numeric(sort(unique(unlist(mylist))) %in% x)})
rownames(mymat) <- sort(unique(unlist(mylist)))
mymat
spl_A spl_B spl_C
spcs1 1 1 0
spcs2 1 0 1
spcs3 0 1 1
spcs4 1 0 0
spcs5 1 0 1
> Thanks for any help,
> Patrick Zimmermann
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] looking at a function's code
Walter Paczkowski wrote:
> Good morning,
>
> I'd like to look at the code for the R function head. When I type just the
> word head, I get back
>
> function(x, ...)
> UseMethod("head")
>
>
>
> I expected to see several lines of R code. Any suggestions?
Have a look at:
https://svn.r-project.org/R/trunk/src/library/utils/R/head.R
Also, ?head shows methods for different types of objects, and you can
see these with
getAnywhere(head.default)
or
utils:::head.default
> Thanks,
>
> Walt Paczkowski
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to write a data.frame into n different csv-files
Dr. Guido Möser wrote:
> Hello,
>
>
> I want to write a data.frame with 10 different entries x (by index y)? into
> ten different csv-files.
> In every output file should be one of the ten numbers of x.
>
>
> I tried this one, but I can't write into ten different files, like test1.csv
> to test10.csv:
>
>
> x <- rnorm(10)
> y <- c(1:10)
> z <- data.frame(y,x)
>
> n <- nrow(z)
> for (i in 1:n) write.csv(z$x[i], file="test[i].csv")
>
>
> Can anyone help me?
Try this:
x <- rnorm(10)
y <- c(1:10)
z <- data.frame(y,x)
n <- nrow(z)
for (i in 1:n) write.csv(z$x[i], file=paste("test", i, ".csv", sep=""))
> Yours,
>
>
>
> Guido
>
>
>
> ---
> Dr. Guido Moeser
> Diplom Volkswirt
> Diplom Sozialwissenschaftler
>
> E-Mail: [EMAIL PROTECTED]
>
>
> Bei AOL gibt's jetzt kostenlos eMail f?r alle. Klicken Sie auf AOL.de um
> heraus zu finden, was es sonst noch kostenlos bei AOL gibt.
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Looping through a series of (csv) files
Yuchen Luo wrote:
> Dear Colleagues.
>
> This should be a very common operation and I believe there should be a nice
> way in R to handle it. I couldn't find it in the manual or by searching
> online. I am wondering if I could ask for some help in this community.
>
> I have 48 csv files; each stores the data for a specific month. The 48
> corresponding months are consecutively from January, 2001 to December, 2004.
> I name the files A200101, A200102,….., A200112, A200201, ……,etc.
>
> I want to process file A2000101 and store the result to a new file named
> B200101, process file A200102 and store the result to a new file named
> B200102… …etc.
>
> I do not want to manually change a little bit of the code to read a
> different file and write to a different file every time. I want the program
> to be able to loop through all the files.
>
> My question is, how to loop through the 48 files?
>
> Your help will be highly appreciated!
>
>
> Best Wishes
>
> Yuchen Luo
This will read each of the 48 csv files, add 50 to the values in each
file, and then write each resulting data frame to its own new file:
for(i in 1:48){
tempdf <- read.csv(paste("c:/myfolder/raw", i, ".csv", sep=""))
tempdf.new <- tempdf + 50
write.table(tempdf.new,
file=paste("c:/myfolder/plus50-", i, ".csv", sep=""),
sep=",", row.names=FALSE)
}
rm(list=c("tempdf", "tempdf.new", "i"))
Of course, you can do something more useful than adding 50.
> [[alternative HTML version deleted]]
>
>
>
>
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to estimate treatment-interaction contrasts
0.1 ' ' 1
summary(lm(score ~ A*B))
Call:
lm(formula = score ~ A * B)
Residuals:
Min 1Q Median 3Q Max
-1.16667 -0.29167 0.0 0.29167 0.7
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.493e+00 1.187e-0154.705 < 2e-16 ***
Aa2 6.533e-01 1.679e-01 3.892 0.000905 ***
B1 -4.000e-01 3.754e-01-1.066 0.299271
B2 -3.815e-16 3.754e-01 -1.02e-15 1.00
B3 -9.000e-01 3.754e-01-2.398 0.026373 *
B4 1.533e+00 3.754e-01 4.085 0.000577 ***
Aa2:B1 -3.667e-01 5.308e-01-0.691 0.497665
Aa2:B2 6.000e-01 5.308e-01 1.130 0.271719
Aa2:B3 1.567e+00 5.308e-01 2.951 0.007893 **
Aa2:B4 -3.667e-01 5.308e-01-0.691 0.497665
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.4597 on 20 degrees of freedom
Multiple R-Squared: 0.8038, Adjusted R-squared: 0.7156
F-statistic: 9.107 on 9 and 20 DF, p-value: 2.324e-05
Because your contrasts are correlated, the order in which they are
entered makes a difference. Thus, the two summaries only agree for the
last interaction contrast (Aa2:B4).
Perhaps more intuitively, you also can see that the tests are
sensitive to order for your contrasts as follows (no output shown):
# Second Contrast First
anova(lm(score ~
A * I(B %in% c("b1","b2")) +
A * I(B %in% "b1") +
A * I(B %in% c("b1","b2","b3")) +
A * I(B %in% c("b1","b2","b3","b4"))), test="F")
# Second Contrast Last
anova(lm(score ~ A * I(B %in% "b1") +
A * I(B %in% c("b1","b2","b3")) +
A * I(B %in% c("b1","b2","b3","b4")) +
A * I(B %in% c("b1","b2"))), test="F")
So if you really do want to evaluate your correlated contrasts
simultaneously, what makes the most sense to me is:
summary(lm(score ~ A*B))
I would be interested in other people's thoughts on this, particularly
how to use something like estimable() in the gmodels package to achieve
these contrasts.
hope this helps,
Chuck
> Quoting Chuck Cleland <[EMAIL PROTECTED]>:
>
>> [EMAIL PROTECTED] wrote:
>>> Hello, R experts,
>>> Sorry for asking this question again again since I really want a help!
>>>
>>> I have a two-factor experiment data and like to calculate estimates of
>>> interation contrasts say factor A has levels of a1, a2, and B has
>>> levels of b1, b2, b3, b4, and b5 with 3 replicates. I am not sure the
>>> constrast estimate I got is right using the script below:
>>>
>>> score<-c(7.2,6.5,6.9,6.4,6.9,6.1,6.9,5.3,7.2,5.7,5.1,5.9,7.6,6.9,6.8,
>>> 7.2,6.6,6.9,6.4,6.0,6.0,6.9,6.9,6.4,7.5,7.7,7.0,8.6,8.8,8.3)
>>>
>>> A <- gl(2, 15, labels=c("a1", "a2"))
>>> B <- rep(gl(5, 3, labels=c("b1", "b2", "b3", "b4", "b5")), 2)
>>>
>>> contrasts(B)<-cbind(c(-4,rep(1,4)),c(rep(-3,2),rep(2,3)),
>>> + c(rep(-2,3),rep(3,2)),c(rep(-1,4), rep(4,1)))
>>> fit1 <- aov(score ~ A*B)
>>> summary(fit1, split=list(B=1:4), expand.split = TRUE)
>>>Df Sum Sq Mean Sq F valuePr(>F)
>>> A1 3.2013 3.2013 15.1483 0.0009054 ***
>>> B4 8.7780 2.1945 10.3841 0.0001019 ***
>>> B: C1 1 0.0301 0.0301 0.1424 0.7099296
>>> B: C2 1 2.0335 2.0335 9.6221 0.0056199 **
>>> B: C3 1 1.2469 1.2469 5.9004 0.0246876 *
>>> B: C4 1 5.4675 5.4675 25.8715 5.637e-05 ***
>>> A:B 4 5.3420 1.3355 6.3194 0.0018616 **
>>> A:B: C11 0.7207 0.7207 3.4105 0.0796342 .
>>> A:B: C21 2.6068 2.6068 12.3350 0.0021927 **
>>> A:B: C31 1.9136 1.9136 9.0549 0.0069317 **
>>> A:B: C41 0.1008 0.1008 0.4771 0.4976647
>>> Residuals 20 4.2267 0.2113
>>> ---
>>> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>>>
>>> Now I like to get interaction contrast estimate for b1 and b2 vs
>>> b3, b4 and b5
>>> cont <- c(1, -1)[A] * c(-3, -3, 2, 2, 2)[B]
>>>
>>> estimat<-sum(cont*score) # value of the contrast estimate for A:B C2
>>>
>>>> estimat
>>> [1] -24.1
>>>
>>> I am not sure the estimate for A:B C2 contrast (-24.1) is correct
>>> because the F value given the output above(1
Re: [R] how to estimate treatment-interaction contrasts
[EMAIL PROTECTED] wrote:
> Hello, R experts,
> Sorry for asking this question again again since I really want a help!
>
> I have a two-factor experiment data and like to calculate estimates of
> interation contrasts say factor A has levels of a1, a2, and B has
> levels of b1, b2, b3, b4, and b5 with 3 replicates. I am not sure the
> constrast estimate I got is right using the script below:
>
> score<-c(7.2,6.5,6.9,6.4,6.9,6.1,6.9,5.3,7.2,5.7,5.1,5.9,7.6,6.9,6.8,
> 7.2,6.6,6.9,6.4,6.0,6.0,6.9,6.9,6.4,7.5,7.7,7.0,8.6,8.8,8.3)
>
> A <- gl(2, 15, labels=c("a1", "a2"))
> B <- rep(gl(5, 3, labels=c("b1", "b2", "b3", "b4", "b5")), 2)
>
> contrasts(B)<-cbind(c(-4,rep(1,4)),c(rep(-3,2),rep(2,3)),
> + c(rep(-2,3),rep(3,2)),c(rep(-1,4), rep(4,1)))
> fit1 <- aov(score ~ A*B)
> summary(fit1, split=list(B=1:4), expand.split = TRUE)
>Df Sum Sq Mean Sq F valuePr(>F)
> A1 3.2013 3.2013 15.1483 0.0009054 ***
> B4 8.7780 2.1945 10.3841 0.0001019 ***
> B: C1 1 0.0301 0.0301 0.1424 0.7099296
> B: C2 1 2.0335 2.0335 9.6221 0.0056199 **
> B: C3 1 1.2469 1.2469 5.9004 0.0246876 *
> B: C4 1 5.4675 5.4675 25.8715 5.637e-05 ***
> A:B 4 5.3420 1.3355 6.3194 0.0018616 **
> A:B: C11 0.7207 0.7207 3.4105 0.0796342 .
> A:B: C21 2.6068 2.6068 12.3350 0.0021927 **
> A:B: C31 1.9136 1.9136 9.0549 0.0069317 **
> A:B: C41 0.1008 0.1008 0.4771 0.4976647
> Residuals 20 4.2267 0.2113
> ---
> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>
> Now I like to get interaction contrast estimate for b1 and b2 vs b3, b4 and b5
> cont <- c(1, -1)[A] * c(-3, -3, 2, 2, 2)[B]
>
> estimat<-sum(cont*score) # value of the contrast estimate for A:B C2
>
>> estimat
> [1] -24.1
>
> I am not sure the estimate for A:B C2 contrast (-24.1) is correct
> because the F value given the output above(12.3350) does not equal to
> those I calculate below (15.2684):
>
> t.stat <- sum(cont*score)/se.contrast(fit1, as.matrix(cont))
>> t.stat^2
> Contrast 1
> 15.2684
>
> Could you please help me calculate the correct the estimate of
> interaction contrast and corresponding F value?
> Thanks in advance!
> Joshua
If the contrasts for B are orthogonal, then you get the result you
expected:
score <- c(7.2,6.5,6.9,6.4,6.9,6.1,6.9,5.3,7.2,5.7,5.1,5.9,7.6,6.9,6.8,
7.2,6.6,6.9,6.4,6.0,6.0,6.9,6.9,6.4,7.5,7.7,7.0,8.6,8.8,8.3)
A <- gl(2, 15, labels=c("a1", "a2"))
B <- rep(gl(5, 3, labels=c("b1", "b2", "b3", "b4", "b5")), 2)
contrasts(B) <- matrix(c(3, -1, 0, 0,
3, 1, 0, 0,
-2, 0, 2, 0,
-2, 0, -1, 1,
-2, 0, -1, -1), ncol=4, byrow=TRUE)
fit1 <- aov(score ~ A*B)
summary(fit1, split=list(B=1:4), expand.split = TRUE)
Df Sum Sq Mean Sq F valuePr(>F)
A1 3.2013 3.2013 15.1483 0.0009054 ***
B4 8.7780 2.1945 10.3841 0.0001019 ***
B: C1 1 1.0427 1.0427 4.9340 0.0380408 *
B: C2 1 1.0208 1.0208 4.8304 0.0399049 *
B: C3 1 1.2469 1.2469 5.9004 0.0246876 *
B: C4 1 5.4675 5.4675 25.8715 5.637e-05 ***
A:B 4 5.3420 1.3355 6.3194 0.0018616 **
A:B: C11 3.2267 3.2267 15.2684 0.0008734 ***
A:B: C21 0.1008 0.1008 0.4771 0.4976647
A:B: C31 1.9136 1.9136 9.0549 0.0069317 **
A:B: C41 0.1008 0.1008 0.4771 0.4976647
Residuals 20 4.2267 0.2113
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Note that I put the contrast of interest for B in the first column of
the contrast matrix.
hope this helps,
Chuck
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] elementary statistics with R (rkward?)
Christopher W. Ryan wrote: > As a fellow beginner, I also found Handbook of Statistical Analyses > Using R, by Brian Everitt, to be a very useful book. There is an > accompanying R package, "HSAUR". > > Also Using R for Introductory Statistics, by John Verzani. There is an > accompanying R package, "UsingR". And for Peter Dalgaard's book there is the ISwR package: http://cran.r-project.org/src/contrib/Descriptions/ISwR.html > Christopher W. Ryan, MD > SUNY Upstate Medical University Clinical Campus at Binghamton > 40 Arch Street, Johnson City, NY 13790 > cryanatbinghamtondotedu > PGP public keys available at http://home.stny.rr.com/ryancw/ > > "If you want to build a ship, don't drum up the men to gather wood, > divide the work and give orders. Instead, teach them to yearn for the > vast and endless sea." [Antoine de St. Exupery] > > Charles Annis, P.E. wrote: >> Face the music and buy the book: _Introductory Statistics with R_ by Peter >> Dalgaard. It's perfect for what you need. It's clear and concise and will >> teach you statistics AND R as painlessly as such a thing can be. It's >> inexpensive and you can get it on Amazon.com and every other major >> bookseller, including the nearest university bookstore. >> >> Charles Annis, P.E. >> >> [EMAIL PROTECTED] >> phone: 561-352-9699 >> eFax: 614-455-3265 >> http://www.StatisticalEngineering.com >> >> -Original Message- >> From: [EMAIL PROTECTED] >> [mailto:[EMAIL PROTECTED] On Behalf Of Donatas G. >> Sent: Wednesday, July 11, 2007 9:27 AM >> To: r-help@stat.math.ethz.ch >> Subject: [R] elementary statistics with R (rkward?) >> >> Hi, I am trying to learn some basic statistics stuff but I cannot find any >> elementary statistics exercises using R language. Using RKward would be even >> better... >> >> I need that in analysing sociological data, obtained through questionnairres >> - >> findind corelations between variables, relations between different types of >> data, etc. >> >> Could anyone recommend simple tutorials/exercises, available on www for me >> to >> work on? >> >> I realize it would be much simple to do this introductory stuff with spss, >> that >> everyone around me is using here in Lithuania, but I'd really like to learn >> to >> do it with R instead... >> > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] a little problem on selecting a subset from dataset A according to dataset B?
zhijie zhang wrote: > Dear Friends, >I want to extract the records from A according to B, but the results are > not correct because R says : > The length of long object is not integer times on the length of short > object. > Anybody have met the same problem? How to do it correctly? > > length(A)=47 > length(B)=6 > > A[A$coords.x1==B$X1,] #the program for the above task. I should get 6 > records, but i only get former 4 records for the above reason. > > Thanks. A[A$coords.x1 %in% B$X1,] coords.x1 coords.x2 1 542250.9 3392404 8 541512.5 3394722 9 541479.3 3394878 10 538903.4 3395943 19 543274.0 3389919 20 543840.8 3392012 ?is.element > The folloing shows dataset A and B. > > >> A >coords.x1 coords.x2 > 0 542250.89 3392404.1 > 1 538813.87 3388339.0 > 2 536049.19 3385821.6 > 3 533659.62 3383194.2 > 4 530642.30 3376834.9 > 5 529573.15 3378177.8 > 6 530853.82 3394838.8 > 7 541512.51 3394721.6 > 8 541479.33 3394877.8 > 9 538903.39 3395942.5 > 10 536019.95 3396286.1 > 11 538675.23 3384213.2 > 12 535127.95 3381255.4 > 13 533852.24 3378660.4 > 14 531360.91 3379273.8 > 15 539289.14 3375759.8 > 16 543410.51 3384353.1 > 17 543089.27 3388170.1 > 18 543274.03 3389919.2 > 19 543840.77 3392012.4 > 20 553383.55 3402401.8 > 21 554621.51 3397938.9 > 22 564096.42 3397524.4 > 23 567529.64 3398702.9 > 24 561798.76 3404864.0 > 25 562868.34 3405502.2 > 26 563145.22 3403192.1 > 27 562419.87 3404090.4 > 28 558321.85 3403879.9 > 29 567050.74 3404973.1 > 30 570609.70 3408742.4 > 31 556777.57 3397858.0 > 32 531353.38 3368596.6 > 33 533513.50 3372749.3 > 34 537543.19 3364284.8 > 35 538779.41 3368224.8 > 36 525930.09 3374067.7 > 37 522990.85 3369213.1 > 38 528826.37 3359019.0 > 39 533865.85 3362595.4 > 40 531200.25 3365053.0 > 41 551054.10 3377181.3 > 42 546974.19 3369284.8 > 43 572315.59 3359541.1 > 44 562703.63 3355173.4 > 45 558959.31 3357804.4 > 46 558531.39 3361741.1 > > >> B > X1X2 > 1 542250.89 3392404.1 > 2 541512.51 3394721.6 > 3 541479.33 3394877.8 > 4 538903.39 3395942.5 > 5 543274.03 3389919.2 > 6 543840.77 3392012.4 -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Levene Test with R
along zeng wrote:
> Hi All,
> is there Levene' test in R ? If not ,Could you give me some
> advice about Levene test with R?
> Thanks a lot! I am waiting for yours.
Did you try RSiteSearch("levene", restrict="function"), which points
to a funtion in the car package?
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] sequences
livia wrote: > Hi, I would like to generate a series in the following form (0.8^1, 0.8^2, > ..., 0.8^600) > Could anyone tell me how can I achieve that? I am really new to R. 8^(1:600) -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to save results from chisq.test or mantelhaen.test to file
[EMAIL PROTECTED] wrote: > Hi, > > I am new to these functions. I'm wondering if there is anyway to save the > entire results (all attributes of the result object) from the chisq.test or > mantelhaen.test functions? For example, from chisq.test function, you will > have statistic, parameter, p.value, expected, etc. in the result list. How > can I save all of them in one shot to, says, a text file or csv file? Thank > you. > > - adschai You could unlist() the result, coerce it to a data frame, then use write.table(). For example, something like this: write.table(as.data.frame(t(unlist(chisq.test(InsectSprays$count > 7, InsectSprays$spray, quote=FALSE) or write.table(as.data.frame(unlist(chisq.test(InsectSprays$count > 7, InsectSprays$spray))), quote=FALSE) > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Meta-Analysis of proportions
Monica Malta wrote:
> Dear colleagues,
>
> I'm conducting a meta-analysis of studies evaluating adherence of
> HIV-positive drug users into AIDS treatment, therefore I'm looking for some
> advice and syntax suggestion for running the meta-regression using
> proportions, not the usual OR/RR frequently used on RCT studies.
>
> Have already searched already several handbooks, R-manuals, mailing lists,
> professors, but... not clue at all...
>
> Does anyone have already tried this? A colleague of mine recently published a
> similar study on JAMA, but he used OpenBUGS - a software I'm not familiar
> with...
>
> If there is any tip/suggestion for a possible syntax, could someone send me?
> I need to finish this paper before my PhD qualify, but I'm completely stuck...
>
> So, any tip will be more than welcome...I will really appreciate it!!!
>
> Thanks in advance and congrats on the amazing mailing-list.
Specifying adherence as the "effect size", if you can also specify a
sampling variance perhaps you can use a tool like mima
(http://www.wvbauer.com/downloads/mima_tutorial.pdf). For example,
working with logits rather than proportions directly, consider the
following example data:
> df
ntot nadhere mod prop logitse
1 100 76 0 0.760 -1.1526795 0.2341463
2 125 98 0 0.784 -1.2891306 0.2173501
350 37 0 0.740 -1.0459686 0.3224127
4 200 159 0 0.795 -1.3553321 0.1751557
5 150 114 0 0.760 -1.1526795 0.1911796
680 56 1 0.700 -0.8472979 0.2439749
7 160 113 1 0.7062500 -0.8772402 0.1735688
8 200 130 1 0.650 -0.6190392 0.1482498
9 120 75 1 0.625 -0.5108256 0.1885618
10 105 78 1 0.7428571 -1.0608720 0.2232879
Then do the meta-regression as follows:
> source("http://www.wvbauer.com/downloads/mima.ssc";)
> with(df, mima(yi=logit, vi=se, mods=mod))
Estimate of (Residual) Heterogeneity: 0
Test for (Residual) Heterogeneity:
QE = 1.2419
df = 8
p-value = 0.9962
Parameter Estimates:
[,1]
intrcpt -1.2161
mods 0.4520
Variance-Covariance Matrix of Parameter Estimates:
intrcptmods
intrcpt 0.0436 -0.0436
mods-0.0436 0.0815
Omnibus Test of all Moderators:
QME = 2.5058
df = 1
p-value = 0.1134
Individual Moderator Tests:
estimate SEzval pvalCI_LCI_U
intrcpt -1.2161 0.2089 -5.8213 0. -1.6256 -0.8067
mods 0.4520 0.2856 1.5830 0.1134 -0.1077 1.0117
> Bests from Rio de Janeiro, Brazil.
>
> Monica
>
> Monica Malta
> Researcher
> Oswaldo Cruz Foundation - FIOCRUZ
> Social Science Department - DCS/ENSP
> Rua Leopoldo Bulhoes, 1480 - room 905
> Manguinhos
> Rio de Janeiro - RJ 21041-210
> Brazil
> phone +55.21.2598-2715
> fax +55.21.2598-2779
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Import from excel
Bruce Willy wrote:
> Hello,
>
> I have imported data from Excel using the command
> cours=read.delim("w:/apprentissage/cours_2.txt")
> after transforming my initial file with tab delimiters
>
> It seemed to work
>
> It is 2-dimensionnal. When I type cours[5,5],
> I get this type of message : "[1] 0,9760942761824 Levels: 0,495628477
> 0,893728761 0,89640592 0,903398558 ... 3,864414217"
> And when I want to multiply it, for example by 2 (cours[5,5]*2), I get :
> "Warning message:* ceci n'est pas pertinent pour des variables facteurs in:
> Ops.factor(cours[5, 5], 2)"
> i.e. approximately "this is not relevant for factor variables.
>
> What can I do if I want to manipulate my variables ?
>
> Thank you very much
You might try the following:
cours <- read.delim("w:/apprentissage/cours_2.txt", dec=",")
See the "dec" argument on the help page for read.delim().
> _
>
> stallées directement dans votre Messenger !
>
> [[alternative HTML version deleted]]
>
>
>
>
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] adding lines to stripchart
James Root wrote:
> I have two points of collection across 20 subjects (pre and post for each),
> so 20 pairs of data points. I would like to plot the actual raw data points
> for each subject for both pre and post and connect lines between these two
> points (20 in all) to depict real change between the two timepoints.
>
> I have tried using stripchart which adequately plots the two lines of
> subject data points. Attempting to use segments however has been
> difficult. It seems that the segments command gives too many coordiate
> points - so where segments has:
>
> x0, y0, x1, y1
>
> I really only have two coordinates for each point -
>
> pre to post
>
> I am sure that I am misusing the command but not sure if I should continue
> to try with segments or if there is another command that would be more
> appropriate.
>
> As always, thanks for any help.
Here is how you might use segments() and stripchart():
df <- data.frame(pre = runif(20), post = runif(20))
stripchart(x = list(df$pre, df$post), vertical=TRUE,
group.names=c("Pre","Post"))
for(i in 1:nrow(df)){segments(1, df$pre[i], 2, df$post[i])}
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] adding lines to stripchart
James Root wrote:
> I have two points of collection across 20 subjects (pre and post for each),
> so 20 pairs of data points. I would like to plot the actual raw data points
> for each subject for both pre and post and connect lines between these two
> points (20 in all) to depict real change between the two timepoints.
>
> I have tried using stripchart which adequately plots the two lines of
> subject data points. Attempting to use segments however has been
> difficult. It seems that the segments command gives too many coordiate
> points - so where segments has:
>
> x0, y0, x1, y1
>
> I really only have two coordinates for each point -
>
> pre to post
>
> I am sure that I am misusing the command but not sure if I should continue
> to try with segments or if there is another command that would be more
> appropriate.
>
> As always, thanks for any help.
How about using matplot() instead? Something like this:
matplot(t(matrix(runif(40), ncol=2)), type="b", col="black", lty=1,
xaxt="n", pch=16)
mtext(side=1, at=c(1,2), c("Pre","Post"))
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] latex of ftable (Hmisc?)
Dieter Menne wrote:
> Dear latexRs,
>
> I tried to make a latex printout of a simple categorial ftable. It should
> look like the output of print.ftable. Any ideas how to get the syntax of
> summary.formula right. Or some alternative? As far I see, xtable does not
> have method for ftable.
>
> Dieter
How about a Sweave approach? Something like this in the *.Rnw file:
\documentclass[letter]{article}
\begin{document}
@
<>=
library(Hmisc)
n=500
sex <- factor(sample(c("m","f"), n, rep=TRUE))
treatment <- factor(sample(c("Drug","Placebo"), n, rep=TRUE))
symptom <- factor(sample( c('H','S','G'), n,TRUE))
@
<>=
ftable(symptom ~ treatment + sex)
@
\end{document}
Then Sweave() the *.Rnw file to produce a *.tex file in the working
directory.
> library(Hmisc)
> n=500
> sex <- factor(sample(c("m","f"), n, rep=TRUE))
> treatment <- factor(sample(c("Drug","Placebo"), n, rep=TRUE))
> symptom <- factor(sample( c('H','S','G'), n,TRUE))
> # I want this output in latex
> ftable(symptom~treatment+sex)
> # No, it's not the same
> ss = summary(symptom~treatment+sex,fun=table)
> latex(ss)
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Extract Df under lme4
George W. Gilchrist wrote: > I need to extract the degrees of freedom and log likelihoods from a > series of mixed models computed using lmer/lme4. If I ask for logLik > (lmer.object), I get something like > > logLik(lmer.object) > 'log Lik.' -177.1000 (df=10) > > Can I easily get that df from there (or elsewhere) into an object? > Thank you for any ideas. The help page for logLik() suggests the following: > library(lme4) > fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy) > logLik(fm1) 'log Lik.' -871.8141 (df=5) > attr(logLik(fm1), "df") [1] 5 > George > > .. > George W. Gilchrist Email: [EMAIL PROTECTED] > Director of Graduate Studies Phone: (757) 221-7751 > Department of Biology, Box 8795Fax: (757) 221-6483 > College of William & Mary > Williamsburg, VA 23187-8795 > http://gwgilc.people.wm.edu/ > > > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Extracting t-tests on coefficients in lm
David C. Howell wrote: > I am writing a resampling program for multiple regression using lm(). I > resample the data 10,000 times, each time extracting the regression > coefficients. At present I extract the individual regression > coefficients using > > brg = lm(Newdv~Teach + Exam + Knowledge + Grade + Enroll) > bcoef[i,] = brg$coef > > This works fine. > > But now I want to extract the t tests on these coefficients. I cannot > find how these coefficients are stored, if at all. When I try > attributes(brg) > I do not find the t values as the attributes of the object. Typing > summary(brg) will PRINT the coefficients, their standard errors, t, and > the associated probability. I would like to type something like > tcoeff[i,] = brg$tvalue > but, of course, tvalue doesn't exist. > > Is there a simple way to extract, or compute if necessary, these values? summary(brg)$coefficients[,3] str(summary(brg)) is sometimes helpful for figuring out how to extract something. Also, you might have a look at John Fox's document on bootstraping regression models if you don't already know about it: http://cran.r-project.org/doc/contrib/Fox-Companion/appendix-bootstrapping.pdf > Thanks, > Dave Howell -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Use of the "by" command (clarification)
Economics Guy wrote:
> Well apparently this has nothing to do with the gini() command.
>
> I cannot get it to work for something as simple as sum()
>
> Here is the little example I am playing with, maybe someone can help me find
> my error:
>
> a<-c("A","B","C","A","B","C","A","A","C","B")
>
> b<-c(23,6534,456,234,7,567,345,9,565,345)
>
> c<-cbind(a,b)
>
> by(c, a, function(x) sum(b))
>
> and I get the output
>
> INDICES: A
> [1] 9085
>
> INDICES: B
> [1] 9085
> --
> INDICES: C
> [1] 9085
>
>
> Same problem as before. It is summing over the whole b vector rather than by
> the groups.
>
> Anybody have any ideas on what I am doing wrong?
Try this:
a <- c("A","B","C","A","B","C","A","A","C","B")
b <- c(23,6534,456,234,7,567,345,9,565,345)
c <- data.frame(a,b)
by(c, a, function(x) sum(x$b))
a: A
[1] 611
-
a: B
[1] 6886
-
a: C
[1] 1588
Also, consider this:
> with(c, tapply(b, list(a), sum))
ABC
611 6886 1588
> Thanks,
>
> EG
>
> On 6/16/07, Economics Guy <[EMAIL PROTECTED]> wrote:
>> I have a data set that contains income data and a group identifier. Sort
>> of like:
>>
>>
>>DATA
>>
>> Group,Income
>> A,2300
>> B,6776
>> A,6668
>> A,6768
>> B,9879
>> C,5577
>> A,7867
>> (etc),(etc)
>>
>> I am trying to compute the gini coefficient for each group.
>>
>> I have tried the following and none seem to do the trick:
>>
>> 1)
>>
>> attach(DATA)
>>
>> by(DATA, group, function(x) gini(income))
>>
>>
>> 2)
>>
>> attach(data)
>>
>> tapply(income, group, function(x) gini(income))
>>
>> Both of these return the same value for all groups. Like:
>>
>> group: A
>> [1] 0.2422496
>>
>> group: B
>> [1] 0.2422496
>>
>> group: C
>> [1] 0.2422496
>> --------
>> group: D
>> [1] 0.2422496
>>
>> Any ideas on how I can make this work? I need the fastest way since I am
>> gonna run a monte carlo based on this routine once I get the basics working.
>>
>>
>> Thanks,
>>
>> EG
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] getting R2.4 (Win)
Jeremy Miles wrote: > Hi, > > I would like to get hold of the R version 2.4.0 windows installer, it > doesn't seem to be available on CRAN (except the source, which needs > compiling). Does anyone know if it's still available anywhere? > > Thanks, > > jeremy Hi Jeremy: The windows binary for 2.4.0 seems to be available here: http://cran.r-project.org/bin/windows/base/old/ hope this helps, Chuck -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] compute new variable
Matthias von Rad wrote:
> Hello,
> maybe my question ist stupid, but I would like to calculate a new
> variable for all cases in my dataset. Inspired by the dialog in Rcmdr
> I tried
> Datenmatrix$cohigha<- with(Datenmatrix,mean (c(M2ORG, M5ORG, M8ORG,
> M11ORG), na.rm = TRUE)
> as output I got the same number for all my cases (possibly the
> overallmean of all cases), instead of a mean for each case.
> Can you help me with this problem?
Datenmatrix$cohigha <- rowMeans(Datenmatrix[,c("M2ORG", "M5ORG",
"M8ORG", "M11ORG")], na.rm=TRUE)
> regards
> Matthias
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Formating the data
A Ezhil wrote: > Hi All, > > I have a vector of length 48, something like: > 0 0 1 1 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 > 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 > > I would like to print (reformat) this vector as: > 00110011 > > by simply removing the spaces between them. I have > been trying with many option but not able to do this > task. > I would greatly appreciate your suggestion on fixing > this simple task. X <- rbinom(n=48, size=1, prob=.3) paste(X, collapse="") [1] "1001001000100011010010010111" print(paste(X, collapse=""), quote=FALSE) [1] 1001001000100011010010010111 > Thanks in advance. > > Kind regards, > Ezhil > > > > > > Bored stiff? Loosen up... > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] rlm results on trellis plot
Alan S Barnett wrote:
> How do I add to a trellis plot the best fit line from a robust fit? I
> can use panel.lm to add a least squares fit, but there is no panel.rlm
> function.
How about using panel.abline() instead of panel.lmline()?
fit1 <- coef(lm(stack.loss ~ Air.Flow, data = stackloss))
fit2 <- coef(rlm(stack.loss ~ Air.Flow, data = stackloss))
xyplot(stack.loss ~ Air.Flow, data=stackloss,
panel = function(x, y, ...){
panel.xyplot(x, y, ...)
panel.abline(fit1, type="l", col="blue")
panel.abline(fit2, type="l", col="red")
}, aspect=1)
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] logical 'or' on list of vectors
Tim Bergsma wrote:
> Suppose I have a list of logicals, such as returned by lapply:
>
> Theoph$Dose[1] <- NA
> Theoph$Time[2] <- NA
> Theoph$conc[3] <- NA
> lapply(Theoph,is.na)
>
> Is there a direct way to execute logical "or" across all vectors? The
> following gives the desired result, but seems unnecessarily complex.
>
> as.logical(apply(do.call("rbind",lapply(Theoph,is.na)),2,"sum"))
Is this what you want?
apply(is.na(Theoph), 1, any)
> Regards,
>
> Tim
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Standard errors of the predicted values from a lme (or lmer)
Dieter Menne wrote: > Fränzi Korner oikostat.ch> writes: > >> sorry for not being more specific. I would like to use R to get a prediction >> (with standard error) of the response in a mixed model at selected values of >> the fixed-effects factors. Hence, in a mixed model, say, for response body >> size with, say, fixed factors sex and age, I would like to get a prediction >> of size for each sex and at selected ages such as 5, 10, 15; and I want a SE >> for that prediction as well. > > > In that case, estimable in gmodels (by Greg Warnes, as also suggested by > Michael > Kubovy) and glht in Thorsten Hothorn's multcomp are best. The first works for > lme out of the box, the very powerful glht can handle lmer(lme4), with strong > support for multiple testing. Too bad the latter does not immediately work > with > lme, but it can be tweaked. > > In both cases, you have to construct the contrast matrix, which can be > error-prone for complex models. To my knowledge (??), there is no > simple-to-handle package that generates this matrix with an intuitive > interface. > > Dieter See contrast() in the contrast package. library(nlme) Orthodont2 <- Orthodont Orthodont2$newAge <- Orthodont$age - 11 fm1Orth.lme2 <- lme(distance ~ Sex*newAge, data = Orthodont2, random = ~ newAge | Subject) contrast( fm1Orth.lme2, a = list(Sex = levels(Orthodont2$Sex)[1], newAge = 8 - 11)) lme model parameter contrast Contrast S.E.LowerUpper t df Pr(>|t|) 22.61563 0.5265041 21.58370 23.64755 42.95 1000 > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Calculating column percentages of a table
Rehceb Rotkiv wrote: > Hello, > > I know, this is a real newbie question, but I can't find anything on this in > the manuals! I know that I get calculate the column totals of a table with > `apply(mytable, 2, sum)'. Now I want each column total to be 100% and > calculate the percentage of each field of the column. How would I do that? > Rcommander, the ultimate newb-tool, has a function `colPercents' which is > exactly what I need, but I do not want to load the whole thing into memory > just for this function. Is there a way to load just this single function > from the Rcmdr package? Does this help? > mymat <- matrix(sample(1:10, 40, replace=TRUE), ncol=4) > mymat [,1] [,2] [,3] [,4] [1,]9258 [2,]8523 [3,]4158 [4,]9867 [5,]7 107 10 [6,]1954 [7,] 10 1013 [8,]5357 [9,]6186 [10,]16 109 > prop.table(mymat, margin=2)*100 [,1] [,2] [,3] [,4] [1,] 15.00 3.636364 9.259259 12.307692 [2,] 13.33 9.090909 3.703704 4.615385 [3,] 6.67 1.818182 9.259259 12.307692 [4,] 15.00 14.545455 11.11 10.769231 [5,] 11.67 18.181818 12.962963 15.384615 [6,] 1.67 16.363636 9.259259 6.153846 [7,] 16.67 18.181818 1.851852 4.615385 [8,] 8.33 5.454545 9.259259 10.769231 [9,] 10.00 1.818182 14.814815 9.230769 [10,] 1.67 10.909091 18.518519 13.846154 > colSums(prop.table(mymat, margin=2)*100) [1] 100 100 100 100 > Many thanks, > Rehceb Rotkiv -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How should i get the quantile 2.5 % and 97.5% in each row of a matrix?
zhijie zhang wrote:
> Dear friends,
> I need the get the 2.5% and 97.5% quantile from each row of a matrix, how
> should i get it?
> BTW, i can get the min/max value from each row of a matrix, using the
> following programs, is there an easy function to do it?
>
> simmin<-matrix(NA,nrow=47,ncol=1)
> for (i in 1:47) {
> simmin[i,]<-min(datas[i,])
> }
>
> Thanks for your help.
mymat <- matrix(rnorm(200), ncol=50)
apply(mymat, 1, quantile, probs = c(0.025,0.975))
[,1] [,2] [,3] [,4]
2.5% -1.690786 -1.691700 -1.678078 -1.564438
97.5% 1.970500 1.954904 2.249030 1.862571
apply(mymat, 1, max)
[1] 2.179982 2.257138 2.772428 2.579247
apply(mymat, 1, min)
[1] -2.830661 -1.989114 -1.710050 -2.316970
?apply
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problem while working with SPSS data
Arun Kumar Saha wrote:
> Dear all R users,
>
> I got a strange problem while working with SPSS data :
>
> I wrote following :
>
> library(foreign)
> data.original = as.data.frame(read.spss(file="c:/Program Files/SPSS/Employee
> data.sav"))
>
> data = as.data.frame(cbind(data.original$MINORITY, data.original$EDUC,
> data.original$PREVEXP, data.original$JOBCAT, data.original$GENDER))
> colnames(data) = c('MINORITY', 'EDUC', 'PREVEXP', 'JOBCAT', 'GENDER')
>
> head( data.original)
>
> ID GENDER BDATE EDUC JOBCAT SALARY SALBEGIN JOBTIME PREVEXP
> MINORITY
> 1 111654150400 15 Manager 5700027000 98 144
> No
> 2 211852956800 16 Clerical 4020018750 98 36
> No
> 3 310943337600 12 Clerical 2145012000 98 381
> No
> 4 4115025184008 Clerical 2190013200 98 190
> No
> 5 511749363200 15 Clerical 4500021000 98 138
> No
> 6 611860819200 15 Clerical 3210013500 98 67
> No
>
> head( data)
> V1 V2 V3 V4 V5
> 1 1 5 144 4 NA
> 2 1 6 36 2 NA
> 3 1 3 381 2 NA
> 4 1 2 190 2 NA
> 5 1 5 138 2 NA
> 6 1 5 67 2 NA
>
>
> here I got the values of variable "V2" as 5,6,3,...etc which should
> be 15,16,12,
> can anyone tell me why I got that?
Your use of cbind() converted the factors to numeric.
> And my second question is that in my "data.original" why I got the values of
> "GENDER" as NA? Is there any way to get the actual values i.e. "m", and "f"?
Gender is of type "string" in the SPSS file, which seems to cause some
problem when you try to use the SPSS value labels. You might set the
use.value.labels argument to FALSE.
df <- read.spss(file="c:/Program Files/SPSS/Employee data.sav",
to.data.frame=TRUE, use.value.labels=FALSE)
summary(df)
IDGENDER BDATEEDUC
Min. : 1.0 f:216 Min. :1.093e+10 Min. : 8.00
1st Qu.:119.3 m:258 1st Qu.:1.153e+10 1st Qu.:12.00
Median :237.5 Median :1.197e+10 Median :12.00
Mean :237.5 Mean :1.180e+10 Mean :13.49
3rd Qu.:355.8 3rd Qu.:1.208e+10 3rd Qu.:15.00
Max. :474.0 Max. :1.225e+10 Max. :21.00
NA's :1.000e+00
JOBCAT SALARY SALBEGINJOBTIME
Min. :1.000 Min. : 15750 Min. : 9000 Min. :63.00
1st Qu.:1.000 1st Qu.: 24000 1st Qu.:12488 1st Qu.:72.00
Median :1.000 Median : 28875 Median :15000 Median :81.00
Mean :1.411 Mean : 34420 Mean :17016 Mean :81.11
3rd Qu.:1.000 3rd Qu.: 36938 3rd Qu.:17490 3rd Qu.:90.00
Max. :3.000 Max. :135000 Max. :79980 Max. :98.00
PREVEXP MINORITY
Min. : 0.00 Min. :0.
1st Qu.: 19.25 1st Qu.:0.
Median : 55.00 Median :0.
Mean : 95.86 Mean :0.2194
3rd Qu.:138.75 3rd Qu.:0.
Max. :476.00 Max. :1.
If you want to retain the labels for all of the variables and get
around the problem with gender, you might do this:
df1 <- read.spss(file="c:/Program Files/SPSS/Employee data.sav",
to.data.frame=TRUE, use.value.labels=TRUE)
df2 <- read.spss(file="c:/Program Files/SPSS/Employee data.sav",
to.data.frame=TRUE, use.value.labels=FALSE)
new.df <- merge(df1[,!names(df1) %in% "GENDER"], df2[,c("ID","GENDER")])
head(new.df)
ID BDATE EDUC JOBCAT SALARY SALBEGIN JOBTIME PREVEXP
1 1 11654150400 15 Manager 5700027000 98 144
2 2 11852956800 16 Clerical 4020018750 98 36
3 3 10943337600 12 Clerical 2145012000 98 381
4 4 115025184008 Clerical 2190013200 98 190
5 5 11749363200 15 Clerical 4500021000 98 138
6 6 11860819200 15 Clerical 3210013500 98 67
MINORITY GENDER
1 No m
2 No m
3 No f
4 No f
5 No m
6 No m
> Thanks
> Arun
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Missing data
John Sorkin wrote:
> Can someone direct me to a package, or packages, that can perform data
> augmentation to deal with missing data?
> Thanks,
> John
RSiteSearch("data augmentation", restrict="function") points to several
functions in packages norm, cat, and mix.
> John Sorkin M.D., Ph.D.
> Chief, Biostatistics and Informatics
> University of Maryland School of Medicine Division of Gerontology
> Baltimore VA Medical Center
> 10 North Greene Street
> GRECC (BT/18/GR)
> Baltimore, MD 21201-1524
> (Phone) 410-605-7119
> (Fax) 410-605-7913 (Please call phone number above prior to faxing)
>
> Confidentiality Statement:
> This email message, including any attachments, is for the so...{{dropped}}
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to select specific rows from a data frame based on values
Arin Basu wrote: > Dear Group: > > I am working with a data frame containing 316 rows of individuals > with 79 variables. Each of these 79 variables have values that range > between -4 to +4, and I want to subset this data frame so that in the > resulting new dataframe, values of _all_ of these variables should > range between -3 and +3. > > Let's say I have the following dataframe (it's a toy example with 4 > individuals and 5 variables): > > subj1 <- cbind(-4, -3, -1, -5, -7) > subj2 <- cbind(-2, -1, -1, -2, +2) > subj3 <- cbind(+2, +1, +2, +1, +2) > subj4 <- cbind(-4, -1, -2, +2, +1, +1) > > mydf <- as.data.frame(rbind(subj1, subj2, subj3, subj4)) > >>From mydf, I want to generate a new dataframe (let's call it mydf1) > which will have records of only subj2 and subj3 in it since only these > two individuals had all values for variables V1 through V5 in mydf to > range between -3 and +3. > > Documentation on subsetting and indexing data frames did not help to > solve this specific problem. There may be an obvious solution to it > but I just cannot seem to get it. > > Would greatly appreciate your inputs. mydf1 <- mydf[apply(mydf >= -3 & mydf <= 3, MARGIN=1, FUN=all),] > [relevant information: R-version: 2.4.1, running on Windows XP] > > /Arin Basu > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] read.table opening a website incl Password
Roland Rau wrote:
> Dear all,
>
> in the past I have been able to access websites with data directly. For
> example the following code works nicely
>
> mydata <-
> read.table("http://www.lifetable.de/data/MPIDR/POL_2004.txt";,
> header=TRUE)
>
> But what happens if I need a username and password (a different site)?
> How do I do that? Or is it not possible to this in R?
> I tried something like this
>
> mydata.frame <- read.table("myusr:[EMAIL PROTECTED]/adir/afile.txt")
>
> but it did not work.
> I'd appreciate any hints.
> My platform is Win32 (and I am actually running R 2.3.1, but I guess
> (hope!) this is not the reason. At least I checked the NEWS file whether
> any changes appeared since 2.3.1 which could affect this behavior).
In what way did it not work? The following seems to work for me:
read.table("ftp://myusr:[EMAIL PROTECTED]/mydir/test.dat")
> Thanks,
> Roland
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] drop a letter
elyakhlifi mustapha wrote:
> hello,
> how can I do to drop C from this character "C325" ?
> x <- "C325"
> substring(x, first=2)
[1] "325"
> gsub("C", "", x)
[1] "325"
> gsub("[A-Z]", "", x)
[1] "325"
?substring
?gsub
>
> _
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] aov problem
Anders Malmendal wrote:
> I am using R to make two-way ANOVA on a number of variables using
>
> g <- aov(var ~ fact1*fact2)
>
> where var is a matrix containing the variables.
> However the outcome seem to be dependent on the order of fact1 and fact2
> (i.e. fact2*fact1) gives a slightly (factor of 1.5) different result.
> Any ideas why this is?
RSiteSearch("sequential sums of squares")
> Thanks for any help
> Anders
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] apply( )
Greg Tarpinian wrote:
> I have a question that must have a simple answer (but eludes me).
> I need a row-by-row logical comparison across three numeric variables
> in
> a data frame: foo$x, foo$y, foo$z. The logic is
>
>if( x < y || x > z ) 1 else 0
>
> for a particular row.
>
> It is simple and very inefficient to use for(i in 1:length(foo$x)){ }
> loops. How can I accomplish this using sappy( ) / lapply( ) / apply( )
> or some other more efficient method?
X <- as.data.frame(matrix(rnorm(30), ncol=3))
X
V1 V2 V3
1 -0.48026236 0.8629789 -1.2600858
2 -1.32408219 -0.5590268 1.1310638
3 0.02717575 -0.5661402 0.7824019
4 0.80783373 0.2300440 -0.4477275
5 1.24518907 -0.3778392 1.7546530
6 -0.39254125 -1.0388962 -0.4436296
7 -1.44473455 1.8606963 0.4253889
8 -0.63543047 -1.6408418 -1.0409473
9 0.81075970 0.3914066 -1.0361739
10 1.66021280 -1.6694101 -0.4810839
with(X, ifelse(V1 < V2 | V1 > V3, 1, 0))
[1] 1 1 0 1 0 1 1 1 1 1
> Thank you in advance,
>
> Greg
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] the ifelse function
Jennifer Dillon wrote: > Hi Everyone, > > I think I found a problem with the ifelse function: If the condition > argument is NA, it treats it as true. Anyone agree or disagree with this? I disagree. > ifelse(c(5,4,NA) == 5, 1, 0) [1] 1 0 NA > Jen -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Package contrast error
Ken Nussear wrote: > Trying to use contrast to look at differences within an lme > > lme.fnl.REML <- lme(Max ~ S + Tr + Yr + Tr:Yr, random = ~1 |TID, > method = "REML") > > I have three levels of Tr I'm trying to contrast among different > years (R, T97, T98), years = 1997-1999, so I'm interested in > contrasts of the interaction term. > > > anova(lme.fnl.REML) > numDF denDF F-value p-value > (Intercept) 1 168 19255.389 <.0001 > S 1 168 5.912 0.0161 > Tr 2 11615.919 <.0001 > Yr 1 16877.837 <.0001 > Tr:Yr 2 16847.584 <.0001 > > summary(lme.fnl.REML) > Linear mixed-effects model fit by REML > Data: NULL > AIC BIClogLik >580.6991 613.5399 -281.3496 > > Random effects: > Formula: ~1 | TID > (Intercept) Residual > StdDev: 0.3697006 0.5316062 > > Fixed effects: Max ~ S + Tr + Yr + Tr:Yr > Value Std.Error DF t-value p-value > (Intercept) -13.5681 113.2623 168 -0.119793 0.9048 > SM 0.21870.0957 168 2.284605 0.0236 > TrT97 1375.5897 164.0060 116 8.387434 0. > TrT98 2890.9462 455.3497 116 6.348848 0. > Yr 0.00990.0567 168 0.174005 0.8621 > TrT97:Yr -0.68830.0821 168 -8.384798 0. > TrT98:Yr -1.44630.2279 168 -6.347310 0. > Correlation: > (Intr) SM TrT97 TrT98 Yr TT97:Y > SM0.067 > TrT97-0.691 -0.049 > TrT98-0.248 -0.001 0.171 > Yr -1.000 -0.067 0.691 0.248 > TrT97:Yr 0.691 0.048 -1.000 -0.171 -0.691 > TrT98:Yr 0.248 0.001 -0.171 -1.000 -0.248 0.171 > > Standardized Within-Group Residuals: > Min Q1 Med Q3 Max > -2.19017911 -0.58108001 -0.04983642 0.57323031 2.39811353 > > Number of Observations: 291 > Number of Groups: 119 > > > > When I try to get the contrast I get one of two errors each time. > > Trying for a paired contrast... > > tst <- contrast(lme.fnl.REML, a=list(Yr=levels(Yr), Tr="R"), b=list > (Yr=levels(Yr, Tr="T97"))) > Error in gendata.default(fit = list(modelStruct = list(reStruct = list > ( : > not enough factors > > Trying to include the other factor to make R happy > > > tst <- contrast(lme.fnl.REML, a=list(Yr=levels(Yr), Tr="R"), b=list > (Yr=levels(Yr, Tr="T97")), c=list(Yr=levels(Yr, Tr="T98")) > + ) > Error in contrastCalc(fit, ...) : argument 4 matches multiple formal > arguments > > Can anyone help with the syntax here? I believe you need to include one or more values S at which to contrast the treatments in each of the lists. So you might try something like this: contrast(lme.fnl.REML, a=list(Yr=levels(Yr), Tr="R", S="M"), b=list(Yr=levels(Yr), Tr="T97", S="M")) or contrast(lme.fnl.REML, a=list(Yr=levels(Yr), Tr="R", S=levels(S)), b=list(Yr=levels(Yr), Tr="T97", S=levels(S))) > Thanks > > Ken > > > > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Single Title for the Multiple plot page
Mohammad Ehsanul Karim wrote:
> Dear List,
>
> In R we can plot multiple graphs in same page using
> par(mfrow = c(*,*)). In each plot we can set title
> using main and sub commands.
>
> However, is there any way that we can place an
> universal title above the set of plots placed in the
> same page (not individual plot titles, all i need is a
> title of the whole graph page) as well as sib-titles?
> Do I need any package to do so?
>
> Thank you for your time.
This is covered in a number of places in the archives and
RSiteSearch("Overall Title") points you to relevant posts. I'm not sure
there is an example of having both a title and subtitle, but that is
easy enough:
par(mfrow=c(1,2), oma=c(2,0,2,0))
plot(1:10)
plot(1:10)
title("Centered Overall Title", outer=TRUE)
mtext(side=1, "Centered Subtitle", outer=TRUE)
You might consider upgrading to a more recent version of R.
> Mohammad Ehsanul Karim (R - 2.3.1 on windows)
> Institute of Statistical Research and Training
> University of Dhaka
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] Intercept Coefficient in a Model with Orthogonal Polynomials
This very likely falls in the category of an unexpected result due to user ignorance. I generated the following data: time <- 0:10 set.seed(4302007) y <- 268 + -9*time + .4*(time^2) + rnorm(11, 0, .1) I then fit models using both orthogonal and raw polynomials: fit1 <- lm(y ~ poly(time, 2)) fit2 <- lm(y ~ poly(time, degree=2, raw=TRUE)) > predict(fit1, data.frame(time = 0:10)) 1234567 268.1339 259.4912 251.6542 244.6230 238.3976 232.9780 228.3642 89 10 11 224.5562 221.5540 219.3575 217.9669 > predict(fit2, data.frame(time = 0:10)) 1234567 268.1339 259.4912 251.6542 244.6230 238.3976 232.9780 228.3642 89 10 11 224.5562 221.5540 219.3575 217.9669 > coef(fit1) (Intercept) poly(time, 2)1 poly(time, 2)2 237.00698 -52.61565 11.80144 > coef(fit2) (Intercept) 268.1339235 poly(time, degree = 2, raw = TRUE)1 -9.0456491 poly(time, degree = 2, raw = TRUE)2 0.4028944 I expected the intercept coefficient in the model with orthogonal polynomials to correspond to the predicted value of y when time=5. Instead, it seems to correspond to y at time between time=4 and time=5. Is there a way of figuring out what time the intercept corresponds to on the original time scale (0:10 here)? Any comments and pointers to references would be greatly appreciated. thanks, Chuck Cleland -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Too slow to execute!
Usman Shehu wrote:
> Greetings,
> I have the following simple function but what worries me is that it takes
> about 5 or more minutes to execute. My machine runs on windows with 1.8GHz
> and 256 Ram.
>> Re=NULL
>> for(i in 1:10){
> + x=rnorm(20)
> + Re[i]=(x-2*10)/20
> + Re
> + }
> I would appreciate any help on how to make it faster.
Why use a loop here at all? It seems like this would do it pretty fast:
Re <- rnorm(n = 2e+06, mean = - 1)
> system.time(Re <- rnorm(n = 2e+06, mean = - 1))
[1] 0.77 0.01 0.78 NA NA
I'm not sure what you expect Re to be (a list with 10 vectors of
length 20?), but you could reshape one long vector into a matrix or
whatever you might need. For example:
Re.mat <- matrix(Re, ncol=20)
> Usman
>
>
>
>
> ___
>
> now.
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Help on 'find.BIB' function
Jason Parcon wrote: > Hello everyone, > > I am trying to use the 'find.BIB' function to construct a balanced > incomplete block design. When I ran the example given in the help file > (find.BIB(10,30,4)), I obtained the following error message: > > Error in optBlock(~., withinData = factor(1:trt), blocksize = rep(k, b)) : > object ".Random.seed" not found > > I investigated what the optBlock function is doing but the manual / help > files did not give me any information regarding the above error. > > I hope somebody can help me regarding this matter. The following seems to work for me: > library(crossdes) Loading required package: AlgDesign Loading required package: gtools Loading required package: MASS > set.seed(671969) > find.BIB(10,30,4) [,1] [,2] [,3] [,4] [1,]457 10 [2,]123 10 [3,]156 10 [4,]289 10 [5,]3567 [6,]349 10 [7,]1589 [8,]1679 [9,]1247 [10,]268 10 [11,]2357 [12,]1679 [13,]267 10 [14,]1239 [15,]2568 [16,]2459 [17,]3468 [18,]158 10 [19,]2478 [20,]369 10 [21,]1246 [22,]378 10 [23,]2359 [24,]145 10 [25,]4689 [26,]479 10 [27,]1378 [28,]3456 [29,]5789 [30,]1348 I get the same error you report if I don't do the set.seed() step. > sessionInfo() R version 2.4.1 Patched (2007-03-31 r41127) i386-pc-mingw32 locale: LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States.1252;LC_MONETARY=English_United States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252 attached base packages: [1] "stats" "graphics" "grDevices" "utils" "datasets" "methods" [7] "base" other attached packages: crossdes MASSgtools AlgDesign "1.0-7" "7.2-33" "2.3.1" "1.0-7" > Best regards, > > Jason Parcon > > > > - > > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to get LSMEANS from linear mixed model?
Yannan Jiang wrote: > Hi there, > > > > I am trying to run simulations using R with linear mixed model (lme). There > are two factors in my fixed effect model, educ (treatment and control) and > mth (visit 1, 2, and 3). What I want to obtain is the estimated treatment > difference (treatment - control) at visit 3, plus the standard error and > p-value. This can be easily obtained in SAS using lsmeans or estimate > statements, but I am not sure how to do this in R. > > > > The fixed effects I obtained are as follows: > > > > Fixed effects: ymth ~ educ * mth - 1 > > Value Std.Error > DFt-value p-value > > educcont 0.14814308 0.006232419 93 > 23.769758 0. > > eductreat 0.13696952 0.006255672 93 > 21.895254 0. > > mthymth2 0.3759 0.006333043 165 > 0.005936 0.9953 > > mthymth3 0.01075489 0.006251328 165 > 1.720416 0.0872 > > eductreat:mthymth2 0.00323847 0.008947291 165 > 0.361950 0.7179 > > eductreat:mthymth3 -0.01246565 0.008941306 165 > -1.394164 0.1651 > > > > > > The estimated treatment difference I am interested are: > > > > a<-0.14814308+0.01075489 > > b<-0.13696952+0.01075489-0.01246565 > >> b-a > > [1] -0.02363921 (treatment effect at visit 3, same as SAS lsmean output) > > > > But I don't know how to get the standard error and corresponding p-value for > this estimate. Any of your helps on that would be greatly appreciated! How about fitting the model this way? df$mth <- relevel(df$mth, ref = "ymth3") lme(ymth ~ educ * mth, random = ~ 1 | id, data = df) The coefficient for educ will contain the simple effect of educ at mth=ymth3, along with a standard error and p-value. > Thanks, > > Yannan > > > > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] robust correlation
Soare Marcian-Alin wrote:
> Hello,
>
> How do I calculate a robust correlation in R?
> I want to compare it to the pearson method.
>
> library(mvoutlier)
> data(chorizon)
> cor(log10(chorizon$Al), log10(chorizon$Na), method=c("pearson"))
How about covRob() in the robust package?
covRob(log10(chorizon[,c("Al","Na")]), corr=TRUE, estim="mcd", quan =
.75, ntrial=1000)
Call:
covRob(data = log10(chorizon[, c("Al", "Na")]), corr = TRUE,
estim = "mcd", quan = 0.75, ntrial = 1000)
Robust Estimate of Correlation:
AlNa
Al 1.000 0.1929950
Na 0.1929950 1.000
Robust Estimate of Location:
Al Na
3.940185 2.185384
> KR,
> Alin Soare
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Using power.t.test over a range of conditions
Inman, Brant A. M.D. wrote:
> R-Helpers:
>
> I would like to perform sample size calculations for an experiment. As
> part of this process, I would like to know how various assumptions
> affect the sample size calculation. For instance, one thing that I
> would like to know is how the calculated sample size changes as I vary
> the difference that I would like to detect. I tried the following
> first, but got the associated error message.
>
> -
>
>> power.t.test(delta=seq(500,2000,100), sd=1000, sig.level=0.05,
> power=0.8,
> + type='two.sample', alt='two.sided')
>
> Error in uniroot(function(n) eval(p.body) - power, c(2, 1e+07)) :
> invalid function value in 'zeroin'
> In addition: Warning message:
> the condition has length > 1 and only the first element will be used in:
> if (f(lower, ...)
> * f(upper, ...) >= 0) stop("f() values at end points not of opposite
> sign")
>
> -
>
>>From the error message I suspected that the function did not handle
> vectors as arguments. I therefore tried the following looping structure
> to solve the problem:
>
> -
>
>
> DELTA <- seq(500,2000,250)
> SD <- seq(1000,2500,250)
> result <- matrix(nrow=length(DELTA), ncol=length(SD))
> colnames(result) <- paste('SD=',SD, sep='')
> rownames(result) <- paste('Delta=',DELTA, sep='')
>
> for(i in 1:length(DELTA)){
> for(j in 1:length(SD)){
> result[i,j] <- power.t.test(delta=DELTA[i], sd=SD[j],
> sig.level=0.05, power=0.8,
> type='two.sample', alt='two.sided')
> }
> }
>
> Error in result[i, j] <- power.t.test(delta = DELTA[i], sd = SD[j],
> sig.level = 0.05, :
> number of items to replace is not a multiple of replacement
> length
>
> -
>
> Can some one tell me what I am doing wrong here?
I think one problem you are having is that power.t.test() returns a
list with multiple components.
Perhaps you could go about it like this:
df <- data.frame(expand.grid(DELTA = seq(500,2000,250),
SD = seq(1000,2500,250)))
df$N <- NA
for(i in 1:dim(df)[1]){
df$N[i] <- power.t.test(delta=df$DELTA[i], sd=df$SD[i],
sig.level=0.05, power=0.8,
type='two.sample', alt='two.sided')$n
}
> Thanks in advance for your help,
>
> Brant Inman
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] general question about plotting multiple regression results
Simon Pickett wrote:
> Hi all,
>
> I have been bumbling around with r for years now and still havent come up
> with a solution for plotting reliable graphs of relationships from a
> linear regression.
>
> Here is an example illustrating my problem
>
> 1.I do a linear regression as follows
>
> summary(lm(n.day13~n.day1+ffemale.yell+fmale.yell+fmale.chroma,data=surv))
>
> which gives some nice sig. results
>
> Coefficients:
> Estimate Std. Error t value Pr(>|t|)
> (Intercept) -0.739170.43742 -1.690 0.093069 .
> n.day11.004600.05369 18.711 < 2e-16 ***
> ffemale.yell 0.224190.06251 3.586 0.000449 ***
> fmale.yell0.258740.06925 3.736 0.000262 ***
> fmale.chroma 0.235250.11633 2.022 0.044868 *
>
> 2. I want to plot the effect of "ffemale.yell", "fmale.yell" and
> "fmale.chroma" on my response variable.
>
> So, I either plot the raw values (which is fine when there is a very
> strong relationship) but what if I want to plot the effects from the
> model?
>
> In this case I would usually plot the fitted values values against the raw
> values of x... Is this the right approach?
>
> fit<-fitted(lm(n.day13~n.day1+ffemale.yell+fmale.yell+fmale.chroma,data=fsurv1))
>
> plot(fit~ffemale.yell)
>
> #make a dummy variable across the range of x
> x<-seq(from=min(fsurv1$ffemale.yell),to=max(fsurv1$ffemale.yell), length=100)
>
> #get the coefficients and draw the line
> co<-coef(lm(fit~ffemale.yell,data=fsurv1))
> y<-(co[2]*x)+co[1]
> lines(x,y, lwd=2)
>
> This often does the trick but for some reason, especially when my model
> has many terms in it or when one of the independent variables is only
> significant when the other independent variables are in the equation, it
> gives me strange lines.
>
> Please can someone show me the light?
Have a close look at John Fox's effects package. For example, the
following seems to be the kind of thing you want:
fm <- lm(n.day13 ~ n.day1 + ffemale.yell + fmale.yell +
fmale.chroma, data=fsurv1)
library(effects)
eff.fm <- effect("ffemale.yell", fm)
plot(eff.fm)
> Thanks in advance,
>
> Simon.
>
> Simon Pickett
> PhD student
> Centre For Ecology and Conservation
> Tremough Campus
> University of Exeter in Cornwall
> TR109EZ
> Tel 01326371852
>
> __________
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting command trouble
Schmitt, Corinna wrote: > Dear R-Experts, > > I have the following command lines: > > windows() > > plot(0:60, 0:0.896, type="n", xlab="Zeit [min]", ylab="Absorptionsmessung >bei 600nm",main="Zellwandstabilität" ) > > dev.off() > > > Can anyone say me why the plot command does not work and how the correct one > should look like? > Important is: x-axis goes from 0 to 60 and the y-axis from 0 to 0.896! It does not work because the lengths of x and y differ. Try evaluating 0:0.896 at the command line. If you are just trying to set up the plotting region, you could do this: par(mar=c(4,6,4,4)) plot(0, 0, type="n", xlab="Zeit [min]", ylab="Absorptionsmessung bei 600nm", main="Zellwandstabilität", xlim=c(0,60), ylim=c(0,0.896)) > Thanks, Corinna > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] A problem about all possible sequences
Paul Smith wrote: > Dear All > > Suppose a sequence of length 10 generated by the following rule: the > first element is 0 with 50% of probability or 1 with the same > probability; the second element likewise; and so on. > > Is there some R command to obtain all possible different sequences > formed by the above rule? I am aware that one could write a small > program to do that, but I am speculating about whether a command is > already existent. levels(interaction(c(0,1), c(0,1), c(0,1), c(0,1), c(0,1), c(0,1), c(0,1), c(0,1), c(0,1), c(0,1), sep="")) > Thanks in advance, > > Paul > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] indexing a subset dataframe
Mandy Barron wrote:
> Hello
> I am having problems indexing a subset dataframe, which was created
> as:
>> waspsNoGV<-subset(wasps,site!="GV")
>
> Fitting a linear model revealed some data points which had high
> leverage, so I attempted to redo the regression without these data
> points:
>> wasps.lm<-lm(r~Nt,data=waspsNoGV[-c(61,69,142),])
> which resulted in a "subscript out of bounds" error.
>
> I'm pretty sure the problem is that the data points identified in the
> regression as having high leverage were the row names carried over from
> the original dataframe which had 150 rows, but when I try to remove data
> point #142 from the subset dataframe this tries to reference by a
> numerical index but there are only 130 data points in the subset
> dataframe hence the "subscript out of bounds" message. So I guess my
> question is how do I reference the data points to drop from the
> regression by name?
Does this do it?
wasps.lm <- lm(r ~ Nt, data = subset(wasps, site != "GV" &
!(rownames(wasps) %in% c(61,69,142)))
> Thanks
> Mandy
>
> WARNING: This email and any attachments may be confidential ...{{dropped}}
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Nonparametric Effect size indices
Martin Plöderl wrote: > Hello! > > For comparing two non-normally distributed samples, Leech (2002) suggested > to report nonparametric effect size indices, such as Vargha & Delaney's A or > Cliff's d. I tried to search the R-sites, but could not find related > procedures or packages that include nonparametric effect sizes. > Thank you for your help! > > Citation: Leech (2002). A call for greater use of nonparametric statistics. > Paper presented at the Annual Meeting of the Mid-South Educational Research > Association, Chattanooga, TN, November 6-8. Based on the description of Cliff's d in http://www.florida-air.org/romano06.pdf, you could do something like the following: g1 <- c(1,1,2,2,2,3,3,3,4,5) g2 <- c(1,2,3,4,4,5) # Dominance matrix sign(outer(g1, g2, FUN="-")) [,1] [,2] [,3] [,4] [,5] [,6] [1,]0 -1 -1 -1 -1 -1 [2,]0 -1 -1 -1 -1 -1 [3,]10 -1 -1 -1 -1 [4,]10 -1 -1 -1 -1 [5,]10 -1 -1 -1 -1 [6,]110 -1 -1 -1 [7,]110 -1 -1 -1 [8,]110 -1 -1 -1 [9,]11100 -1 [10,]111110 mean(rowMeans(sign(outer(g1, g2, FUN="-" [1] -0.25 If you can point us to a description of Vargha & Delaney's A, someone can likely suggest a way of obtaining that too. > Regards, > > Martin Plöderl PhD > Suicide Prevention Research Program, Institute of Public Health > Paracelsus Private Medical University > Dept. of Suicide Prevention, University Clinic of Psychiatry and > Psychotherapy > Christian - Doppler - Klinik > Ignaz Harrerstrasse 79 > A-5020 Salzburg > AUSTRIA > Phone: +43-662-4483-4345 > Fax: +43-662-4483-4344 > E-mail: [EMAIL PROTECTED] > > __ > [EMAIL PROTECTED] mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Labelling boxplot with fivenumber summary
Daniel Siddle wrote:
> I am very new to R so forgive me if this seems basic but I have done
> extensive searching and failed to come up with the answer for myself.
>
> I am trying to label a boxplot I have created with the values for the median,
> upper and lower quartiles and max and min values. I have been unable to do
> this or find anything on the net to say how it might be done. Is this
> possible and if so how? Regards,
Here is another idea:
fn <- boxplot(ToothGrowth$len, plot=FALSE)$stats
par(mar=c(4,6,4,2))
boxplot(ToothGrowth$len, ylab="Length", at=.80)
text(1.15, fn[1], paste("Minimum Value =", fn[1]), adj=0, cex=.7)
text(1.15, fn[2], paste("Lower Quartile =", fn[2]), adj=0, cex=.7)
text(1.15, fn[3], paste("Median =", fn[3]), adj=0, cex=.7)
text(1.15, fn[4], paste("Upper Quartile =", fn[4]), adj=0, cex=.7)
text(1.15, fn[5], paste("Maximum Value =", fn[5]), adj=0, cex=.7)
arrows(1.14, fn[1], 1.02, fn[1])
arrows(1.14, fn[2], 1.02, fn[2])
arrows(1.14, fn[3], 1.02, fn[3])
arrows(1.14, fn[4], 1.02, fn[4])
arrows(1.14, fn[5], 1.02, fn[5])
title("Annotated Boxplot of Tooth Growth")
> Daniel Siddle
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Labelling boxplot with fivenumber summary
Daniel Siddle wrote: > I am very new to R so forgive me if this seems basic but I have done > extensive searching and failed to come up with the answer for myself. > > I am trying to label a boxplot I have created with the values for the median, > upper and lower quartiles and max and min values. I have been unable to do > this or find anything on the net to say how it might be done. Is this > possible and if so how? Regards, This message from back in 2002 gives a function called bp.example(), which shows how a boxplot might be annotated: http://tolstoy.newcastle.edu.au/R/help/02a/1515.html You could easily modify it into a stripped down version that does what you want. > Daniel Siddle > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting multilevel / lme lines
Rense Nieuwenhuis wrote:
> Dear expeRts,
>
> I am trying to plot a lme-object {package nlme) in such a way, that
> on a selected level the x-axis represents the value on a selected
> predictor and the y-axis represents the predicted-outcome variable.
> The graphs would than consist of several lines that each represent
> one group. I can't find such a plotting function.
>
> I could write such a function myself, based on ranef() and fixef(),
> but it would be a waste of time if such a function would already exist.
>
> Does any of you such a function?
I don't know of a single function with an lme object as argument, but
for what I think you have in mind, here is how you might go about it:
library(nlme)
fm2 <- lme(distance ~ poly(age, 2) * Sex,
data = Orthodont, random = ~ 1)
newdat <- expand.grid(age = 8:14, Sex = c("Male","Female"))
newdat$PREDDIST <- predict(fm2, newdat, level = 0)
library(lattice)
xyplot(PREDDIST ~ age, groups=Sex, ylab="Model Predicted Distance",
data = newdat, xlab="Age",
panel = function(x, y, ...){
panel.grid(h=6,v=6)
panel.superpose(x, y, type="l", ...)},
main="Orthodont Growth Model",
key = simpleKey(levels(newdat$Sex),
lines=TRUE, points=FALSE)
)
> Regards,
>
> Rense Nieuwenhuis
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] composed matrices
Eva Ubl wrote: > I have just started programming in R and so my question might be basic but > I cant find it in the manual > I have four matrices A,B,C and D and want them to compose in a big matrix X > with A in the upper left corner, B in the right upper corner, c below A > and D below B. A <- matrix(1:4, ncol=2) B <- matrix(5:8, ncol=2) C <- matrix(9:12, ncol=2) D <- matrix(13:16, ncol=2) rbind(cbind(A,B), cbind(C,D)) [,1] [,2] [,3] [,4] [1,]1357 [2,]2468 [3,]9 11 13 15 [4,] 10 12 14 16 ?rbind ?cbind > Thanks for your help > eva > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] multi-level modeling & R?
John Kane wrote:
> A colleague was asking me if R does multi-level
> modelling as opposed to multiple regression. Since I
> have no knowledge of multi-level modelling (except 5
> minutes googling ) I thought that I would as here.
>
> Does are offer any multi-level modeling packages? It
> looked like arm might be one but I was not sure.
RSiteSearch("multilevel model") and a search for "multilevel" on CRAN
point to a number of other relevant packages and docs, including:
http://finzi.psych.upenn.edu/R/library/nlme/html/lme.html
http://cran.r-project.org/src/contrib/Descriptions/lme4.html
http://cran.r-project.org/doc/packages/mlmRev.pdf
http://cran.r-project.org/src/contrib/Descriptions/multilevel.html
http://cran.r-project.org/doc/contrib/Bliese_Multilevel.pdf
> Thanks
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting symbol
Peter Dalgaard wrote: > Chuck Cleland wrote: >> Federico Calboli wrote: >> >>> Hi All, >>> >>> can I have a plot where the symbol for the dots is smaller than pch >>> =20 but bigger than pch = '.'? >>> >> >> Have you considered using the cex argument to reduce the size of >> pch=20? >> >> X <- rnorm(20) >> >> par(mfrow=c(2,2)) >> plot(X, pch=20, cex=1.0) >> plot(X, pch=20, cex=0.6) >> plot(X, pch=20, cex=1.5) >> plot(X, pch=".") >> >> > Notice though, that it depends on the device. On X11 (presumably Windows > too) at "normal" resolutions, you really have only one sice which is > less than the default. Try e.g. > > plot(1:10, rep(0,10), pch=20,cex=seq(1,.1,,10)) Peter: Thanks for making that point, which certainly did not occur to me. I see the same thing with the windows device with pointsize at the default of 12 (win.graph too). With pointsize=15 I can get two sizes smaller than the default. And on the pdf device I get many sizes less than the default. >>> Best, >>> >>> Fede >>> >>> -- >>> Federico C. F. Calboli >>> Department of Epidemiology and Public Health >>> Imperial College, St. Mary's Campus >>> Norfolk Place, London W2 1PG >>> >>> Tel +44 (0)20 75941602 Fax +44 (0)20 75943193 >>> >>> f.calboli [.a.t] imperial.ac.uk >>> f.calboli [.a.t] gmail.com >>> >>> __ >>> R-help@stat.math.ethz.ch mailing list >>> https://stat.ethz.ch/mailman/listinfo/r-help >>> PLEASE do read the posting guide >>> http://www.R-project.org/posting-guide.html >>> and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting symbol
Federico Calboli wrote: > Hi All, > > can I have a plot where the symbol for the dots is smaller than pch > =20 but bigger than pch = '.'? Have you considered using the cex argument to reduce the size of pch=20? X <- rnorm(20) par(mfrow=c(2,2)) plot(X, pch=20, cex=1.0) plot(X, pch=20, cex=0.6) plot(X, pch=20, cex=1.5) plot(X, pch=".") > Best, > > Fede > > -- > Federico C. F. Calboli > Department of Epidemiology and Public Health > Imperial College, St. Mary's Campus > Norfolk Place, London W2 1PG > > Tel +44 (0)20 75941602 Fax +44 (0)20 75943193 > > f.calboli [.a.t] imperial.ac.uk > f.calboli [.a.t] gmail.com > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] export table
Sergio Della Franca wrote: > Dear R-Helpers, > > I have a problem. > > I want to export from R to .txt my data set(y): > > YEARS PRODUCTS > 1990 10 > 1995 15 > 1997 26 > 1998 29 > 2000 34 > > > I used this code: > > write.table(y,"D:/Desktop/export_table.txt"). > > This procedure run correctly, but i acquired this result: > > YEARS PRODUCTS > 1 1990 10 > 2 1995 15 > 3 1997 26 > 4 1998 29 > 5 2000 34 > > The prolem is that R add a column in the export procedure, but it doesn't > give a name at column then this new column get the name of the first column > of my data set. > > There is a command that a must add to my export procedure to not export this > column or to give at this a name? write.table(y,"D:/Desktop/export_table.txt", row.names=FALSE) The row.names and col.names arguments are described in the help page for write.table(). > Thank you in advance. > > > Sergio. > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to get "lsmeans"?
Liaw, Andy wrote:
> I verified the result from the following with output from JMP 6 on the
> same data (don't have SAS: don't need it):
>
> set.seed(631)
> n <- 100
> dat <- data.frame(y=rnorm(n), A=factor(sample(1:2, n, replace=TRUE)),
> B=factor(sample(1:2, n, replace=TRUE)),
> C=factor(sample(1:2, n, replace=TRUE)),
> d=rnorm(n))
> fm <- lm(y ~ A + B + C + d, dat)
> ## Form a data frame of points to predict: all combinations of the
> ## three factors and the mean of the covariate.
> p <- data.frame(expand.grid(A=1:2, B=1:2, C=1:2))
> p[] <- lapply(p, factor)
> p <- cbind(p, d=mean(dat$d))
> p <- cbind(yhat=predict(fm, p), p)
> ## lsmeans for the three factors:
> with(p, tapply(yhat, A, mean))
> with(p, tapply(yhat, B, mean))
> with(p, tapply(yhat, C, mean))
Using Andy's example data, these are the LSMEANS and intervals I get
from SAS:
Ay LSMEAN 95% Confidence Limits
1 -0.071847 -0.387507 0.243813
2 -0.029621 -0.342358 0.283117
By LSMEAN 95% Confidence Limits
1 -0.104859 -0.397935 0.188216
20.003391 -0.333476 0.340258
Cy LSMEAN 95% Confidence Limits
1 -0.084679 -0.392343 0.222986
2 -0.016789 -0.336374 0.302795
One way of reproducing the LSMEANS and intervals from SAS using
predict() seems to be the following:
> dat.lm <- lm(y ~ A + as.numeric(B) + as.numeric(C) + d, data = dat)
> newdat <- expand.grid(A=factor(c(1,2)),B=1.5,C=1.5,d=mean(dat$d))
> cbind(newdat, predict(dat.lm, newdat, interval="confidence"))
A B C d fitlwr upr
1 1 1.5 1.5 0.09838595 -0.07184709 -0.3875070 0.2438128
2 2 1.5 1.5 0.09838595 -0.02962086 -0.3423582 0.2831165
However, another possibility seems to be:
> dat.lm <- lm(y ~ A + as.numeric(B) + as.numeric(C) + d, data = dat)
> newdat <-
expand.grid(A=factor(c(1,2)),B=mean(as.numeric(dat$B)),C=mean(as.numeric(dat$C)),d=mean(dat$d))
> cbind(newdat, predict(dat.lm, newdat, interval="confidence"))
ABC d fitlwr upr
1 1 1.43 1.48 0.09838595 -0.08078243 -0.3964661 0.2349012
2 2 1.43 1.48 0.09838595 -0.03855619 -0.3479589 0.2708465
The predictions directly above match what effect() in the effects
package by John Fox returns:
library(effects)
> effect("A", fm, xlevels=list(d = mean(dat$D)))
A effect
A
1 2
-0.08078243 -0.03855619
But for some reason the predict() and effect() intervals are a little
different:
> effect("A", fm, xlevels=list(d = mean(dat$D)))$lower
[,1]
101 -0.3924451
102 -0.3440179
> effect("A", fm, xlevels=list(d = mean(dat$D)))$upper
[,1]
101 0.2308802
102 0.2669055
I would be interested in any comments on these different approaches
and on the difference in intervals returned by predict() and effect().
hope this helps,
Chuck
> Andy
>
> From: Xingwang Ye
>> Dear all,
>>
>> I search the mail list about this topic and learn that no
>> simple way is available to get "lsmeans" in R as in SAS.
>> Dr.John Fox and Dr.Frank E Harrell have given very useful
>> information about "lsmeans" topic.
>> Dr. Frank E Harrell suggests not to think about lsmeans,
>> just to think about what predicted values wanted
>> and to use the predict function. However, after reading
>> the R help file for a whole day, I am still unclear how to do it.
>> Could some one give me a hand?
>>
>> for example:
>>
>> A,B and C are binomial variables(factors); d is a continuous
>> variable ; The response variable Y is a continuous variable too.
>>
>> To get lsmeans of Y according to A,B and C, respectively, in
>> SAS, I tried proc glm data=a; class A B C; model Y=A B C d;
>> lsmeans A B C/cl; run;
>>
>> In R, I tried this:
>> library(Design)
>> ddist<-datadist(a)
>> options(datadist="ddist")
>> f<-ols(Y~A+B+C+D,data=a,x=TRUE,y=TRUE,se.fit=TRUE)
>>
>> then how to get the "lsmeans" for A, B, and C, respectively
>> with predict function?
>>
>>
>>
>> Best wishes
>> yours, sincerely
>> Xingwang Ye
>> PhD candidate
>> Research Group of Nutrition Related Cancers and Other Chronic
>> Diseases
>> Institute for Nutritional Sciences,
>> Shanghai Institutes of Biological Sciences,
>> Chinese Academy of Sciences
>> P.O.Box 32
>> 294 Taiyuan Road
>> Shanghai 200031
>> P.R.CHINA
>>
Re: [R] multcomp
Nair, Murlidharan T wrote:
> Thanks for the hint.
> Here is what I had done earlier according to the old package
>
> mult.comp<-simint(a~b, data=z, conf.level=0.99, type=c("Tukey") )
> isoforms<-as.vector(rownames(mult.comp$estimate))
> estimate<-as.vector(mult.comp$estimate)
> lower<-as.vector(mult.comp$conf.int[,1])
> upper<-as.vector(mult.comp$conf.int[,2])
> p.val.raw<-as.vector(mult.comp$p.value.raw)
>
> Here is how I modified the above so that I can use the new package. Can
> someone comment if what I have done is correct or not?
> amod<-aov(a~b, data=z)
> mult.comp<-glht(amod, linfct=mcp(b="Tukey"))
> estimate<-confint(mult.comp)[1,1]
> lower<-confint(mult.comp)[1,2]
> upper<-confint(mult.comp)[1,3]
> summary(mult.comp) does give the adjusted pvalues. I tried all the
> subscripts to extract it to store in a variable but could not, can any one
> help me here. Also, how do I specify the conf.level to be equal to 0.99. The
> default is 0.95
To specify a confidence level other than the default, see the help
page for confint().
library(multcomp)
amod <- aov(breaks ~ tension, data = warpbreaks)
mc <- glht(amod, linfct = mcp(tension = "Tukey"))
confint(mc, level = .99)
Simultaneous Confidence Intervals for General Linear Hypotheses
Multiple Comparisons of Means: Tukey Contrasts
Fit: aov(formula = breaks ~ tension, data = warpbreaks)
Estimated Quantile = 3.0494
Linear Hypotheses:
Estimate lwr upr
M - L == 0 -10. -22.0764 2.0764
H - L == 0 -14.7222 -26.7986 -2.6458
H - M == 0 -4.7222 -16.7986 7.3542
99% family-wise confidence level
To see how to extract the p-values, str() is useful:
str(summary(mc))
summary(mc)$test$pvalues
[1] 0.038569618 0.001440005 0.463049425
attr(,"error")
[1] 0.000226511
> Thanks for your help,
> Cheers../Murli
>
>
>
>
> From: talepanda [mailto:[EMAIL PROTECTED]
> Sent: Sun 3/18/2007 11:07 PM
> To: Nair, Murlidharan T
> Cc: r-help@stat.math.ethz.ch
> Subject: Re: [R] multcomp
>
>
>
> ?glht says
>> with 'print', 'summary', 'confint', 'coef' and 'vcov' methods
>> being available.
>
> try:
>
> example(glht)
> summary(glht(amod, linfct = mcp(tension = "Tukey")))
> confint(glht(amod, linfct = mcp(tension = "Tukey")))
>
> On 3/19/07, Nair, Murlidharan T <[EMAIL PROTECTED]> wrote:
>> I used the multcomp package sometime back for doing multiple
>> comparisons. I see that it has been updated and the methods like simint
>> are no longer supported. When I run the program it prompts to me to use
>> glht. How do I get the lower and upper conf int and the pValues using
>> glht? Does anyone have an example?
>>
>> Thanks ../Murli
>>
>>
>>
>>
>> [[alternative HTML version deleted]]
>>
>> __
>> R-help@stat.math.ethz.ch mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] order of values in vector
Tord Snäll wrote: > Dear all, > I would like to get the order of the values in a vector. I have tried > rank(), order() and searched the archive, though without success. > > Here is an example of a try > x= c(20,30,50,40,60,10) > cbind(sort.list(x),x) > x > [1,] 6 20 > [2,] 1 30 > [3,] 2 50 > [4,] 4 40 > [5,] 3 60 > [6,] 5 10 > but I was hoping to get this: > x > [1,] 2 20 > [2,] 3 30 > [3,] 5 50 > [4,] 4 40 > [5,] 6 60 > [6,] 1 10 > > I'm most grateful for a tip! cbind(rank(x), x) x [1,] 2 20 [2,] 3 30 [3,] 5 50 [4,] 4 40 [5,] 6 60 [6,] 1 10 > cheers, > Tord -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] dataframe layout
Robert Baer wrote: > Can someone remind me how to change the columns in df.a into a two column > df.b that contains one column of data and another column of the original > column headings as levels. > > Example: > a=1:3 > b=4:6 > c=7:9 > df.a=data.frame(a,b,c) > > Should become in df.b: > dat lev > 1 a > 2 a > 3 a > 4 b > 5 b > 6 b > 7 c > 8 c > 9 c Here are a couple of different approaches: df.b <- data.frame(dat = unlist(df.a), lev = rep(names(df.a), each = dim(df.a)[1])) df.b dat lev a1 1 a a2 2 a a3 3 a b1 4 b b2 5 b b3 6 b c1 7 c c2 8 c c3 9 c library(reshape) melt(df.a, measure.var = names(df.a), variable_name = "lev") lev value 1 a 1 2 a 2 3 a 3 4 b 4 5 b 5 6 b 6 7 c 7 8 c 8 9 c 9 > Thanks. > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] about bootstrapping
Ray Cheung wrote: > Dear All, > > I've a 10 by 5 data frame like this > > set.seed(1001) > a <- rnorm(10) > b <- rnorm(10) > c <- rnorm(10) > d <- rnorm(10) > e <- rnorm(10) > A <- cbind(a,b,c,d,e) > > Each row is one datum. I want to resample within A[,5]. Fit the regression > lines lm(A[,5] ~ A[,1] + A[,2]) and lm(A[,5] ~ A[,3] + A[,4]) to each > bootstrap sample. I don't know how to write the statistics part. Can anybody > help me? Thank you very much. You might have a look at John Fox's Appendix to An R and S-PLUS Companion to Applied Regression on this topic: http://cran.r-project.org/doc/contrib/Fox-Companion/appendix-bootstrapping.pdf > Regards, > > Ray > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] selecting rows with more than x occurrences in a given column (data type is names)
Mike Jasper wrote:
> Despite a long search on the archives, I couldn't find how to do this.
> Thanks in advance for what is likely a simple issue.
>
> I have a data set where the first column is name (i.e., 'Joe Smith',
> 'Jane Doe', etc). The following columns are data associated with that
> person. I have many people with multiple rows. What I want is to get a
> new data frame out with only the people who have more than x
> occurrences in the first column.
>
> Here's what I've done, that's not working:
>
> Let's call my old data.frame "all.data"
>
> table(all.data$names)>10
>
> I get a list of names and TRUE/FALSE values. I then want to make a
> list of the TRUEs and pass that to some subset type command like
>
> dup.names=table(all.data$names)>10
>
> new.data=(all.data[all.data$names==dup.names,])
>
> That's not working because the dimensions are wrong (I think). But
> even when I tried to do part of it manually (to troubleshoot) like
> this
>
> dup.names=c('Joe Smith','Jane Doe','etc')
>
> I got warnings and it didn't work correctly. There must be a simple
> way to do this that I'm just not seeing. Thanks.
Does this help?
df <- data.frame(PERSON = rep(c("John","Tom","Sara","Mary"),
c(5,4,5,4)),
Y = runif(18))
subset(df, PERSON %in% names(which(table(PERSON) >= 5)))
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] where can I find Durbin-Waston test tables for Confidence Level 2.5% or 0.5%?
Michael wrote:
> Hi all,
>
> I am doing a two-sided DW test:
>
> H0: rho = 0
> H1: rho =/= 0
>
> My understanding is that most test statistics tables are one-sided. It's the
> way they created the table.
>
> So from online, by doing Googling, I found a bunch of DW tables for
> Confidence Level 5%.
>
> Those tables can answer my two-sided question at 5x2 = 10% confidence level.
>
> But what if I want two-sided test at 1% and 5% confidence level?
>
> Then I need 0.5% and 2.5% tables on those one-sided table.
>
> My sample size is 278, and the number of parameters is 2, these adds to the
> hardship of finding a good table...
>
> Could anybody give me some pointers of two-sided DW tables or 1-sided DW
> table with 0.5% and 2.5% confidence levels?
>
> Thanks a lot!
>
> Moreover, I appreciate any pointers about electronic tables that I can use
> in programs, I want to implement DW test myself, but let the program
> searching a table automatically...
Are you aware of the implementations of this test in the car, lmtest,
and fMultivar packages?
RSiteSearch("Durbin-Watson", restrict="function") finds those functions.
> Thanks a lot!
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Identifying points in a plot that have duplicate values
David Lloyd wrote:
> I have code like this: -
>
> #---
> --
>
> x=scan()
> 0 0 0 0 0 1 2 3 4
>
> y=scan()
> 1 1 1 2 2 1 3 4 5
>
> plot(x,y)
>
> identify(0,1,3) #Allows me to select manually to identify co-ordinate
> (0,1) as being duplicated 3 times
> identify(0,2,2) #Allows me to select manually to identify co-ordinate
> (0,2) as being duplicated 2 times
> #---
> --
>
> Is there not a way I can automatically display if points are duplicated
> and by how many times?
>
> I thought if I 'jittered' the points ever so slightly I could get an
> idea of how many duplicates there are but with >100 points the graph
> looks very messy.
You might consider using alpha transparency - the more times a point
is duplicated the darker it will be. For example:
df <- data.frame(x=c(0, 0, 0, 0, 0, 1, 2, 3, 4),
y=c(1, 1, 1, 2, 2, 1, 3, 4, 5))
pdf("alphaExample.pdf", version = "1.4", width = 6, height = 6)
with(df, plot(x,y, col=rgb(1,0,0,.3), pch=16))
dev.off()
RSiteSearch("alpha transparency")
> Regards
> DaveL
>
>
>
>
>
>
>
>
> Click for free info on getting an MBA and make $200K/ year
>
>
>
>
> Need cash? Click to get a payday loan
> <http://tagline.bidsystem.com/fc/CAaCDCZ60nyjrrOboFeUJgRjigwgNftK/>
>
>
>
> style="font-size:13.5px">___Get
> the Free email that has everyone talking at href=http://www.mail2world.com target=new>http://www.mail2world.com
> Unlimited Email Storage – POP3 – Calendar
> – SMS – Translator – Much More!
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] format of summary.lm for 2-way ANOVA
Patrick Drechsler wrote: > Hi, > > I am performing a two-way ANOVA (2 factors with 4 and 5 levels, > respectively). If I'm interpreting the output of summary correctly, > then the interaction between both factors is significant: > > , > | ## Two-way ANOVA with possible interaction: > | > model1 <- aov(log(y) ~ xForce*xVel, data=mydataset) > | > | > summary(model1) > | Df Sum Sq Mean Sq F valuePr(>F) > | xForce3 16.640 5.547 19.0191 1.708e-11 *** > | xVel 4 96.391 24.098 82.6312 < 2.2e-16 *** > | xForce:xVel 12 10.037 0.836 2.8681 0.0008528 *** > | Residuals 371 108.194 0.292 > | --- > | Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 > | 3 observations deleted due to missingness > ` > > To see the interactions in detail I call summary.lm: > > , > | > summary.lm(model1) > | > | Call: > | aov(formula = log(y) ~ xForce * xVel, data = mydataset) > | > | Residuals: > | Min 1Q Median 3Q Max > | -2.04830 -0.32420 -0.04653 0.34928 1.46755 > | > | Coefficients: > | Estimate Std. Error t value Pr(>|t|) > | (Intercept) -1.663406 0.027335 -60.853 < 2e-16 *** > | xForce.L 0.408977 0.054726 7.473 5.68e-13 *** > | xForce.Q 0.101240 0.054670 1.852 0.0648 . > | xForce.C-0.068068 0.054613 -1.246 0.2134 > | xVel.L 1.079042 0.061859 17.444 < 2e-16 *** > | xVel.Q 0.339802 0.061439 5.531 6.03e-08 *** > | xVel.C 0.015422 0.060751 0.254 0.7997 > | xVel^4 -0.044399 0.060430 -0.735 0.4630 > | xForce.L:xVel.L 0.622060 0.123966 5.018 8.12e-07 *** > | xForce.Q:xVel.L 0.034298 0.123718 0.277 0.7818 > | xForce.C:xVel.L -0.114776 0.123470 -0.930 0.3532 > | xForce.L:xVel.Q 0.309293 0.123057 2.513 0.0124 * > | xForce.Q:xVel.Q 0.054798 0.122879 0.446 0.6559 > | xForce.C:xVel.Q -0.144219 0.122700 -1.175 0.2406 > | xForce.L:xVel.C 0.110588 0.121565 0.910 0.3636 > | xForce.Q:xVel.C -0.001929 0.121502 -0.016 0.9873 > | xForce.C:xVel.C -0.039477 0.121438 -0.325 0.7453 > | xForce.L:xVel^4 0.090491 0.120870 0.749 0.4545 > | xForce.Q:xVel^4 -0.002762 0.120861 -0.023 0.9818 > | xForce.C:xVel^4 -0.028836 0.120852 -0.239 0.8115 > | --- > | Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 > | > | Residual standard error: 0.54 on 371 degrees of freedom > | (3 observations deleted due to missingness) > | Multiple R-Squared: 0.5322, Adjusted R-squared: 0.5082 > | F-statistic: 22.21 on 19 and 371 DF, p-value: < 2.2e-16 > ` > > I am wondering what the logic is behind the formatting of > rownames. What do the strings "L", "Q", "C" and "^4" mean? The default contrast for ordered factors is contr.poly(), and "L", "Q", "C" and "^4" refer to linear, quadratic, cubic, and quartic, respectively. You might look at some plots if you have not already. For example: with(mydataset, interaction.plot(xForce, xVel, log(y))) library(lattice) bwplot(log(y) ~ xForce | xVel, data = mydataset) > Apologies in case I missed something obvious in the relevant > documentations/archives. > > Thanks for any pointers, > > Patrick > > PS: Here are some details about the dataset: > > , > | > summary(mydataset) > |yxForce xVel > | Min. :0.03662 0.01:97 10 :79 > | 1st Qu.:0.10376 0.1 :98 50 :80 > | Median :0.16314 1 :98 100 :80 > | Mean :0.26592 2 :98 500 :80 > | 3rd Qu.:0.28077 NA's: 3 5000:72 > | Max. :2.39490 NA's: 3 > | NA's :3.0 > | > | > str(mydataset) > | 'data.frame': 394 obs. of 3 variables: > | $ y : num 0.167 0.158 0.152 0.158 0.131 ... > | $ xForce: Ord.factor w/ 4 levels "0.01"<"0.1"<"1"<..: 1 2 3 4 1 2 3 4 1 2 > ... > | $ xVel : Ord.factor w/ 5 levels "10"<"50"<"100"<..: 3 3 3 3 1 1 1 1 2 2 > ... > ` > > I am using > > platform i486-pc-linux-gnu > arch i486 > os linux-gnu > system i486, linux-gnu > status > major 2 > minor 4.1 > year 2006 > month 12
Re: [R] sampling random groups with all observations in the group
Chuck Cleland wrote: > Wadud, Zia wrote: >> Hi >> I have a panel dataset with large number of groups and differing number >> of observations for each group. I want to randomly select say, 20% of >> the groups or 200 groups, but along with all observations from the >> selcted groups (with the corresponding data). >> I guess it is possible to generate a random sample from the groups ids >> and then match that with the entire dataset to have the intended >> dataset, but it sounds cumbersome and possibly there is an easier way to >> do this? checked the package 'sampling' or command 'sample', but they >> cant do exactly the same thing. >> I was wondering if someone on this list will be able to share his/her >> knowldege? > > How about something like this? > > df <- data.frame(GROUP = rep(1:5, c(2,3,4,2,2)), Y = runif(13)) > > # Sample Two of the Five Groups > > subset(df, GROUP %in% with(df, sample(unique(GROUP), 2))) The with() part can be dropped too. subset(df, GROUP %in% sample(unique(GROUP), 2)) >> Thanks in advance, >> Zia >> ** >> Zia Wadud >> PhD Student >> Centre for Transport Studies >> Department of Civil and Environmental Engineering >> Imperial College London >> London SW7 2AZ >> Tel +44 (0) 207 594 6055 >> >> >> [[alternative HTML version deleted]] >> >> __ >> R-help@stat.math.ethz.ch mailing list >> https://stat.ethz.ch/mailman/listinfo/r-help >> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html >> and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] sampling random groups with all observations in the group
Wadud, Zia wrote: > Hi > I have a panel dataset with large number of groups and differing number > of observations for each group. I want to randomly select say, 20% of > the groups or 200 groups, but along with all observations from the > selcted groups (with the corresponding data). > I guess it is possible to generate a random sample from the groups ids > and then match that with the entire dataset to have the intended > dataset, but it sounds cumbersome and possibly there is an easier way to > do this? checked the package 'sampling' or command 'sample', but they > cant do exactly the same thing. > I was wondering if someone on this list will be able to share his/her > knowldege? How about something like this? df <- data.frame(GROUP = rep(1:5, c(2,3,4,2,2)), Y = runif(13)) # Sample Two of the Five Groups subset(df, GROUP %in% with(df, sample(unique(GROUP), 2))) > Thanks in advance, > Zia > ** > Zia Wadud > PhD Student > Centre for Transport Studies > Department of Civil and Environmental Engineering > Imperial College London > London SW7 2AZ > Tel +44 (0) 207 594 6055 > > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] count the # of appearances...
bunny , lautloscrew.com wrote: > Hi there, > > is there a possibility to count the number of appearances of an > element in a vector ? > i mean of any given element.. deliver all elements which are exactly > xtimes in this vector ? > > thx in advance !! > myvec <- c(0, rep(1:6, each = 2), 7, 7, 7) > table(myvec) myvec 0 1 2 3 4 5 6 7 1 2 2 2 2 2 2 3 # Extract elements appearing twice > myvec[myvec %in% names(which(table(myvec) == 2))] [1] 1 1 2 2 3 3 4 4 5 5 6 6 > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to apply the function cut( ) to many columns in a data.frame?
ahimsa campos-arceiz wrote:
> Dear useRs,
>
> In a data.frame (df) I have several columns (x1, x2, x3xn) containing
> data as a continuous numerical response:
>
> df
> var x1x2 x3
> 1143 147 137
> 2 9393 117
> 316439 101
> 4123 11897
> 5 63 125 97
> 612983 124
> 712393 136
> 812380 79
> 9 89 107 150
> 10 7895121
>
> I want to classify the values in the columns x1, x2, etc, into bins of fix
> margins (0-5, 5-10, ). For one vector I can do it easily with the
> function cut:
>
>> df$x1 <- cut(df$x1, br=5*(0:40), labels=5*(1:40))
>> df$x1
> [1] 145 95 165 125 65 130 125 125 90 80
> 40 Levels: 5 10 15 20 25 30 35 40 45 50 55 60 65 70 75 80 85 90 95 100 ...
> 200
>
> However if I try to use a subset of my data.frame:
>
> df[,3:4] <- cut(df[,3:4], br=5*(0:40), labels=5*(1:40))
>
> Error in cut.default(df[, 3:4], br = 5 * (0:40), labels = 5 * (1:40)) :
> 'x' must be numeric
>
>
> How can I make this work with data frames in which I want to apply the
> function cut( ) to many columns in a data.frame?
You have an answer within your question - use one of the various
"apply" functions. For example:
lapply(df[,3:4], function(x){cut(x, br=5*(0:40), labels=5*(1:40))})
?lapply
?sapply
?apply
> I guess that I might have to use something like for ( ) (which I'm not
> familiar with), but maybe you know a straight forward method to use with
> data.frames.
>
>
> Thanks a lot!
>
> Ahimsa
>
> *
>
> # data
> var <- 1:10
> x1 <- rnorm(10, mean=100, sd=25)
> x2 <- rnorm(10, mean=100, sd=25)
> x3 <- rnorm(10, mean=100, sd=25)
> df <- data.frame(var,x1,x2,x3)
> df
>
> # classifying the values of the vector df$x1 into bins of width 5
> df$x1 <- cut(df$x1, br=5*(0:40), labels=5*(1:40))
> df$x1
>
> # trying it a subset of the data.frame
> df[,3:4] <- cut(df[,3:4], br=5*(0:40), labels=5*(1:40))
> df[,3:4]
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Adding duplicates by rows
Serguei Kaniovski wrote:
> Hi,
>
> I am trying to add duplicates of matrix "mat" by row. Commands
>
> subset(mat,duplicated(rownames(mat)))
>
> or
>
> mat[which(duplicated(rownames(mat))),]
>
> return only half of the required indices. How can I find the remaining
> ones, ie the matches, so that I can add them up?
mat <- matrix(runif(70), ncol=5)
rownames(mat) <- c("Z", rep(LETTERS[1:6], each=2), "G")
There is probably a more elegant way, but this seems to do what you want:
mat[rownames(mat) %in% names(which(table(rownames(mat)) > 1)),]
Also, have you considered aggregate()?
aggregate(mat, list(ROW = rownames(mat)), sum)
> Thanks,
> Serguei
>
> ___
>
> Austrian Institute of Economic Research (WIFO)
>
> Name: Serguei Kaniovski P.O.Box 91
> Tel.: +43-1-7982601-231 Arsenal Objekt 20
> Fax: +43-1-7989386 1103 Vienna, Austria
> Mail: [EMAIL PROTECTED] A-1030 Wien
>
> http://www.wifo.ac.at/Serguei.Kaniovski
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] rotating labels
Gitta Lubke wrote:
> I'd like to rotate x-axis labels to get a diagonal orientation (rather
> than horizontal or vertical), this would result in an easier read, any
> ideas?
> thanks, Gitta
This is FAQ 7.27
http://cran.r-project.org/doc/FAQ/R-FAQ.html#How-can-I-create-rotated-axis-labels_003f
And for more examples in the archives, RSiteSearch("rotate label")
>
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] TRUE/FALSE as numeric values
Thomas Preuth wrote: > Hello, > > I want to select in a column of a dataframe all numbers smaller than a > value x > but when I type in test<-(RSF_EU$AREA<=x) I receiv as answer: > > test > [1] TRUE FALSE FALSE TRUE TRUE TRUE FALSE FALSE TRUE TRUE TRUE > FALSE TRUE TRUE TRUE TRUE TRUE > [18] TRUE TRUE TRUE TRUE FALSE FALSE TRUE TRUE TRUE TRUE TRUE > TRUE TRUE FALSE TRUE TRUE TRUE > [35] FALSE TRUE TRUE TRUE TRUE FALSE FALSE FALSE TRUE TRUE FALSE > TRUE TRUE FALSE FALSE TRUE FALSE > [52] TRUE TRUE TRUE TRUE TRUE FALSE TRUE TRUE FALSE TRUE > > How can i get the values smaller than x and not the TRUE/FALSE reply? This is covered in http://cran.r-project.org/doc/manuals/R-intro.html (Section 2.7; Index vectors; selecting and modifying subsets of a data set) RSF_EU$AREA[RSF_EU$AREA <= x] or to retain other columns in the data frame subset(RSF_EU, AREA <= x) > Thanks in advance, > Thomas > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Principal Component Analysis & explained variance
Milton Cezar Ribeiro wrote: > Hi there, > > How can I know the explaned variance of a PC axis generated by prcomp()? From the standard deviations of each component, you could do something like this maybe: prcomp(USArrests, scale = TRUE)$sdev^2 / ncol(USArrests) [1] 0.62006039 0.24744129 0.08914080 0.04335752 > Kind regards, > > miltinho > Brazil > > __ > > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] bootstrapping Levene's test
[EMAIL PROTECTED] wrote:
> Hello all,
>
> I am low down on the learning curve of R but so far I have had little
> trouble using most of the packages. However, recently I have run into
> a wall when it comes to bootstrapping a Levene's test (from the car
> package) and thought you might be able to help. I have not been able
> to find R examples for the "boot" package where the test statistic
> specifically uses a grouping variable (or at least a simple example
> with this condition). I would like to do a non-parametric bootstrap
> to eventually get 95% confidence intervals using the boot.ci command.
> I have included the coding I have tried on a simple data set below. If
> anyone could provide some help, specifically with regards to how the
> "statistic" arguement should be set up in the boot package, it would
> be greatly appreciated.
>
>> library(boot)
>> library(car)
>> data<-c(2,45,555,1,77,1,2,1,2,1)
>> group<-c(1,1,1,1,1,2,2,2,2,2)
>> levene.test(data,group)
> Levene's Test for Homogeneity of Variance
>Df F value Pr(>F)
> group 1 1.6929 0.2294
> 8
>> stat<-function(a){levene.test(a,group)}
>> trial1<-boot(data,statistic,100)
> Error in statistic(data, original, ...) : unused argument(s) ( ...)
One problem is that levene.test() returns an ANOVA table, which is not
a statistic. Looking inside levene.test() to see what might be used as
a statistic, for the two-group case it seems like you could do something
like this:
library(boot)
library(car)
set.seed(671969)
df <- data.frame(Y = c(rnorm(100,0,1), rnorm(100,0,1)),
G = rep(c(1,2), each = 100))
boot.levene <- function(data,indices){
levene.diff <- diff(tapply(df$Y[indices],
list(df$G[indices]), mad))
return(levene.diff)
}
first.try <- boot(df, boot.levene, 1000)
first.try
ORDINARY NONPARAMETRIC BOOTSTRAP
Call:
boot(data = df, statistic = boot.levene, R = 1000)
Bootstrap Statistics :
original biasstd. error
t1* -0.1944924 0.02439172 0.1753512
boot.ci(first.try, index=1, type="bca")
BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
Based on 1000 bootstrap replicates
CALL :
boot.ci(boot.out = first.try, type = "bca", index = 1)
Intervals :
Level BCa
95% (-0.5772, 0.1159 )
Calculations and Intervals on Original Scale
> Best regards,
> Kevin
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Data frame: how to create list of row max?
Johannes Graumann wrote:
> Dear all,
>
> Can anyone please shed some light onto how to do this?
>
> This will give me all "intensity" columsn in my data frame:
> intensityindeces <- grep("^Intensity",names(dataframe),value=TRUE)
>
> This will give me the maximum intensity for the first row:
> intensityone <- max(dataframe[1,intensityindeces])
>
> What I'm now looking for is how to dfo this for the whole data frame. Should
> yield a list of maximum intensities of all rows.
> Can't figure it out ... please nudge me where I need to go.
If you want the values themselves:
apply(dataframe[,intensityindeces], 1, max)
If you want the column in which the max appears:
apply(dataframe[,intensityindeces], 1, which.max)
?apply
> Thanks, Joh
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] How is the line in av.plot calculated?
Tom La Bone wrote: > I suspect that the line in the added variable plot (library car) is a SLR of > the residuals, but I can't seem to find this written anywhere. Can someone > confirm this? Thanks. See the Fox(1997) reference on the help page (pages 281-283). Note that both the x-axis and the y-axis are residuals. For example, you could do the av.plot for income in the Duncan data as follows: library(car) res1 <- residuals(lm(prestige ~ education+type, data=Duncan)) res2 <- residuals(lm(income ~ education+type, data=Duncan)) plot(res2, res1) abline(lm(res1 ~ res2)) > Tom > > > > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] R implementations of scatterplot/map labeling algorithims?
Michael Friendly wrote: > Dear R community > > In a current paper, I'm (briefly) considering the topic of producing > scatterplots or maps with point labels positioned in such a way as to > minimize label overlap and occlusion. This is a topic with a large, but > scattered literature. In CS, it is considered NP-hard, but there are > a variety of approximate solutions. The most complete bibliography I've > found is > the Map-Labeling Bibliography, > http://liinwww.ira.uka.de/bibliography/Theory/map.labeling.html > > AFAIK, the only concrete and published implementation is a Fortran > program published by Noma (below), and then adapted by Warren Kuhfeld > at SAS in PROC PLOT, and used in his %plotit macro. > > Journal Title - Psychometrika > Article Title - Heuristic method for label placement in scatterplots > Volume - Volume 52 > Issue - 3 > First Page - 463 > Last Page - 468 > Issue Cover Date - 1987-09-27 > Author - Elliot Noma > DOI - 10.1007/BF02294366 > Link - http://www.springerlink.com/content/c4k6205r83156165 > > Does anyone know of a related R (or other) public implementation of a > solution to this problem? One thing that comes to mind is labcurve() in Frank Harrell's Hmisc package. "Optionally draws a set of curves then labels the curves. A variety of methods for drawing labels are implemented, ranging from positioning using the mouse to automatic labeling to automatic placement of key symbols with manual placement of key legends to automatic placement of legends. For automatic positioning of labels or keys, a curve is labeled at a point that is maximally separated from all of the other curves. Gaps occurring when curves do not start or end at the same x-coordinates are given preference for positioning labels. If labels are offset from the curves (the default behaviour), if the closest curve to curve i is above curve i, curve i is labeled below its line. If the closest curve is below curve i, curve i is labeled above its line. These directions are reversed if the resulting labels would appear outside the plot region." > thanks, > -Michael -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] missing -> nonmissing levels
Jon Minton wrote:
> Hi,
>
> I expect this is simple but haven’t found an answer looking on the
> archives...
>
> I want to convert ‘NA’ (missing) to particular levels (nonmissing) in factor
> vectors.
>
> e.g. I know
>
>> X <- c(1, 2, 3)
>
>> summary(X)
>
>Min. 1st Qu. MedianMean 3rd Qu.Max.
>
> 1.0 1.5 2.0 2.0 2.5 3.0
>
>> X <- as.factor(X)
>
>> summary(X)
>
> 1 2 3
>
> 1 1 1
>
>> levels(X)
>
> [1] "1" "2" "3"
>
>> levels(X) <- c("A", NA, "B")
>
>> summary(X)
>
>AB NA's
>
>111
>
> But what if I want to turn the NA back into a level?
>
> How do I do this?
One way is to use recode() in the car package by John Fox. For example:
library(car)
X <- factor(c(1, NA, 3))
X
[1] 1 3
Levels: 1 3
recode(X, "NA='J'")
[1] 1 J 3
Levels: 1 3 J
> Thanks, Jon
>
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Generating MVN Data
Rauf Ahmad wrote: > Dear All > > I want to generate multivariate normal data in R for a given covariance > matrix, i.e. my generated data must have the given covariance matrix. I > know the rmvnorm command is to be used but may be I am failing to > properly assign the covariance matrix. > > Any help will be greatly appreciated library(MASS) mat <- mvrnorm(100, mu=rep(0,4), Sigma = diag(4), empirical=TRUE) cov(mat) [,1] [,2] [,3] [,4] [1,] 1.00e+00 -1.634448e-16 -1.004223e-16 -2.015521e-16 [2,] -1.634448e-16 1.00e+00 4.244391e-17 1.544399e-17 [3,] -1.004223e-16 4.244391e-17 1.00e+00 -1.921951e-16 [4,] -2.015521e-16 1.544399e-17 -1.921951e-16 1.00e+00 ?mvrnorm > thanks. > > M. R. Ahmad > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Help in using multcomp.
A Ezhil wrote:
> Hi Chuck,
>
> Thank you very for this help. I am able to store the
> results. Now, I am facing the following problems:
>
> 1. I am trying to extract only the p-values from
>
> dunres <- lapply(amod, function(x){summary(glht(x,
> linfct=mcp(f = contr)))})
>
> but I am stuck. If I use dunres[[1]], it displays the
> results. But I don't know how to extract the p value
> from this.
Ezhil:
str(summary(glht(amod[[1]], linfct=mcp(f = contr suggests a way to
extract just the p values. Try this:
library(multcomp)
dat <- matrix(rnorm(45), nrow=5, ncol=9)
f <- gl(3,3,9, label=c("C", "Tl", "T2"))
aof <- function(x) {
m <- data.frame(f, x);
aov(x ~ f, m)
}
amod <- apply(dat,1,aof)
my.pvals <- sapply(amod, function(x){summary(glht(x, linfct=mcp(f =
contr)))$test$pvalues})
rownames(my.pvals) <- rownames(contr)
colnames(my.pvals) <- paste("amod", 1:5, sep="")
my.pvals
amod1 amod2 amod3 amod4 amod5
C - T10.03348242 0.3581771 0.9873633 0.9764219 0.9225445
C - T20.67794496 0.7138491 0.2183949 0.9962458 0.8439224
C - All T 0.10621039 0.4344881 0.4913578 0.9970689 0.8519888
> 2. If I want to get raw pvalues instead of adjusted
> ones, what should I do in summary(glht)?
See the test argument of summary.glht(). You could do something like
this:
my.pvals <- sapply(amod, function(x){summary(glht(x, linfct=mcp(f =
contr)), test = adjusted("none"))$test$pvalues})
rownames(my.pvals) <- rownames(contr)
colnames(my.pvals) <- paste("amod", 1:5, sep="")
my.pvals
amod1 amod2 amod3 amod4 amod5
C - T10.01770489 0.2197886 0.8985143 0.8611198 0.740
C - T20.46651736 0.4994115 0.1278132 0.9448342 0.6364631
C - All T 0.05953715 0.2733017 0.3149223 0.9512189 0.6464735
hope this helps,
Chuck Cleland
> Thanks again for your help. I look forward to your
> reply.
>
> Kind regards,
> Ezhil
>
>
> --- Chuck Cleland <[EMAIL PROTECTED]> wrote:
>
>> A Ezhil wrote:
>>> Hi All,
>>>
>>> I am trying use 'multcomp' for multiple
>> comparisons
>>> after my ANOVA analysis. I have used the following
>>> code to do ANOVA:
>>>
>>> dat <- matrix(rnorm(45), nrow=5, ncol=9)
>>> f <- gl(3,3,9, label=c("C", "Tl", "T2"))
>>>
>>> aof <- function(x) {
>>> m <- data.frame(f, x);
>>> aov(x ~ f, m)
>>> }
>>> amod <- apply(dat,1,aof)
>>>
>>> Now, how can I use 'glht' for the above amod. I
>> know
>>> that I cannot use simply
>>>
>>> glht(amod, linfct = mcp(f = "Dunnett")).
>> Since amod is a list of models rather than one
>> model, do you want
>> something like this?
>>
>> lapply(amod, function(x){summary(glht(x, linfct =
>> mcp(f = "Dunnett")))})
>>
>>> Also, if I want to use Dunnett for comparing C vs
>> (T1
>>> and T2), how can I specify this in the glht
>> function.
>>
>> How about doing it with user-defined contrasts?
>>
>> contr <- rbind("C - T1 " = c(-1, 1, 0),
>>"C - T2 " = c(-1, 0, 1),
>>"C - All T" = c(-1,.5,.5))
>>
>> lapply(amod, function(x){summary(glht(x, linfct =
>> mcp(f = contr)))})
>>
>>> Thanks in advance.
>>> Regards,
>>> Ezhil
>>>
>>> __________
>>> R-help@stat.math.ethz.ch mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>>> and provide commented, minimal, self-contained,
>> reproducible code.
>>
>> --
>> Chuck Cleland, Ph.D.
>> NDRI, Inc.
>> 71 West 23rd Street, 8th floor
>> New York, NY 10010
>> tel: (212) 845-4495 (Tu, Th)
>> tel: (732) 512-0171 (M, W, F)
>> fax: (917) 438-0894
>>
>
>
>
>
>
> Need a quick answer? Get one in minutes from people who know.
> Ask your question on www.Answers.yahoo.com
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
