Re: [R] (coxph, se) Obtaining standard errors of coefficients from coxph to store
David, It would be helpful to give an example of what you would like to extract. I guess you know how to extract elements from vectors and lists. However, sometimes the objects returned by functions can be rather complex (output of coxph() is...) A general method to capture printed output is via capture.output(). Maybe not fast, but if you have no other solution... Joris a - rnorm(10,1,1) b - rnorm(10,1,1) mod - lm(a~b) smod - summary(mod) smod Call: lm(formula = a ~ b) Residuals: Min 1Q Median 3Q Max -1.7482 -0.5991 0.1211 0.8341 1.4975 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) 1.6210 0.5332 3.040 0.0161 * b-0.7667 0.5037 -1.522 0.1664 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 1.142 on 8 degrees of freedom Multiple R-Squared: 0.2246, Adjusted R-squared: 0.1277 F-statistic: 2.317 on 1 and 8 DF, p-value: 0.1664 output - capture.output(print(smod)) output [1] [2] Call: [3] lm(formula = a ~ b) [4] [5] Residuals: [6] Min 1Q Median 3Q Max [7] -1.7482 -0.5991 0.1211 0.8341 1.4975 [8] [9] Coefficients: [10] Estimate Std. Error t value Pr(|t|) [11] (Intercept) 1.6210 0.5332 3.040 0.0161 * [12] b-0.7667 0.5037 -1.522 0.1664 [13] --- [14] Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 [15] [16] Residual standard error: 1.142 on 8 degrees of freedom [17] Multiple R-Squared: 0.2246,\tAdjusted R-squared: 0.1277 [18] F-statistic: 2.317 on 1 and 8 DF, p-value: 0.1664 [19] David Lloyd [EMAIL PROTECTED] lloyd.com To Sent by: r-help@stat.math.ethz.ch [EMAIL PROTECTED] cc at.math.ethz.ch Subject [R] (coxph, se) Obtaining standard 16/08/2007 11:31 errors of coefficients from coxph to store Hi all, I'm wanting to be able to find and store the z-score of coxph below: - modz=coxph(Surv(TSURV,STATUS)~RAGE+DAGE+REG_WTIME_M+CLD_ISCH+POLY_VS, data=kidneyT,method=breslow) I know summary(modz) will give me this, but how do i extract the standard error or z-score values in a similar way to obtaining the coefficients by coef(modz) ? I think it must be something to do with modz$var but I'm having a complete mental blank. I need this info so I can write a function to use within a bootstrap so I can record the number of times (proportion) each variable in the Cox PH model is actually significant over all the bootstrap resamples. Any assistance is greatly appreciated DL Click to find local singles for dating, romance and fun. http://tagline.bidsystem.com/fc/Ioyw36XJJVs581mfqGSywy0Z69Mq8VM03oVytPu 8otqP84CBZmNX2G/ span id=m2wTlpfont face=Arial, Helvetica, sans-serif size=2 style=font-size:13.5 px___BRGet the Free email that has everyone talking at a href=http://www.mail2world.com target=newhttp://www.mail2world.com/abr font color=#99Unlimited Email Storage #150; POP3 #150; Calendar #150; SMS #150; Translator #150; Much More!/font/font/span [[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.
Re: [R] Custom axis
x - 1:10 y - rnorm(10,10,1) x2 - 3*x + 2 plot(y ~ x, xaxt = n) axis(side=1,at = x, labels = x2) Joris Florent Bresson [EMAIL PROTECTED] frTo Sent by: r-help@stat.math.ethz.ch [EMAIL PROTECTED] cc at.math.ethz.ch Subject [R] Custom axis 01/08/2007 11:49 Dear R users, I would like to draw a plot with a custom scale for the axis. More precisely, instead of plotting y on x, I want to plot y on a monotone function of x (for instance a*x+b). Which command and/or package should I use in order to get this result? Thanks Florent Bresson _ __ 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.
[R] Re : Custom axis
What about the other way round? x - 1:10 y - rnorm(10,10,1) x2 - 3*x + 2 plot(y ~ x2, xaxt = n) axis(side=1,at = x2, labels = x) Florent Bresson [EMAIL PROTECTED] frTo [EMAIL PROTECTED] 01/08/2007 12:14 cc r-help@stat.math.ethz.ch, [EMAIL PROTECTED] Subject Re : [R] Custom axis Maybe I do not explain well what I would like to do. I do not want to change the labels of the axis, but the scale. What I want is a general procedure for changing the scale. Its like using a logarithmic scale on a plot. Labels are the same, but the increases of x along the x-axis are defined by a known monotone and continuous function. Florent Bresson - Message d'origine De : [EMAIL PROTECTED] [EMAIL PROTECTED] À : Florent Bresson [EMAIL PROTECTED] Cc : r-help@stat.math.ethz.ch; [EMAIL PROTECTED] Envoyé le : Mercredi, 1 Août 2007, 12h01mn 58s Objet : Re: [R] Custom axis x - 1:10 y - rnorm(10,10,1) x2 - 3*x + 2 plot(y ~ x, xaxt = n) axis(side=1,at = x, labels = x2) Joris Florent Bresson [EMAIL PROTECTED] frTo Sent by: r-help@stat.math.ethz.ch [EMAIL PROTECTED] cc at.math.ethz.ch Subject [R] Custom axis 01/08/2007 11:49 Dear R users, I would like to draw a plot with a custom scale for the axis. More precisely, instead of plotting y on x, I want to plot y on a monotone function of x (for instance a*x+b). Which command and/or package should I use in order to get this result? Thanks Florent Bresson _ __ 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. _ Ne gardez plus qu'une seule adresse mail ! Copiez vos mails vers Yahoo! Mail __ 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] lme or gls prediction intervals
Martin, Have you checked ?intervals.gls This intervals are approximate, but this would be the obvious starting point to me. Joris Martin Henry H. Stevens [EMAIL PROTECTED] To edu R-Help r-help@stat.math.ethz.ch Sent by: cc [EMAIL PROTECTED] at.math.ethz.ch Subject [R] lme or gls prediction intervals 24/07/2007 19:22 Hi folks, I am trying to generate 95% confidence intervals for a gls model using predict.nlme with R version 2.5.1 (2007-06-27) . nlme: Linear and Nonlinear Mixed Effects Models. R package version 3.1-83. I have looked in help, and I can do it for lm and glm models, and I can generate simple predictions for lme models with various levels -- I am familiar with the basics. Is there a way to get prediction intervals for gls models? My best model uses varPower(), so I am reluctant to fall back on lm predictions. Thank you, Hank Dr. Hank Stevens, Associate Professor 338 Pearson Hall Botany Department Miami University Oxford, OH 45056 Office: (513) 529-4206 Lab: (513) 529-4262 FAX: (513) 529-4243 http://www.cas.muohio.edu/~stevenmh/ http://www.muohio.edu/ecology/ http://www.muohio.edu/botany/ E Pluribus Unum __ 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.
Re: [R] The 'REP' term in AMMI{agricolae}
The AMMI senso strictu part starts from the corrected data for genotype and environment. In most cases where AMMI is applied (maybe also in the agricolae package, I haven't checked), starts from the interaction effects obtained through a general linear model anova. It should be possible to replace the input by blups from your mixed model. R news 7/1 (p14-19) gave some basic code for AMMI. There you should be able to get ideas how to modify the code. Joris CG Pettersson [EMAIL PROTECTED] e.slu.se To Sent by: r-help@stat.math.ethz.ch [EMAIL PROTECTED] cc at.math.ethz.ch Subject [R] The 'REP' term in 23/07/2007 19:59 AMMI{agricolae} Dear all, W2k, R 2.5.1 I am trying out the AMMI function in the agricolae package, to analyse the dependence of environment for a certain cultivar. The function responds to four basic variables: ENV Environment GEN Genotype REP Replication Y Response My question concerns the 'REP' term. When I normally do an analysis of variance, the replication would refer to repeated observations within the each field trial. In the case I am analysing right now, I have five years of observations for each 'ENV', in such a case it feels natural to use year as the most important replication of the environment- especially as I am interested in long term trends. This approach also seems to work allright. But the field trials are also structured in three main blocks, subdivided into five 'lattice' blocks, a structure that is powerful in the analysis of variance. (I use a random call in lme{nlme}). Is it possible to use the block structure also in the AMMI analysis? If that is possible, how should I code? I have tried to find out in the documentation, but if it is stated there I do not understand it. Thank you /CG -- CG Pettersson, PhD Swedish University of Agricultural Sciences (SLU) Dept. of Crop Production Ecology. Box 7043. SE-750 07 Uppsala, Sweden [EMAIL PROTECTED] __ 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.
Re: [R] extract from anova
Mustapha, You should tell us what fm1 actually is (which class), before we can give you a correct answer. It is helpfull to use str() in this case. str(fm1) will give you the structure of fm1, and with that information you extract any value you want. I suppose that this is what you need: summary(fm1)[[1]][Residuals,Mean Sq] elyakhlifi mustapha elyakhlifi_musta To [EMAIL PROTECTED] R-help@stat.math.ethz.ch Sent by: cc [EMAIL PROTECTED] at.math.ethz.ch Subject [R] extract from anova 16/07/2007 15:01 Hello, I would like extract the Mean Sq of the Residuals how should I do please? thanks. summary(fm1) Df Sum Sq Mean Sq F valuePr(F) groupe 20 300.987 15.049 22.853 3.369e-16 *** Residuals 41 27.000 0.659 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 _ [[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.
Re: [R] A More efficient method?
or Cat - c('a','a','a','b','b','b','a','a','b') C1 - (Cat=='a')*1 ONKELINX, Thierry Thierry.ONKELINX To @inbo.be Keith Alan Chamberlain Sent by: [EMAIL PROTECTED], [EMAIL PROTECTED] r-help@stat.math.ethz.ch at.math.ethz.chcc Subject 04/07/2007 17:17 Re: [R] A More efficient method? Cat - c('a','a','a','b','b','b','a','a','b') C1 - ifelse(Cat == 'a', -1, 1) ir. Thierry Onkelinx Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest Cel biometrie, methodologie en kwaliteitszorg / Section biometrics, methodology and quality assurance Gaverstraat 4 9500 Geraardsbergen Belgium tel. + 32 54/436 185 [EMAIL PROTECTED] www.inbo.be Do not put your faith in what statistics say until you have carefully considered what they do not say. ~William W. Watt A statistical analysis, properly conducted, is a delicate dissection of uncertainties, a surgery of suppositions. ~M.J.Moroney -Oorspronkelijk bericht- Van: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Namens Keith Alan Chamberlain Verzonden: woensdag 4 juli 2007 15:45 Aan: r-help@stat.math.ethz.ch Onderwerp: [R] A More efficient method? Dear Rhelpers, Is there a faster way than below to set a vector based on values from another vector? I'd like to call a pre-existing function for this, but one which can also handle an arbitrarily large number of categories. Any ideas? Cat=c('a','a','a','b','b','b','a','a','b') # Categorical variable C1=vector(length=length(Cat))# New vector for numeric values # Cycle through each column and set C1 to corresponding value of Cat. for(i in 1:length(C1)){ if(Cat[i]=='a') C1[i]=-1 else C1[i]=1 } C1 [1] -1 -1 -1 1 1 1 -1 -1 1 Cat [1] a a a b b b a a b Sincerely, KeithC. Psych Undergrad, CU Boulder (US) RE McNair Scholar __ 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. __ 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] unequal variance assumption for lme (mixed effect model)
gls() from the package nlme is similar to lm but is meant for models without random effects. Joris Spencer Graves [EMAIL PROTECTED] df.comTo Sent by: shirley zhang [EMAIL PROTECTED] [EMAIL PROTECTED] at.math.ethz.chcc R-help@stat.math.ethz.ch Subject 01/07/2007 23:30 Re: [R] unequal variance assumption for lme (mixed effect model) The 'weights' argument on 'lm' is assumed to identify a vector of the same length as the response, giving numbers that are inversely proportional to the variance for each observation. However, 'lm' provides no capability to estimate weights. If you want to do that, the varFunc capabilities in the 'nlme' package is the best tool I know for that purpose. If someone thinks there are better tools available for estimating heterscedasticity, I hope s/he will enlighten us both. Hope this helps. Spencer Graves shirley zhang wrote: Thanks for Spencer and Simon's help. I've got very interesting results based on your suggestions. One more question, how to handle unequal variance problme in lm()? Isn't the weights option also, which means weighted least squares, right? Can you give me an example of setting this parameter in lm() to account for different variance assumption in each group? Thanks again, Shirley On 6/29/07, Spencer Graves [EMAIL PROTECTED] wrote: comments in line shirley zhang wrote: Hi Simon, Thanks for your reply. Your reply reminds me that book. I've read it long time ago, but haven't try the weights option in my projects yet:) Is the heteroscedastic test always less powerful because we have to estimate the within group variance from the given data? SG: In general, I suspect we generally lose power when we estimate more parameters. SG: You can check this using the 'simulate.lme' function, whose use is illustrated in the seminal work reported in sect. 2.4 of Pinheiro and Bates (2000) Mixed-Effects Models in S and S-Plus (Springer). Should we check whether each group has equal variance before using weights=varIdent()? If we should, what is the function for linear mixed model? SG: The general advice I've seen is to avoid excessive overparameterization of heterscedasticity and correlations. However, parsimonious correlation had heterscedasticity models would likely be wise. Years ago, George Box expressed concern about people worrying too much about outliers, which are often fairly obvious and relatively easy to detect, while they worried too little, he thought, about dependence, especially serial dependence, which is generally more difficult to detect and creates bigger problems in inference than outliers. He wrote, Why worry about mice when there are tigers about? SG: Issues of this type can be fairly easily evaluated using 'simulate.lme'. Hope this helps. Spencer Graves Thanks, Shirley On 6/27/07, Simon Blomberg [EMAIL PROTECTED] wrote: The default settings for lme do assume equal variances within groups. You can change that by using the various varClasses. see ?varClasses. A simple example would be to allow unequal variances across groups. So if your call to lme was: lme(...,random=~1|group,...) then to allow each group to have its own variance, use: lme(...,random=~1|group, weights=varIdent(form=~1|group),...) You really really should read Pinheiro Bates (2000). It's all there. HTH, Simon. , On Wed, 2007-06-27 at 21:55 -0400, shirley zhang wrote: Dear Douglas and R-help, Does lme assume normal distribution AND equal variance among groups like anova() does? If it does, is there any method like unequal variance T-test (Welch T) in lme when each group has unequal variance in my data? Thanks, Shirley __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE
Re: [R] Linear Mixed Models with nlme, more than one random effect
Guilio, Have a look at Rnew volume 5/1 (http://cran.r-project.org/doc/Rnews/) where Doug Bates explains this nicely. Condider using lme4 for your purpose. Joris Giulio Di Giovanni [EMAIL PROTECTED] To otmail.com r-help@stat.math.ethz.ch Sent by: cc [EMAIL PROTECTED] at.math.ethz.ch Subject [R] Linear Mixed Models with nlme, more than one random effect 20/06/2007 14:09 Hi, I' trying to learn how to use lme for Linear Mixed Models and I have a problem when I have to include more than one random effect in my model. I know that this could be a stupid question to ask, but I'm not able to solve it by myself... One example: if my model is response = operator + block + day with operator and block as fixed effects and day as random effect, I use res.lme - lme(resp ~ oper + block , random=~1|day) If I want to include also another random effect, as experiment, what I should do ? This effect doesn't have to be nested, at the and I would like to have the COV matrix using (if I'm not wrong) getVarCov function. Thanks in advance for any help or suggestions, I'm a beginner on this field... Davide _ Cinema, Tv, Gossip e Orsoscopo…Tutto su MSN intrattenimento! __ 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.
Re: [R] Inverse BoxCox transformation
to backtransform 'estimate': if (lambda == 0 ) { log(estimate) } else { estimate^(1/lambda) } Des Callaghan [EMAIL PROTECTED] er.co.uk To Sent by: r-help@stat.math.ethz.ch [EMAIL PROTECTED] cc at.math.ethz.ch Subject [R] Inverse BoxCox transformation 18/06/2007 09:32 Hi, I can't seem to find a function in R that will reverse a BoxCox transformation. Can somebody help me locate one please? Thanks in advance. Best wishes, Des [[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.
Re: [R] str() to extract components
mod - lmer(resp3~b$age+b$size+b$pcfat+(1|sex), data=b) coef(mod)[1]$Subject[1,1] [EMAIL PROTECTED] wrote on 27/02/2007 16:52:54: Hi, I have been dabbling with str() to extract values from outputs such as lmer etc and have found it very helpful sometimes. but only seem to manage to extract the values when the output is one simple table, any more complicated and I'm stumped :-( take this example of the extracted coeficients from a lmer analysis... using str(coef(lmer(resp3~b$age+b$size+b$pcfat+(1|sex), data=b))) yields Formal class 'lmer.coef' [package Matrix] with 3 slots ..@ .Data :List of 1 .. ..$ :`data.frame': 2 obs. of 4 variables: .. .. ..$ (Intercept): num [1:2] 1.07 1.13 .. .. ..$ b$age : num [1:2] 0.00702 0.00702 .. .. ..$ b$size : num [1:2] 0.0343 0.0343 .. .. ..$ b$pcfat: num [1:2] 0.0451 0.0451 ..@ varFac: list() ..@ stdErr: num(0) how do I get inside the first table to get the value 1.07 for instance? Any help much appreciated. 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. __ 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] Hierarchical ANOVA
Romain, Look for info on mixed models. In R you do this either with the library nlme or lme4. A good starting point is an article by Doug Bates in Rnews http://cran.r-project.org/doc/Rnews/Rnews_2005-1.pdf Joris [EMAIL PROTECTED] wrote on 13/02/2007 17:16:41: Hello, Does somebody could help me in the computation(formulation) of a hierarchical ANOVA using linear model in R? I'm working in a population biology study of an endangered species. My aim is to see if I have effects of Density of individuals/m2 on several measured plants fitness traits. As independent variables I have: Humidity (continuous) Nutritive substance (continuous) LUX (continuous) Density of individuals/m2 (continuous) Population (categorical) Year (categorical) I want to test the 4 continous variables with Population as error term. And the Population, the Year and interaction term, with the error of Population x Year. Thank you for any help, best regards. Romain Mayor __ 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.
Re: [R] changes in the structure of mer objects?
Thanks Martin, R --vanilla did the trick. By the way, is there a way to check which version of a lme4 or Matrix I am using in a certain instance of R? We are running multiple versions of R on the same server until we are sure all our operational code is behaving well under a new version of R or of a package. Now I start doubting if we ever have been testing what we were intendng to test... Joris [EMAIL PROTECTED] wrote on 15/01/2007 14:00:12: Hi Joris, I suspect you somehow load an older version lme4 or Matrix than you think you are loading. Or then you have an lmer() function or a class definition {from a saved workspace } in your work space. example(lmer) *must* run correctly for the 'lme4' package to get onto CRAN at all, hence it must be unique to your setup. Maybe as a first step, run R as R --vanilla ? Regards, Martin Maechler, ETH Zurich joris == joris dewolf [EMAIL PROTECTED] on Sun, 14 Jan 2007 19:15:43 +0100 writes: joris Dear all, joris I try to run the example of lmer and get the following error message. library(lme4) example(lmer) lmer (fm1 - lmer(Reaction ~ Days + (Days | Subject), sleepstudy)) joris [[1]] joris Error in get(x, envir, mode, inherits) : variable as. dpoMatrix was not joris found joris This error message is similar to what I get with joris other models. It looks like the mer class has a joris slightly different structure. Anybody an idea how to joris solve this? joris I am using R 2.4.1 under linux and the latest releases of lme4 and Matrix joris lme4_0.9975-10 joris Matrix_0.9975-8 version joris _ joris platform x86_64-unknown-linux-gnu joris arch x86_64 joris os linux-gnu joris system x86_64, linux-gnu joris status joris major 2 joris minor 4.1 joris year 2006 joris month 12 joris day18 joris svn rev40228 joris language R joris version.string R version 2.4.1 (2006-12-18) __ 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.
[R] changes in the structure of mer objects?
Dear all, I try to run the example of lmer and get the following error message. library(lme4) example(lmer) lmer (fm1 - lmer(Reaction ~ Days + (Days | Subject), sleepstudy)) [[1]] Error in get(x, envir, mode, inherits) : variable as.dpoMatrix was not found This error message is similar to what I get with other models. It looks like the mer class has a slightly different structure. Anybody an idea how to solve this? I am using R 2.4.1 under linux and the latest releases of lme4 and Matrix lme4_0.9975-10 Matrix_0.9975-8 version _ platform x86_64-unknown-linux-gnu arch x86_64 os linux-gnu system x86_64, linux-gnu status major 2 minor 4.1 year 2006 month 12 day18 svn rev40228 language R version.string R version 2.4.1 (2006-12-18) Thanks Joris De Wolf Phone: +32 9 2429155, E-Mail: [EMAIL PROTECTED] Postal Address: CropDesign N.V. Technologiepark 3, 9052 Gent Belgium __ 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 extract the variance componets from lme
I advice you strongly to use VarCorr(), but if you insist tmp - as.matrix(m1$modelStruct$reStruct$MI) c(sqrt(diag(tmp)), Residual = 1) * m1$sigma Joris [EMAIL PROTECTED] wrote on 02/01/2007 12:50:05: Here is a piece of code fitting a model to a (part) of a dataset, just for illustration. I can extract the random interaction and the residual variance in group meth==1 using VarCorr, but how do I get the other residual variance? Is there any way to get the other variances in numerical form directly - it seems a litte contraintuitive to use as.numeric when extracting estimates, it's a bit like good old days writing SAS-programs that reads the SAS output files... Bendix --- library( nlme ) dfr - structure(list(meth = structure(c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2), .Label = c(CO, pulse), class = factor), item = c(1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4, 4, 1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4, 4), repl = c(1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3), y = c(78, 76.4, 77.2, 68.7, 67.6, 68.3, 82.9, 80.1, 80.7, 62.3, 65.8, 67.5, 71, 72, 73, 68, 67, 68, 82, 77, 77, 43, 69, 77)), .Names = c(meth, item, repl, y), row.names = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 184, 185, 186, 187, 188, 189, 190, 191, 192, 193, 194, 195), class = data.frame) m1 - lme( y ~ factor( meth ) + factor( item ), random = ~1 | MI, weights = varIdent( form = ~1 | meth ), method =REML, data = cbind( dfr, MI=interaction( dfr$meth, dfr$item ) ) ) m1 # The MI std and the residual std for meth==1 as.numeric(VarCorr(m1)[,2]) __ Bendix Carstensen Senior Statistician Steno Diabetes Center Niels Steensens Vej 2-4 DK-2820 Gentofte Denmark +45 44 43 87 38 (direct) +45 30 75 87 38 (mobile) +45 44 43 73 13 (fax) [EMAIL PROTECTED] http://www.biostat.ku.dk/~bxc __ 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.
[R] contrasts involving random terms in lme
Hello, Has somebody an idea of how to fit contrast involving random terms and obtain their standard errors with lme? I am referring to Welham, S.; Cullis, B.; Gogel, B.; Gilmour, A. and Thompson, R. (2004). Prediction in linear mixed models. Australian and New Zealand Journal of Statistics 46, 325-347 who develop the ideas of McLean, R.A.; Sanders,W.L.; and Stroup,W.W. (1991). A unified approach to linear mixed models. American Statistician 45, 5464. Any hint would be welcome. Joris -- == Joris De Wolf CropDesign N.V. Plant Evaluation Group Technologiepark 3 B-9052 Zwijnaarde Belgium Tel. : +32 9 242 91 55 Fax : +32 9 241 91 73 == confidentiality notice: The information contained in this e-mail is confidential and...{{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
[R] Lattice: how to index in a custom panel function?
Hi, I have a lattice xyplot that contains panels according to FactorA, and curves for the 2 levels of Factor B within a panel. I try to add text in the panels of a lattice graph. I suppose I have to write a custom function (panel.txt). What I really would like is to adapt the text in the panel according to the levels of FactorA. In the manuals, I find examples for the strips using which.given and which.panel. But this does not work for the panels... Can Anybody give a hint? Thanks Joris I work with R 1.9.1 under Linux panel.txt = function(sometxt,x,y,...){ grid.text(sometxt,x,y) } xyplot(data = Pdata, P ~ DAS | FactorA, groups = FactorB, type =s, col = c(red,blue), panel = function(x,y,...){ panel.abline(h = 100, lty = 5, lwd =0.5, col = darkgrey) panel.txt(mytext, 0.2, 0.8) } -- == Joris De Wolf CropDesign N.V. Plant Evaluation Group Technologiepark 3 B-9052 Zwijnaarde Belgium Tel. : +32 9 242 91 55 Fax : +32 9 241 91 73 == confidentiality notice: The information contained in this e-mail is confidential and...{{dropped}} __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] difference in order() between Linux and Windows with mixtures of caps and normal letters
Could anybody explain this difference in the function order() between R under windows and R under Linux ? It seems that under Linux the order in of character dsitinguishes between caps and normal letters and sorts them starting with the capitals (see below). How can I avoid this (I mean: get the same results as I get under Windows)? Thanks Joris Under Windows R 1.6.2 mstring - c(b,c,a,t,d) mstring[order(mstring)] [1] a b c d t Mstring - c(b,c,a,T,D) Mstring[order(Mstring)] [1] a b c D T the same under Linux R.1.6.2 mstring - c(b,c,a,t,d) mstring[order(mstring)] [1] a b c d t Mstring - c(b,c,a,T,D) Mstring[order(Mstring)] [1] D T a b c -- == Joris De Wolf CropDesign N.V. Plant Evaluation Group Technologiepark 3 B-9052 Zwijnaarde Belgium Tel. : +32 9 242 91 55 Fax : +32 9 241 91 73 == confidentiality notice: The information contained in this e-mail is confidential and...{{dropped}} __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Problems with postscript output
This is not a PS problem. The barplot you made has its horizontl axis suppressed by default. See argument axis.lty in ?barplot Joris Frank Gerrit Zoellner wrote: Hi all! I have a little problem with saving plots to file. I use the command postscript() followed by the plotting command and a dev.off(). When I then look at the resulting image saved to disk, some of the axis labels are missing (see attached image). Is there a way to fix this. Yours, __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- == Joris De Wolf CropDesign N.V. Plant Evaluation Group Technologiepark 3 B-9052 Zwijnaarde Belgium Tel. : +32 9 242 91 55 Fax : +32 9 241 91 73 == confidentiality notice: The information contained in this e-mail is confidential and...{{dropped}} __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Adding data.frames together
check rbind Joris John Sweval wrote: I have a series of data frames that are identical structurally, i.e. - made with the same code, but I need to add them together so that they become one, longer, data frame, i.e. - each of the slot vectors are increased in length by the length of the added data frame vectors. So if I have df1 with a slot A so that length(df1$A) = 100 and I have df2 with a slot A so that length(df2$A)=200 then I need a method to create df3 its slot A is the df1$A plus df2$A such that length(df3$A) = 300. It does not appear that if you use data.frame to join two data frames it just adds the slots of both sources to the destination data frame and that is not what I want. In my finally solution, I need to do this with multiple data.frames that slot-wise are identical, but each slot length is different between data frames. Seems like there should be an easy solution, but I just have not stumbled across it in the documentation. Thanks, John C. Sweval Database Architect Illumigen Biosciences, Inc. Email: [EMAIL PROTECTED] Phone: 206-378-0400 Fax: 206-378-0408 [[alternative HTML version deleted]] __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- == Joris De Wolf CropDesign N.V. Plant Evaluation Group Technologiepark 3 B-9052 Zwijnaarde Belgium Tel. : +32 9 242 91 55 Fax : +32 9 241 91 73 == confidentiality notice: The information contained in this e-mail is confidential and...{{dropped}} __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html