[R] Unexpected result of read.dbf
Hi there, Unexpected result is given with read.dbf from foreign 0.8-9, which is reproducible in R 2.1.1 with the following sample data: test.dbf KEYCODE,N,19,0 422010010 42201002101 42201002102 42201002103 42201002104 422010060 422010071 422010072 42201008001 42201008002 422010100 42201011001 42201011002 422010130 422010140 42201015001 42201015002 42201015003 422010180 422010190 In the example, the name of column is KEYCODE with numeric type, 19 digits, and no decimal. The data was saved as CSV format for confirmation. dbf.test - read.dbf(test.dbf) dbf.test[1:20,] [1] 422010010NANANANA 422010060 422010071 [8] 422010072NANA 422010100NANA 422010130 [15] 422010140NANANA 422010180 422010190 csv.test - read.csv(test.csv) csv.test[1:20,] [1] 422010010 42201002101 42201002102 42201002103 42201002104 422010060 [7] 422010071 422010072 42201008001 42201008002 422010100 42201011001 [13] 42201011002 422010130 422010140 42201015001 42201015002 42201015003 [19] 422010180 422010190 I wonder why I get NA from test.dbf. I have this kind of trouble when handling ESRI Shape files with maptools. There is no choice to use read.csv instead of read.dbf at executing read.shape. -- Susumu Tanimura __ 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] Problem with building R packages under Windows
Dear all, I am coming to the guru for advise here. I am a native user of Windows, and S-plus, and the process of migrating my libraries from S to R has been more than a painstaking task... I am currently using R version 2.1.1 in Windows XP SP2. I have read the Writing R extensions, the FAQ in R, and your valuable document R for Windows Users, but still unable to compile a package in R using Rcmd (likely stupidity on my part). The followings are what I have done: 1. Downloaded and setup ActivePerl and the path (in control panel) of c:\Perl\bin has been added automatically. I have downloaded version (5.8.7). 2. Downloaded and setup tools bundle from www.stats.ox.ac.uk/pub/Rtools/tools.zip and unzip them into c:\bin and added the path (in control panel) of c:\bin 3. Downloaded and setup the entire package of MikTex at http://www.miktex.org/ and added the path (in control panel) of C:\texmf\miktex\bin 4. Added the path (in control panel) referencing R (this version being R2.0patched) which is C:\Program Files\R\rw2000pat\bin 5. Since my codes are all in R scripts without compiled codes, I did not install MinGW. 6. At DOS prompt at the preamble of the directory called foo (which contained all subdirectories as created in R using package.skeleton(name = foo, list=cdgam.func)), I typed the following: C:\Rcmd build --binary foo Can't locate R/Dcf.pm in @INC @INC contains: c:\share\perl c:/Perl/lib c:/Perl/site/lib . at C:/bin/build line 29. BEGIN failed--compilation aborted at C:/bin/build line 29. I'd be very grateful if anyone could kindly advise me on this matter. Thanks. Lin __ 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
Re: [R] trouble with wilcox.test
On Wed, 17 Aug 2005, Greg Hather wrote: I'm having trouble with the wilcox.test command in R. Are you sure it is not the concepts that are giving 'trouble'? What real problem are you trying to solve here? To demonstrate the anomalous behavior of wilcox.test, consider wilcox.test(c(1.5,5.5), c(1:1), exact = F)$p.value [1] 0.01438390 wilcox.test(c(1.5,5.5), c(1:1), exact = T)$p.value [1] 6.39808e-07 (this calculation takes noticeably longer). wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value (R closes/crashes) I believe that wilcox.test(c(1.5,5.5), c(1:1), exact = F)$p.value yields a bad result because of the normal approximation which R uses when exact = F. Expecting an approximation to be good in the tail for m=2 is pretty unrealistic. But then so is believing the null hypothesis of a common *continuous* distribution. Why worry about the distribution under a hypothesis that is patently false? People often refer to this class of tests as `distribution-free', but they are not. The Wilcoxon test is designed for power against shift alternatives, but here there appears to be a very large difference in spread. So wilcox.test(5000+c(1.5,5.5), c(1:1), exact = T)$p.value [1] 0.9989005 even though the two samples differ in important ways. Any suggestions for how to compute wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value? I get (current R 2.1.1 on Linux) wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value [1] 1.59976e-07 and no crash. So the suggestion is to use a machine adequate to the task, and that probably means an OS with adequate stack size. [[alternative HTML version deleted]] PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Please do heed it. What version of R and what machine is this? And do take note of the request about HTML mail. -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ 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
Re: [R] Problem with building R packages under Windows
Dr L. Y Hin wrote: Dear all, I am coming to the guru for advise here. I am a native user of Windows, and S-plus, and the process of migrating my libraries from S to R has been more than a painstaking task... I am currently using R version 2.1.1 in Windows XP SP2. I have read the Writing R extensions, the FAQ in R, and your valuable document R for Windows Users, but still unable to compile a package in R using Rcmd (likely stupidity on my part). The followings are what I have done: 1. Downloaded and setup ActivePerl and the path (in control panel) of c:\Perl\bin has been added automatically. I have downloaded version (5.8.7). You need the native port, not the cygwin one. 2. Downloaded and setup tools bundle from www.stats.ox.ac.uk/pub/Rtools/tools.zip and unzip them into c:\bin and added the path (in control panel) of c:\bin Hmm, the tools are no longer available from this location! The current link is: http://www.murdoch-sutherland.com/Rtools/tools.zip For the current R version, the procedure is documented in the manual R Installation and Administration. 3. Downloaded and setup the entire package of MikTex at http://www.miktex.org/ and added the path (in control panel) of C:\texmf\miktex\bin 4. Added the path (in control panel) referencing R (this version being R2.0patched) which is C:\Program Files\R\rw2000pat\bin Above you told us about R-2.1.1 .! Please, really use a recent version. Uwe Ligges 5. Since my codes are all in R scripts without compiled codes, I did not install MinGW. 6. At DOS prompt at the preamble of the directory called foo (which contained all subdirectories as created in R using package.skeleton(name = foo, list=cdgam.func)), I typed the following: C:\Rcmd build --binary foo Can't locate R/Dcf.pm in @INC @INC contains: c:\share\perl c:/Perl/lib c:/Perl/site/lib . at C:/bin/build line 29. BEGIN failed--compilation aborted at C:/bin/build line 29. I'd be very grateful if anyone could kindly advise me on this matter. Thanks. Lin __ 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-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
Re: [R] trouble with wilcox.test
Prof Brian Ripley wrote: On Wed, 17 Aug 2005, Greg Hather wrote: I'm having trouble with the wilcox.test command in R. Are you sure it is not the concepts that are giving 'trouble'? What real problem are you trying to solve here? To demonstrate the anomalous behavior of wilcox.test, consider wilcox.test(c(1.5,5.5), c(1:1), exact = F)$p.value [1] 0.01438390 wilcox.test(c(1.5,5.5), c(1:1), exact = T)$p.value [1] 6.39808e-07 (this calculation takes noticeably longer). wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value (R closes/crashes) I believe that wilcox.test(c(1.5,5.5), c(1:1), exact = F)$p.value yields a bad result because of the normal approximation which R uses when exact = F. Expecting an approximation to be good in the tail for m=2 is pretty unrealistic. But then so is believing the null hypothesis of a common *continuous* distribution. Why worry about the distribution under a hypothesis that is patently false? People often refer to this class of tests as `distribution-free', but they are not. The Wilcoxon test is designed for power against shift alternatives, but here there appears to be a very large difference in spread. So wilcox.test(5000+c(1.5,5.5), c(1:1), exact = T)$p.value [1] 0.9989005 even though the two samples differ in important ways. Any suggestions for how to compute wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value? I get (current R 2.1.1 on Linux) wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value [1] 1.59976e-07 and no crash. So the suggestion is to use a machine adequate to the task, and that probably means an OS with adequate stack size. [[alternative HTML version deleted]] PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Please do heed it. What version of R and what machine is this? And do take note of the request about HTML mail. One could also try wilcox.exact() in package exactRankTests (0.8-11) which also gives (with suitable patience) [1] 1.59976e-07 even on my puny 256M Windows laptop. Still, it might be worthwhile adding a don't do something this silly error message to wilcox.test() rather than having it crash R. Low priority, IMHO. Windows XP SP2 R version 2.1.1, 2005-08-11 Peter Ehlers __ 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] Binary kernel discrimination
Hello, Could you tell me if a package exists to perform a binary kernel discrimination using a training set compose of molecules represented by binary fingerprint. This method was first described by Harper in J. Chem. Inf. Comput. Sci 2001 41 1295 and is also described in recent papers published in the same journal by Hert Jerome. I have attached the page describing the BKD method used in the last paper to give all the details of what I would like to try with R . Thanks for your help Fred BKD method.pdf __ 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
Re: [R] line/bar through median in lattice's bwplot?
M. K. mk_lists at yahoo.ca writes: Is there a way to render a line through the median point in the boxplot generated by the Lattice command bwplot? The line basically bisects the bar at the median point... bwplot(height~voice.part , pch='|', data=singer, xlab=Height (inches)) How to find this (haven't checked, maybe it's documented) library(lattice) panel.bwplot # no ()! Check the code for something like points (which is the default). Find a mysterious '|' and pch Dieter __ 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
Re: [R] accesing slots of S4 class in C code
I am trying to use a custom S4 object in my C code and I cannot get the
access to its slots working.
The following is a toy example, but even this crashes.
In R I have:
setClass(pd, representation(data=numeric))
x - new(pd, data=1:5)
test - function(pd.obj) {
res - .C(TestFun, pd.obj)
.Call(TestFun, pd.obj)
should work.
Torsten
res}
test(x)
(Of couse I load the DLL as well.)
The corresponding C file has:
SEXP TestFun(SEXP pd_obj)
{
SEXP t=R_NilValue;
PROTECT(t = GET_SLOT(pd_obj, install(data)));
UNPROTECT(1);
return(t);
}
What I hope to get as a result is the (1:5) vector.
In the long term, the vector will be a multi-dimensional array and I
will want to do calculations using its contents in the C program.
Thanks for any help,
Aniko
Huntsman Cancer Institute wishes to promote open communication while
protecting confidential and/or privileged information. If you have received
this message in error, please inform the sender and delete all copies.
[[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
__
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] Error messages using LMER
Dear All, After playing with lmer for couple of days, I have to say that I am amazed! I've been using quite some multilevel/mixed modeling packages, lme4 is a strong candidate for the overall winner, especially for multilevel generzlized linear models. Now go back to my two-level poisson model with cross-classified model. I've been testing various different model specificatios for the past couple of days. Here are the models I tried: 1) Two level random intercept model with level-1 covariates only m1 - lmer(.D ~ offset(log(.Y)) + (1|provn) +(1|bcohort) + x1 + x2 , data, poisson, method=Laplace) 2) Two-level random intercept model with both level-1 and level-2 covariates, but no cross-level interactions: m2 - lmer(.D ~ offset(log(.Y)) + (1|provn) +(1|bcohort) + x1 + x2 + z1 + z2, data, poisson, method=Laplace) 3) Two-level random intercept with cross-level interaction m3 - lmer(.D ~ offset(log(.Y)) + (1|provn) +(1|bcohort) + x1 + x2 + z1 + z2 + x1:z1 + x2:z2, data, poisson, method=Laplace) Both model 1 and 2 run fine. For model 3, I got error message: -- Error in fn(par, ...) : Unable to invert singular factor of downdated X'X In addition: Warning messages: 1: optim or nlminb returned message ERROR: ABNORMAL_TERMINATION_IN_LNSRCH in: LMEopt(x = mer, value = cv) 2: Leading minor of size 1 of downdated X'X is indefinite -- What is going on here? Any workarounds? Thanks! Best, Shige __ 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
Re: [R] trouble with wilcox.test
P Ehlers wrote: Prof Brian Ripley wrote: On Wed, 17 Aug 2005, Greg Hather wrote: I'm having trouble with the wilcox.test command in R. Are you sure it is not the concepts that are giving 'trouble'? What real problem are you trying to solve here? To demonstrate the anomalous behavior of wilcox.test, consider wilcox.test(c(1.5,5.5), c(1:1), exact = F)$p.value [1] 0.01438390 wilcox.test(c(1.5,5.5), c(1:1), exact = T)$p.value [1] 6.39808e-07 (this calculation takes noticeably longer). wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value (R closes/crashes) I believe that wilcox.test(c(1.5,5.5), c(1:1), exact = F)$p.value yields a bad result because of the normal approximation which R uses when exact = F. Expecting an approximation to be good in the tail for m=2 is pretty unrealistic. But then so is believing the null hypothesis of a common *continuous* distribution. Why worry about the distribution under a hypothesis that is patently false? People often refer to this class of tests as `distribution-free', but they are not. The Wilcoxon test is designed for power against shift alternatives, but here there appears to be a very large difference in spread. So wilcox.test(5000+c(1.5,5.5), c(1:1), exact = T)$p.value [1] 0.9989005 even though the two samples differ in important ways. Any suggestions for how to compute wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value? I get (current R 2.1.1 on Linux) wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value [1] 1.59976e-07 and no crash. So the suggestion is to use a machine adequate to the task, and that probably means an OS with adequate stack size. [[alternative HTML version deleted]] PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Please do heed it. What version of R and what machine is this? And do take note of the request about HTML mail. One could also try wilcox.exact() in package exactRankTests (0.8-11) which also gives (with suitable patience) [1] 1.59976e-07 even on my puny 256M Windows laptop. Still, it might be worthwhile adding a don't do something this silly error message to wilcox.test() rather than having it crash R. Low priority, IMHO. Windows XP SP2 R version 2.1.1, 2005-08-11 Peter Ehlers I should also mention package coin's wilcox_test() which does the job in about a quarter of the time used by exactRankTests. Peter Ehlers __ 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
Re: [R] newbie- how to add a row to a matrix
thank you all rbind solved my problems!!! simone Il giorno 18/ago/05, alle ore 00:34, Martin Lam ha scritto: # row bind a - matrix(1:5) a a - rbind(a, 6) a # column bind b - matrix(1:5) b b - cbind(b, 6:12) b b - cbind(b, 13) b Hope this helps, Martin --- Simone Gabbriellini [EMAIL PROTECTED] wrote: dear list, if I have a matrix s-matrix(1:5, ncol=5) how can I add another row with other data to that matrix? thank you, simone __ 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 __ Yahoo! Mail Stay connected, organized, and protected. Take the tour: http://tour.mail.yahoo.com/mailtour.html __ 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] matrix indexing
Dear R-users, I was wondering for the following: Let 'x' be a matrix and 'ind' and indicator matrix, i.e., x - array(1:20, dim = c(4, 5)) ind - array(c(1:3, 3:1), dim = c(3, 2)) I'd like to get (as a vector) the elements of 'x' which are not indexed by 'ind'. Since negative indices are not allowed in index matrices I thought of using something like: x[ind] - NA x[!is.na(x)] but are there any more elegant solutions. Thanks in advance, toka __ 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] GLMM - Am I trying the impossible?
Dear all, I have tried to calculate a GLMM fit with lmer (lme4) and glmmPQL (MASS), I also used glm for comparison. I am getting very different results from different functions, and I suspect that the problem is with our dataset rather than the functions, but I would appreciate help in deciding whether my suspicions are right. If indeed we are attempting the wrong type of analysis, some guidance about what would be the right thing to do would be greatly appreciated. The details: The data: The data are from the end-point of a survival experiment with fish. The design of the experiment is a 2 x 2 factorial, with each factor (Bacteria and Parasite) at two levels (yes and no). There were 16 fish in each tank, and the treatment was applied to the whole tank. There were in all 10 tanks (160 fish), with 2 tanks for controls (no/no), 2 tanks for (Parasite:yes/Bacteria:no) and 3 tanks for each of the other 2 treatments. A dead fish was considered a success, and a binomial family with the default logit link was used in the fits. No fish died in the control treatment (Is this the problem?). Overall probabilities as dead/total for the four treatments were: Paras Bact Prob nono0 yes no0.0625 noyes 0.208 yes yes 0.458 We are interested in testing main effects and interaction, but the interaction is of special interest to us. The data for dead are coded as 0/1 with 1 indicating a dead fish, and the file has one row per fish. Some results: lme4 (ver 0.98-1, R 2.1.1, Windows XP) fish1.lmerPQL - lmer(dead ~ Parasite * Bacteria + (1|Tank), data=fish.data, family=binomial) Error in lmer(dead ~ Parasite * Bacteria + (1 | Tank), data = fish.data, : Unable to invert singular factor of downdated X'X In addition: Warning message: Leading minor of size 4 of downdated X'X is indefinite without the interaction: fish3.lmerPQL - lmer(dead ~ Parasite + Bacteria + (1|Tank), data=fish.data, family=binomial) anova(fish3.lmerPQL) Analysis of Variance Table Df Sum Sq Mean Sq Denom F value Pr(F) Parasite 1 0.012 0.012 157.000 0.0103 0.9192 Bacteria 1 0.058 0.058 157.000 0.0524 0.8192 summary(fish3.lmerPQL) Generalized linear mixed model fit using PQL Formula: dead ~ Parasite + Bacteria + (1 | Tank) Data: fish.data Family: binomial(logit link) AIC BIClogLik deviance 141.3818 156.7577 -65.69091 131.3818 Random effects: GroupsNameVarianceStd.Dev. Tank (Intercept) 5e-10 2.2361e-05 # of obs: 160, groups: Tank, 10 Estimated scale (compare to 1) 0.9318747 Fixed effects: Estimate Std. Error z value Pr(|z|) (Intercept) -4.243800.79924 -5.3098 1.098e-07 *** Parasiteyes 1.264070.44694 2.8283 0.0046801 ** Bacteriayes 2.850620.75970 3.7523 0.0001752 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Correlation of Fixed Effects: (Intr) Prstys Parasiteyes -0.429 Bacteriayes -0.898 0.093 Very different P-values from anova and summary. MASS: (ver 7.2-17, R 2.1.1, Windows XP) ~ summary(fish1.glmmpql) Linear mixed-effects model fit by maximum likelihood Data: fish.data AIC BIC logLik 1236.652 1255.103 -612.326 Random effects: Formula: ~1 | Tank (Intercept) Residual StdDev: 0.02001341 0.8944214 Variance function: Structure: fixed weights Formula: ~invwt Fixed effects: dead ~ Parasite * Bacteria Value Std.Error DF t-value p-value (Intercept) -18.56607 1044.451 150 -0.01777591 0.9858 Parasiteyes 15.85802 1044.451 6 0.01518311 0.9884 Bacteriayes 17.23107 1044.451 6 0.01649772 0.9874 Parasiteyes:Bacteriayes -14.69007 1044.452 6 -0.01406487 0.9892 Correlation: (Intr) Prstys Bctrys Parasiteyes -1 Bacteriayes -1 1 Parasiteyes:Bacteriayes 1 -1 -1 Standardized Within-Group Residuals: MinQ1 MedQ3 Max -1.028634e+00 -5.734674e-01 -2.886770e-01 -4.224474e-14 4.330155e+00 Number of Observations: 160 Number of Groups: 10 anova(fish1.glmmpql) numDF denDF F-value p-value (Intercept) 1 150 17.452150 .0001 Parasite 1 6 4.136142 0.0882 Bacteria 1 6 12.740212 0.0118 Parasite:Bacteria 1 6 0.000198 0.9892 anova(glmmPQL(dead~Bacteria*Parasite, random=~1|Tank, family=binomial, data=fish.data)) iteration 1 numDF denDF F-value p-value (Intercept) 1 150 17.452150 .0001 Bacteria 1 6 8.980833 0.0241 Parasite 1 6 7.895521 0.0308 Bacteria:Parasite 1 6 0.000198 0.9892 Now anova indicates significance, but summary gives huge P-values. I have looked in MASS, ISwR, Fox's and Crawley's book, for hints, but I have probably missed the right
[R] R: stepwise procedures
hi all i found step(), stepAIC() and some other functions in leaps and subselect. is there a package/function that does the traditional stepwise regression procedures using F in and F out values? i know that step does the selctions based on AIC. / allan__ 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
Re: [R] trouble with wilcox.test
If this is stack overflow (and I don't know that yet: when I tried this on Windows the traceback was clearly corrupt, referring to bratio), the issue is that it is impossible to catch such an error, and it is not even AFAIK portably possible to find the stack size limit (or even the current usage) to do some estimates. (The amount of RAM is not relevant.) On Unix-alikes the stack size limit can be controlled from the shell used to launch R so we don't have any a priori knowledge. The underlying code could be rewritten not to use recursion, but that seems not worth the effort involved. All I can see we can do it to put a warning in the help file. On Thu, 18 Aug 2005, P Ehlers wrote: Prof Brian Ripley wrote: On Wed, 17 Aug 2005, Greg Hather wrote: I'm having trouble with the wilcox.test command in R. Are you sure it is not the concepts that are giving 'trouble'? What real problem are you trying to solve here? To demonstrate the anomalous behavior of wilcox.test, consider wilcox.test(c(1.5,5.5), c(1:1), exact = F)$p.value [1] 0.01438390 wilcox.test(c(1.5,5.5), c(1:1), exact = T)$p.value [1] 6.39808e-07 (this calculation takes noticeably longer). wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value (R closes/crashes) I believe that wilcox.test(c(1.5,5.5), c(1:1), exact = F)$p.value yields a bad result because of the normal approximation which R uses when exact = F. Expecting an approximation to be good in the tail for m=2 is pretty unrealistic. But then so is believing the null hypothesis of a common *continuous* distribution. Why worry about the distribution under a hypothesis that is patently false? People often refer to this class of tests as `distribution-free', but they are not. The Wilcoxon test is designed for power against shift alternatives, but here there appears to be a very large difference in spread. So wilcox.test(5000+c(1.5,5.5), c(1:1), exact = T)$p.value [1] 0.9989005 even though the two samples differ in important ways. Any suggestions for how to compute wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value? I get (current R 2.1.1 on Linux) wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value [1] 1.59976e-07 and no crash. So the suggestion is to use a machine adequate to the task, and that probably means an OS with adequate stack size. [[alternative HTML version deleted]] PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Please do heed it. What version of R and what machine is this? And do take note of the request about HTML mail. One could also try wilcox.exact() in package exactRankTests (0.8-11) which also gives (with suitable patience) [1] 1.59976e-07 even on my puny 256M Windows laptop. Still, it might be worthwhile adding a don't do something this silly error message to wilcox.test() rather than having it crash R. Low priority, IMHO. Windows XP SP2 R version 2.1.1, 2005-08-11 Peter Ehlers -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ 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
Re: [R] RMySQL not loading on Mac OS X
Bill, thanks a lot for your answer. I did not know about the sig-Mac list. I will post there next time if I do not find a solution. Concerning your suggestion: The system default compiler is gcc 4.0, but RMySQL seems to be built using gcc-3.3 regardless if I switch to 3.3. or not. Would it be a solution to force RMySQL to use gcc 4.0 during built? (Might be a naive idea, I am quite new to this). And if yes, how could I do this? Best, Georg On 18 Aug 2005, at 04:21, Bill Northcott wrote: On 11/08/2005, at 8:00 PM, Georg Otto wrote: I have a problem loading RMySQL 0.5-5 on Mac OS 10.4.2 running R 2.1.1. I installed RMySQL using: export PKG_CPPFLAGS=-I/usr/local/mysql/include export PKG_LIBS=-L/usr/local/mysql/lib -lmysqlclient R CMD INSTALL /Users/gwo/Desktop/RMySQL_0.5-5.tar.gz The installation seemed to work ok, but when I load RMySQL in R I get an error message: library(RMySQL) Error in dyn.load(x, as.logical(local), as.logical(now)) : unable to load shared library '/Library/Frameworks/ R.framework/Resources/library/RMySQL/libs/RMySQL.so': dlopen(/Library/Frameworks/R.framework/Resources/library/RMySQL/ libs/RMySQL.so, 6): Symbol not found: _printf$LDBLStub Referenced from: /Library/Frameworks/R.framework/Resources/ library/ RMySQL/libs/RMySQL.so Expected in: flat namespace Error in library(RMySQL) : .First.lib failed for 'RMySQL' Any hint will be highly appreciated! I notice this question never got a reply. It would have been better asked on the sig-Mac list, but here are some pointers. The libs/RMySQL.so, 6): Symbol not found: _printf$LDBLStub is caused by not linking all the necessary system libraries. Either 1. there is an attempt to link objects and/or static libraries built with gcc-3.x/g77 with objects produced by gcc-4.x/gfortran. This will not work. or 2. there is an attempt to link using ld or libtool rather than the gcc compiler driver which will ensure that appropriate system libraries are used. FWIW I had no problem building it, but I was using an R package which I built from source. So I know the same compiler was used throughout. If you are using the R binary distribution, make sure you have run 'sudo gcc_select 3.3' to get the right default compiler. Bill Northcott __ 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
Re: [R] matrix indexing
Hi, you might not consider this more elegant, but how about this x[-apply(ind, 1, function(i) (i[1]-1)*nrow(x) + i[2])] Ben On 8/18/05, toka tokas [EMAIL PROTECTED] wrote: Dear R-users, I was wondering for the following: Let 'x' be a matrix and 'ind' and indicator matrix, i.e., x - array(1:20, dim = c(4, 5)) ind - array(c(1:3, 3:1), dim = c(3, 2)) I'd like to get (as a vector) the elements of 'x' which are not indexed by 'ind'. Since negative indices are not allowed in index matrices I thought of using something like: x[ind] - NA x[!is.na(x)] but are there any more elegant solutions. Thanks in advance, toka __ 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-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
Re: [R] How to assess significance of random effect in lme4
Yes, it is a different issue. ranef() extracts the empirical Bayes estimates,
which are the empirical posterior modes. The bVar slot holds the corresponding
posterior variances of these modes.
Technically, (according to D. Bates) the values in the bVar slot are the the
diagonal elements of
(Z'Z+\Omega)^{-1}.
The original post was asking how to test and compare a specific random effect,
not a general assessment of how much information is provided by the data via
LRT.
Shige asked how to test whether a specific EB estimate is different than some
other value.
LRT doesn't answer this question, but the values in the bVar slot do.
-Original Message-
From: Spencer Graves [mailto:[EMAIL PROTECTED]
Sent: Wed 8/17/2005 10:08 PM
To: Doran, Harold
Cc: Shige Song; r-help@stat.math.ethz.ch
Subject:Re: [R] How to assess significance of random effect in lme4
Is there some reason you are NOT using anova, as in Examples
section of ?lmer?
Permit me to summarize what I know about this, and I'll be pleased if
someone else who thinks they know different would kindly enlighten me
and others who might otherwise be misled if anything I say is
inconsistent with the best literature available at the moment:
1. Doug Bates in his PhD dissertation and later in his book with Don
Watts (1988) Nonlinear Regression Analysis and Its Applications (Wiley)
split approximation errors in nonlinear least squares into intrinsic
curvature and parameter effects curvature. He quantified these two
problems in the context of roughly three dozen published examples, if my
memory is correct, and found that in not quite all cases, the parameter
effects were at least an order of magnitude greater than the intrinsic
curvature.
2. In nonnormal situations, maximum likelihood is subject to more
approximation error -- intrinsic curvature -- than simple nonlinear
least squares. However, I would expect this comparison to still be
fairly accurate, even if the differences may not be quite as stark.
3. The traditional use of standard errors to judge statistical
significance is subject to both intrinsic and parameter effects errors,
while likelihood ratio procedures such as anova are subject only to the
intrinsic curvature (assuming there are no substantive problems with
nonconvergence). Consequently, to judge statistical significance of an
effect, anova is usually substantially better than the so-called Wald
procedure using approximate standard errors, and is almost never worse.
If anyone knows of a case where this is NOT true, I'd like to know.
4. With parameters at a boundary as with variance components, the
best procedure seems to double the p-value from a nested anova (unless
the reported p-value is already large). This is because the
2*log(likelihood ratio) in such cases is roughly a 50-50 mixture of 0
and chi-square(1) [if testing only 1 variance component parameter].
This is supported by a substantial amount of research, including
simulations discussed in a chapter in Pinheiro and Bates (2000)
Mixed-Effects Models in S and S-Plus (Springer). The may be more
accurate procedures available in the literature, but none so simple as
this as far as I know.
Comments?
spencer graves
p.s. It looks like [EMAIL PROTECTED] is a list containing vectors of length 29
and 6 in your example. I don't know what they are, but I don't see how
they can be standard errors in the usual sense.
Doran, Harold wrote:
These are the posterior variances of the random effects (I think more
properly termed empirical posteriors). Your model apparently includes
three levels of random variation (commu, bcohort, residual). The first
are the variances associated with your commu random effect and the
second are the variances associated with the bcohort random effect.
Accessing either one would require
[EMAIL PROTECTED] or [EMAIL PROTECTED]
Obviously, replace fm with the name of your fitted model.
-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Shige Song
Sent: Wednesday, August 17, 2005 7:50 AM
To: r-help@stat.math.ethz.ch
Subject: Re: [R] How to assess significance of random effect in lme4
Hi Harold,
Thanks for the reply. I looked at my outputs using str() as you
suggested, here is the part you mentioned:
..@ bVar :List of 2
.. ..$ commu : num [1, 1, 1:29] 5e-10 5e-10 5e-10 5e-10 5e-10 ...
.. ..$ bcohort: num [1, 1, 1:6] 1.05e-05 7.45e-06 6.53e-06 8.25e-06
7.11e-06 ...
where commu and bcohort are the two second-level units. Are these
standard errors? Why the second vector contains a series of different
numbers?
Thanks!
Shige
On 8/17/05, Doran, Harold [EMAIL PROTECTED] wrote:
You can extract the posterior variance of the random effect from the
bVar slot of the fitted lmer model. It is not a hidden option, but a
part of
Re: [R] GLMM - Am I trying the impossible?
It is not supported to call anova() on a glmmPQL fit. For the glmmPQL fit you show, you have very large parameter estimates for a logistic and have partial separation (as you comment on for the control group): in that case PQL is not a reasonable method. Try fit - glm(dead ~ Parasite * Bacteria, data=fish.data, family=binomial) summary(fit) anova(fit, test=Chisq) fitted(fit) and you have fitted values of zero (up to numerical tolerances). This *is* discussed in MASS, around pp.198-9. So there is little point in adding random efects for that model. Now try fit2 - glmmPQL(dead ~ Parasite + Bacteria, random=~1|Tank, family=binomial, data=fish.data) summary(fit2) Fixed effects: dead ~ Bacteria + Parasite Value Std.Error DF t-value p-value (Intercept) -4.243838 0.7519194 150 -5.644007 0. Parasiteyes 1.264102 0.4205313 7 3.005964 0.0198 Bacteriayes 2.850640 0.7147180 7 3.988483 0.0053 which is pretty similar to the lmer fit you show. I don't know what anova is doing for your lmer fit, but I do know that it should not be working with sums of squares as is being reported. On Thu, 18 Aug 2005, Pedro J. Aphalo wrote: Dear all, I have tried to calculate a GLMM fit with lmer (lme4) and glmmPQL (MASS), I also used glm for comparison. I think you missed what glm was trying to tell you. I am getting very different results from different functions, and I suspect that the problem is with our dataset rather than the functions, but I would appreciate help in deciding whether my suspicions are right. If indeed we are attempting the wrong type of analysis, some guidance about what would be the right thing to do would be greatly appreciated. The details: The data: The data are from the end-point of a survival experiment with fish. The design of the experiment is a 2 x 2 factorial, with each factor (Bacteria and Parasite) at two levels (yes and no). There were 16 fish in each tank, and the treatment was applied to the whole tank. There were in all 10 tanks (160 fish), with 2 tanks for controls (no/no), 2 tanks for (Parasite:yes/Bacteria:no) and 3 tanks for each of the other 2 treatments. A dead fish was considered a success, and a binomial family with the default logit link was used in the fits. No fish died in the control treatment (Is this the problem?). Overall probabilities as dead/total for the four treatments were: Paras Bact Prob nono0 yes no0.0625 noyes 0.208 yes yes 0.458 We are interested in testing main effects and interaction, but the interaction is of special interest to us. The data for dead are coded as 0/1 with 1 indicating a dead fish, and the file has one row per fish. Some results: lme4 (ver 0.98-1, R 2.1.1, Windows XP) fish1.lmerPQL - lmer(dead ~ Parasite * Bacteria + (1|Tank), data=fish.data, family=binomial) Error in lmer(dead ~ Parasite * Bacteria + (1 | Tank), data = fish.data, : Unable to invert singular factor of downdated X'X In addition: Warning message: Leading minor of size 4 of downdated X'X is indefinite without the interaction: fish3.lmerPQL - lmer(dead ~ Parasite + Bacteria + (1|Tank), data=fish.data, family=binomial) anova(fish3.lmerPQL) Analysis of Variance Table Df Sum Sq Mean Sq Denom F value Pr(F) Parasite 1 0.012 0.012 157.000 0.0103 0.9192 Bacteria 1 0.058 0.058 157.000 0.0524 0.8192 summary(fish3.lmerPQL) Generalized linear mixed model fit using PQL Formula: dead ~ Parasite + Bacteria + (1 | Tank) Data: fish.data Family: binomial(logit link) AIC BIClogLik deviance 141.3818 156.7577 -65.69091 131.3818 Random effects: GroupsNameVarianceStd.Dev. Tank (Intercept) 5e-10 2.2361e-05 # of obs: 160, groups: Tank, 10 Estimated scale (compare to 1) 0.9318747 Fixed effects: Estimate Std. Error z value Pr(|z|) (Intercept) -4.243800.79924 -5.3098 1.098e-07 *** Parasiteyes 1.264070.44694 2.8283 0.0046801 ** Bacteriayes 2.850620.75970 3.7523 0.0001752 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Correlation of Fixed Effects: (Intr) Prstys Parasiteyes -0.429 Bacteriayes -0.898 0.093 Very different P-values from anova and summary. MASS: (ver 7.2-17, R 2.1.1, Windows XP) ~ summary(fish1.glmmpql) Linear mixed-effects model fit by maximum likelihood Data: fish.data AIC BIC logLik 1236.652 1255.103 -612.326 Random effects: Formula: ~1 | Tank (Intercept) Residual StdDev: 0.02001341 0.8944214 Variance function: Structure: fixed weights Formula: ~invwt Fixed effects: dead ~ Parasite * Bacteria Value Std.Error DF t-value p-value (Intercept) -18.56607 1044.451 150 -0.01777591 0.9858 Parasiteyes 15.85802 1044.451 6 0.01518311 0.9884 Bacteriayes 17.23107 1044.451 6 0.01649772 0.9874 Parasiteyes:Bacteriayes -14.69007
Re: [R] RMySQL not loading on Mac OS X
On 18/08/2005, at 7:46 PM, Georg Otto wrote: Concerning your suggestion: The system default compiler is gcc 4.0, but RMySQL seems to be built using gcc-3.3 regardless if I switch to 3.3. or not. Would it be a solution to force RMySQL to use gcc 4.0 during built? (Might be a naive idea, I am quite new to this). And if yes, how could I do this? When R is built it stores the configuration including the compilers which were used. When you try to build a source package, it uses the same compile commands. However, this is a less than perfect mechanism and requires that the packagers of binaries understand what is going on. For instance on Panther 'gcc' by default means gcc-3.3, whereas on Tiger by default it means gcc-4.0. So a package built with 'gcc' and 'g77' will work on Panther but not Tiger. See the problem. OTOH 'gcc-3.3' means gcc-3.3 on either Panther or Tiger, but 'gcc-4.0' won't work on Panther at all. The best solution is to build everything from source with the same compilers. That way you cannot have a compatibility problems. Hopefully soon the transition to gcc-4 will be complete. The change to Intel guarantees that, because the Apple gcc-3.3 compilers don't do x86 code. In a perfect future Apple will include Fortran in their Developer tools distribution, but for now they want gcc-4 and gfortran is not quite ready for the big time. Bill Northcott __ 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
Re: [R] power of a matrix
Hi, Gabor: Thanks. spencer graves Gabor Grothendieck wrote: Its expm. On 8/17/05, Spencer Graves [EMAIL PROTECTED] wrote: Hi, Peter: I couldn't find mexp in the Matrix package, but I did find it in fMultivar and in Lindsey's rmutil. These are different functions, but produced essentially the same answer for mexp(array(1:4, dim=c(2,2))). While hunting for that, I also also found reference by Doug Bates in a previous interchange on r-help to a classic paper ... I would recommend reading: Moler C., van Loan C., (2003); _Nineteen dubious ways to compute the exponential of a matrix, twenty-five years later_, SIAM Review 45, 3-49. This paper was cited in the help page for mexp in fMultivar but not in rmutil. spencer graves Peter Dalgaard wrote: Rau, Roland [EMAIL PROTECTED] writes: Thank you very much! Thanks also to the authors of this function, Vincente Canto Cassola and Martin Maechler! This is exactly what I hoped for. look at function ?mtx.exp() in the Malmig package, e.g., Also, there is mexp() in the Matrix package. I'm not sure about the relative merits. mexp() is one of the less dubious implementations of matrix exponentials, but it does require to and from class Matrix. mtx.exp is a bit unfortunately named as it appears to calculate matrix *powers* (which in this case is what you need). -- Spencer Graves, PhD Senior Development Engineer PDF Solutions, Inc. 333 West San Carlos Street Suite 700 San Jose, CA 95110, USA [EMAIL PROTECTED] www.pdf.com http://www.pdf.com Tel: 408-938-4420 Fax: 408-280-7915 __ 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 -- Spencer Graves, PhD Senior Development Engineer PDF Solutions, Inc. 333 West San Carlos Street Suite 700 San Jose, CA 95110, USA [EMAIL PROTECTED] www.pdf.com http://www.pdf.com Tel: 408-938-4420 Fax: 408-280-7915 __ 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
Re: [R] How to assess significance of random effect in lme4
Hi, Harold: Thanks for the clarification. I thought I had read the
original post. Obviously, I had misread it. Thanks again. spencer graves
Doran, Harold wrote:
Yes, it is a different issue. ranef() extracts the empirical Bayes
estimates, which are the empirical posterior modes. The bVar slot holds
the corresponding posterior variances of these modes.
Technically, (according to D. Bates) the values in the bVar slot are the
the diagonal elements of
(Z'Z+\Omega)^{-1}.
The original post was asking how to test and compare a specific random
effect, not a general assessment of how much information is provided by
the data via LRT.
Shige asked how to test whether a specific EB estimate is different than
some other value.
LRT doesn't answer this question, but the values in the bVar slot do.
-Original Message-
From: Spencer Graves [mailto:[EMAIL PROTECTED]
Sent: Wed 8/17/2005 10:08 PM
To: Doran, Harold
Cc: Shige Song; r-help@stat.math.ethz.ch
Subject:Re: [R] How to assess significance of random effect in lme4
Is there some reason you are NOT using anova, as in Examples
section of ?lmer?
Permit me to summarize what I know about this, and I'll be
pleased if
someone else who thinks they know different would kindly enlighten me
and others who might otherwise be misled if anything I say is
inconsistent with the best literature available at the moment:
1. Doug Bates in his PhD dissertation and later in his book
with Don
Watts (1988) Nonlinear Regression Analysis and Its Applications (Wiley)
split approximation errors in nonlinear least squares into intrinsic
curvature and parameter effects curvature. He quantified these two
problems in the context of roughly three dozen published examples, if my
memory is correct, and found that in not quite all cases, the parameter
effects were at least an order of magnitude greater than the intrinsic
curvature.
2. In nonnormal situations, maximum likelihood is subject to more
approximation error -- intrinsic curvature -- than simple nonlinear
least squares. However, I would expect this comparison to still be
fairly accurate, even if the differences may not be quite as stark.
3. The traditional use of standard errors to judge statistical
significance is subject to both intrinsic and parameter effects errors,
while likelihood ratio procedures such as anova are subject only to the
intrinsic curvature (assuming there are no substantive problems with
nonconvergence). Consequently, to judge statistical significance of an
effect, anova is usually substantially better than the so-called Wald
procedure using approximate standard errors, and is almost never worse.
If anyone knows of a case where this is NOT true, I'd like to know.
4. With parameters at a boundary as with variance components, the
best procedure seems to double the p-value from a nested anova (unless
the reported p-value is already large). This is because the
2*log(likelihood ratio) in such cases is roughly a 50-50 mixture of 0
and chi-square(1) [if testing only 1 variance component parameter].
This is supported by a substantial amount of research, including
simulations discussed in a chapter in Pinheiro and Bates (2000)
Mixed-Effects Models in S and S-Plus (Springer). The may be more
accurate procedures available in the literature, but none so simple as
this as far as I know.
Comments?
spencer graves
p.s. It looks like [EMAIL PROTECTED] is a list containing vectors of length
29
and 6 in your example. I don't know what they are, but I don't see how
they can be standard errors in the usual sense.
Doran, Harold wrote:
These are the posterior variances of the random effects (I think more
properly termed empirical posteriors). Your model apparently includes
three levels of random variation (commu, bcohort, residual). The first
are the variances associated with your commu random effect and the
second are the variances associated with the bcohort random effect.
Accessing either one would require
[EMAIL PROTECTED] or [EMAIL PROTECTED]
Obviously, replace fm with the name of your fitted model.
-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Shige Song
Sent: Wednesday, August 17, 2005 7:50 AM
To: r-help@stat.math.ethz.ch
Subject: Re: [R] How to assess significance of random effect in lme4
Hi Harold,
Thanks for the reply. I looked at my outputs using str() as you
suggested, here is the part you mentioned:
..@ bVar :List of 2
.. ..$ commu : num [1, 1, 1:29] 5e-10 5e-10 5e-10 5e-10 5e-10 ...
.. ..$ bcohort: num [1, 1, 1:6] 1.05e-05 7.45e-06 6.53e-06 8.25e-06
7.11e-06 ...
where commu and bcohort are the two second-level units. Are these
standard errors? Why the second vector
Re: [R] matrix indexing
x - array(1:20, dim = c(4, 5)) ind - array(c(1:3, 3:1), dim = c(3, 2)) # instead of using ind (pairs of coordinates) for getting the items in the matrix, you can convert it to a list of single coordinates to point to the item in the matrix: # t = transpose # nrow = get the number of rows indices - t((ind[,2]-1) * nrow(x) + ind[,1]) # remove the ones indicated by indices # to get an item from y you can do y[value] y - x[-indices] # you could transpose y, but I don't know if that's what you want y - t(y) Hope this helps, Martin --- toka tokas [EMAIL PROTECTED] wrote: Dear R-users, I was wondering for the following: Let 'x' be a matrix and 'ind' and indicator matrix, i.e., x - array(1:20, dim = c(4, 5)) ind - array(c(1:3, 3:1), dim = c(3, 2)) I'd like to get (as a vector) the elements of 'x' which are not indexed by 'ind'. Since negative indices are not allowed in index matrices I thought of using something like: x[ind] - NA x[!is.na(x)] but are there any more elegant solutions. Thanks in advance, toka __ 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
Re: [R] [SPAM] - Re: How to assess significance of random effect in lme4 - Bayesian Filter detected spam
Actually, I re-read the post and think it needs clarification. We may
both be right. If the question is I am building a model and want to
know if I should retain this random effect? (or something like that)
then the LRT should be used to compare the fitted model against another
model. This would be accomplished via anova().
In other multilevel programs, the variance components are often
associated with a chi-square statistic and a test statistic associated
with the variance. lmer() does not report this test statistic.
Now, if the question is something like I want to know if a specific
realization of the random variable (i.e., a specific empirical Bayes
estimate) is different from a population value? then one would need the
posterior means.
So, Shige, I hope this hasn't been confusing, but there are many things
happening in these models and it is easy to get confused. Maybe if you
could clarify your question.
-Original Message-
From: Spencer Graves [mailto:[EMAIL PROTECTED]
Sent: Thursday, August 18, 2005 8:26 AM
To: Doran, Harold
Cc: Shige Song; r-help@stat.math.ethz.ch
Subject: [SPAM] - Re: [R] How to assess significance of random effect in
lme4 - Bayesian Filter detected spam
Hi, Harold: Thanks for the clarification. I thought I had read the
original post. Obviously, I had misread it. Thanks again. spencer
graves
Doran, Harold wrote:
Yes, it is a different issue. ranef() extracts the empirical Bayes
estimates, which are the empirical posterior modes. The bVar slot
holds the corresponding posterior variances of these modes.
Technically, (according to D. Bates) the values in the bVar slot are
the the diagonal elements of (Z'Z+\Omega)^{-1}.
The original post was asking how to test and compare a specific random
effect, not a general assessment of how much information is provided
by the data via LRT.
Shige asked how to test whether a specific EB estimate is different
than some other value.
LRT doesn't answer this question, but the values in the bVar slot do.
-Original Message-
From: Spencer Graves [mailto:[EMAIL PROTECTED]
Sent: Wed 8/17/2005 10:08 PM
To: Doran, Harold
Cc: Shige Song; r-help@stat.math.ethz.ch
Subject:Re: [R] How to assess significance of random effect in
lme4
Is there some reason you are NOT using anova, as in
Examples
section of ?lmer?
Permit me to summarize what I know about this, and I'll be
pleased if someone else who thinks they know different would kindly
enlighten me and others who might otherwise be misled if anything I
say is inconsistent with the best literature available at the moment:
1. Doug Bates in his PhD dissertation and later in his book
with Don Watts (1988) Nonlinear Regression Analysis and Its
Applications (Wiley) split approximation errors in nonlinear least
squares into intrinsic curvature and parameter effects curvature.
He quantified these two problems in the context of roughly three dozen
published examples, if my memory is correct, and found that in not
quite all cases, the parameter effects were at least an order of
magnitude greater than the intrinsic curvature.
2. In nonnormal situations, maximum likelihood is subject
to more approximation error -- intrinsic curvature -- than simple
nonlinear least squares. However, I would expect this comparison to
still be fairly accurate, even if the differences may not be quite as
stark.
3. The traditional use of standard errors to judge
statistical significance is subject to both intrinsic and parameter
effects errors, while likelihood ratio procedures such as anova are
subject only to the intrinsic curvature (assuming there are no
substantive problems with nonconvergence). Consequently, to judge
statistical significance of an effect, anova is usually substantially
better than the so-called Wald procedure using approximate standard
errors, and is almost never worse.
If anyone knows of a case where this is NOT true, I'd like to know.
4. With parameters at a boundary as with variance
components, the best procedure seems to double the p-value from a
nested anova (unless the reported p-value is already large). This is
because the 2*log(likelihood ratio) in such cases is roughly a 50-50
mixture of 0 and chi-square(1) [if testing only 1 variance component
parameter].
This is supported by a substantial amount of research, including
simulations discussed in a chapter in Pinheiro and Bates (2000)
Mixed-Effects Models in S and S-Plus (Springer). The may be more
accurate procedures available in the literature, but none so simple as
this as far as I know.
Comments?
spencer graves
p.s. It looks like [EMAIL PROTECTED] is a list containing vectors of length
29
and 6 in your example. I don't know what they are, but I don't see
how they can be standard errors in the usual sense.
[R] lme model: Error in MEEM
Hi, We have data of two groups of subjects: 32 elderly, 14 young adults. for each subject we have 15 observations, each observation consisting of a reaction-time measure (RT) and an activation maesure (betadlpcv). since we want to analyze the influence of (age-)group and RT on the activation, we call: lme(betadlpcv ~ RT*group, data=our.data, random=~ RT |subject) this yields: Error in MEEM(object, conLin, control$niterEM) : Singularity in backsolve at level 0, block 1 In addition: Warning message: Fewer observations than random effects in all level 1 groups in: lme.formula(betadlpcv ~ RT * group, data = patrizia.data, random = ~RT | what's the problem here? thanks for your kind help christoph -- Lust, ein paar Euro nebenbei zu verdienen? Ohne Kosten, ohne Risiko! Satte Provisionen für GMX Partner: http://www.gmx.net/de/go/partner __ 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] select previous date
Hi everybody, could anyone help me in finding a way for selecting from a dataframe each row that comes before another that matches a condition? EXAMPLE: df raw.number name date temp 1aaa 2001-04-15 15 2aaa 2001-01-15 12 3aaa 2001-03-15 13 ... i-1 bbb 2002-04-15 15 ibbb 2002-03-15 14 the condition is: df$temp==15 matching raws are 1 and i-1: I need something to select (only) rows where date=one month before the matching raws, so raws 3 and i. (the variable name has more than one level) 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
Re: [R] lme model: Error in MEEM
On 8/18/05, Christoph Lehmann [EMAIL PROTECTED] wrote: Hi, We have data of two groups of subjects: 32 elderly, 14 young adults. for each subject we have 15 observations, each observation consisting of a reaction-time measure (RT) and an activation maesure (betadlpcv). since we want to analyze the influence of (age-)group and RT on the activation, we call: lme(betadlpcv ~ RT*group, data=our.data, random=~ RT |subject) this yields: Error in MEEM(object, conLin, control$niterEM) : Singularity in backsolve at level 0, block 1 In addition: Warning message: Fewer observations than random effects in all level 1 groups in: lme.formula(betadlpcv ~ RT * group, data = patrizia.data, random = ~RT | what's the problem here? It seems that you only have one observation per subject and you are trying to estimate a model with two random effects per subject plus the per-observation noise term. These terms are completely confounded. thanks for your kind help christoph -- Lust, ein paar Euro nebenbei zu verdienen? Ohne Kosten, ohne Risiko! Satte Provisionen für GMX Partner: http://www.gmx.net/de/go/partner __ 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-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] problems when installing R in Fedora core 4
Hi, I got a problem when installing R in Fedora core 4. When I ran .configure, it gave the following error message: configure: error: --with-x=yes (default) and X11 headers/libs are not available Could anyone tell me what's wrong? Am I missing some package in Fedora? Thanks a lot for your help. Peter __ 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] help with unknown function
Hello I am working on radio tracking data, with a short programme someone gave me and ... which should, supposedly, work ... In this programme, there is the function : getareahr(kern, levels = 95). But i cannot find any 'getareahr' in R ... could anyone help me? thanks! Agnès --- Agnès GAULT graduate student UMR 5173 MNHN-CNRS-UPMC case postale 50 Species Conservation, Restoration and Population Survey (CERSP) 61 rue Buffon, 1er étage 75005 PARIS FRANCE Tel: 33 (0)1 40 79 57 64 Fax: 33 (0)1 40 79 38 35 Email: [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
Re: [R] accesing slots of S4 class in C code
On 8/17/05, Aniko Szabo [EMAIL PROTECTED] wrote:
I am trying to use a custom S4 object in my C code and I cannot get the
access to its slots working.
The following is a toy example, but even this crashes.
In R I have:
setClass(pd, representation(data=numeric))
x - new(pd, data=1:5)
test - function(pd.obj) {
res - .C(TestFun, pd.obj)
res}
test(x)
(Of couse I load the DLL as well.)
The corresponding C file has:
SEXP TestFun(SEXP pd_obj)
{
SEXP t=R_NilValue;
PROTECT(t = GET_SLOT(pd_obj, install(data)));
UNPROTECT(1);
return(t);
}
Your C code is set for the .Call interface, not the .C interface.
Change the R code to
res - .Call(TestFun, pd.obj)
What I hope to get as a result is the (1:5) vector.
In the long term, the vector will be a multi-dimensional array and I
will want to do calculations using its contents in the C program.
Thanks for any help,
Aniko
Huntsman Cancer Institute wishes to promote open communication while
protecting confidential and/or privileged information. If you have received
this message in error, please inform the sender and delete all copies.
[[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
__
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] display of a loess fitted surface
Good morning, I am Marta Colombo,student at Politecnico,Milan. I am studying local regression models and I am using loess function. My problem is that when I have a loess object I don't know how to display the fitted surface; in fact, while in S when you have a loess object you can see it writing plot(object), in R this dosen't work. Also I'd like to know if there is something like the S function pointwise that computes upper and lower confidence intervals. Thank you very much for your attention. Marta Colombo __ 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] display of a loess fitted surface
Good morning, I am Marta Colombo,student at Politecnico,Milan. I am studying local regression models and I am using loess function. My problem is that when I have a loess object I don't know how to display the fitted surface; in fact, while in S when you have a loess object you can see it writing plot(object), in R this dosen't work. Also I'd like to know if there is something like the S function pointwise that computes upper and lower confidence intervals. Thank you very much for your attention. Marta Colombo __ 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
Re: [R] problems when installing R in Fedora core 4
You are probably missing the 'devel' package for x11 which includes header files. -roger Peter Yang wrote: Hi, I got a problem when installing R in Fedora core 4. When I ran .configure, it gave the following error message: configure: error: --with-x=yes (default) and X11 headers/libs are not available Could anyone tell me what's wrong? Am I missing some package in Fedora? Thanks a lot for your help. Peter __ 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 -- Roger D. Peng http://www.biostat.jhsph.edu/~rpeng/ __ 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] how to draw an ellipse
That's probably a stupid question, but I'm looking for a low-level command which plots ellipse, specifying only center and deformation axes. The purpose is to illustrate bivariates gaussians with 2D .95 confident levels without using any specific library. Thanxs for your help, Guillaume __ 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
Re: [R] accesing slots of S4 class in C code
Hi Aniko, On 17 Aug 2005, [EMAIL PROTECTED] wrote: I am trying to use a custom S4 object in my C code and I cannot get the access to its slots working. The following is a toy example, but even this crashes. res - .C(TestFun, pd.obj) I'm pretty sure you want to use the .Call interface for this sort of thing, not .C. Other than referring you to the Writing R Extensions manual, you might want to look at the Ruuid package in Bioconductor which is quite simple, but demonstrates accessing S4 classes from C. See http://bioconductor.org/ to grab a source package of Ruuid. Best, + seth PS: Questions of this nature (C code type stuff) are better suited to the R-devel mail list. __ 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
Re: [R] problems when installing R in Fedora core 4
On 08/18/05 09:31, Peter Yang wrote: Hi, I got a problem when installing R in Fedora core 4. When I ran .configure, it gave the following error message: configure: error: --with-x=yes (default) and X11 headers/libs are not available Could anyone tell me what's wrong? Am I missing some package in Fedora? Yes, probably xorg-x11-devel, so do yum install xorg-x11-devel and you will probably get rid of THAT error message. I'm treading into deeper water in what follows, but installation depends on several other packages as well. I discovered most of these by looking at the error messages and then finding what rpm provided the files needed. I did this either with rpm -q --whatprovides [the missing file] on some other installation that happened to work, or by searching http://rpm.pbone.net. Whether you have the needed files already depends on what kind of installation you did. Some of the packages are devel and others are compat. Here is my list of compat rpms that I have installed, and I think I installed all of these just to get R to build: compat-libf2c-32-3.2.3-47.fc4 compat-libstdc++-296-2.96-132.fc4 compat-readline43-4.3-2 compat-gcc-32-3.2.3-47.fc4 My hunch is that I still do not have the optimal installation, but it is possible that the newest versions of gcc have solved some of the problems with the ones that originally came with FC4. I've seen some discussion suggesting that the way to go is to use an older version of gcc, but I did not search for it just now. Jon -- Jonathan Baron, Professor of Psychology, University of Pennsylvania Home page: http://www.sas.upenn.edu/~baron __ 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
Re: [R] Error messages using LMER
On 8/18/05, Shige Song [EMAIL PROTECTED] wrote: Dear All, After playing with lmer for couple of days, I have to say that I am amazed! I've been using quite some multilevel/mixed modeling packages, lme4 is a strong candidate for the overall winner, especially for multilevel generzlized linear models. Now go back to my two-level poisson model with cross-classified model. I've been testing various different model specificatios for the past couple of days. Here are the models I tried: 1) Two level random intercept model with level-1 covariates only m1 - lmer(.D ~ offset(log(.Y)) + (1|provn) +(1|bcohort) + x1 + x2 , data, poisson, method=Laplace) 2) Two-level random intercept model with both level-1 and level-2 covariates, but no cross-level interactions: m2 - lmer(.D ~ offset(log(.Y)) + (1|provn) +(1|bcohort) + x1 + x2 + z1 + z2, data, poisson, method=Laplace) 3) Two-level random intercept with cross-level interaction m3 - lmer(.D ~ offset(log(.Y)) + (1|provn) +(1|bcohort) + x1 + x2 + z1 + z2 + x1:z1 + x2:z2, data, poisson, method=Laplace) Both model 1 and 2 run fine. For model 3, I got error message: -- Error in fn(par, ...) : Unable to invert singular factor of downdated X'X In addition: Warning messages: 1: optim or nlminb returned message ERROR: ABNORMAL_TERMINATION_IN_LNSRCH in: LMEopt(x = mer, value = cv) 2: Leading minor of size 1 of downdated X'X is indefinite -- What is going on here? Any workarounds? Thanks! The first thing I would try is set the EMverbose and msVerbose flags in the control list to see what occurs within the optimization. That is append the argument control = list(EMverbose = TRUE, msVerbose = TRUE) to your call to lmer(). You may also want to try the call in a recently compiled R-devel, which will be released as R-2.2.0 in October. You will notice that the first warning message reads optim or nlminb. In R-2.1.1 lmer uses optim for the optimization. Starting with R-2.2.0 the default is to use nlminb. Test compilations of R-devel for Windows are available from CRAN. __ 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
Re: [R] problems when installing R in Fedora core 4
On Thu, 2005-08-18 at 09:31 -0400, Peter Yang wrote: Hi, I got a problem when installing R in Fedora core 4. When I ran .configure, it gave the following error message: configure: error: --with-x=yes (default) and X11 headers/libs are not available Could anyone tell me what's wrong? Am I missing some package in Fedora? Thanks a lot for your help. Peter Yes, the development headers for X11. If you use YUM, do the following in a shell: su -c yum install xorg-x11-devel Enter your root password and it will download and install the relevant headers. As an aside, are you aware of Martyn Plummer's R rpm binary? You can get it from here: http://www.stats.bris.ac.uk/R/bin/linux/redhat/fc4/ HTH Gav __ 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 -- %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~% Gavin Simpson [T] +44 (0)20 7679 5522 ENSIS Research Fellow [F] +44 (0)20 7679 7565 ENSIS Ltd. ECRC [E] gavin.simpsonATNOSPAMucl.ac.uk UCL Department of Geography [W] http://www.ucl.ac.uk/~ucfagls/cv/ 26 Bedford Way[W] http://www.ucl.ac.uk/~ucfagls/ London. WC1H 0AP. %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~% __ 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
Re: [R] how to draw an ellipse
Guillaume Allain guillaume.allain at cbconseil.com writes: That's probably a stupid question, but I'm looking for a low-level command which plots ellipse, specifying only center and deformation axes. The purpose is to illustrate bivariates gaussians with 2D .95 confident levels without using any specific library. Thanxs for your help, Guillaume I don't know exactly why you want to avoid using any specific library, but the ellipse package (sic) would seem to do what you want pretty conveniently ... cheers Ben Bolker __ 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
Re: [R] lme model: Error in MEEM
sorry, RT had an error in raw data and was treated as a factor. after correction of the raw data (RT is numeric) now it works fine. thanks a lot christoph --- Ursprüngliche Nachricht --- Von: Douglas Bates [EMAIL PROTECTED] An: Christoph Lehmann [EMAIL PROTECTED] Kopie: r-help@stat.math.ethz.ch Betreff: Re: [R] lme model: Error in MEEM Datum: Thu, 18 Aug 2005 08:29:52 -0500 On 8/18/05, Christoph Lehmann [EMAIL PROTECTED] wrote: Hi, We have data of two groups of subjects: 32 elderly, 14 young adults. for each subject we have 15 observations, each observation consisting of a reaction-time measure (RT) and an activation maesure (betadlpcv). since we want to analyze the influence of (age-)group and RT on the activation, we call: lme(betadlpcv ~ RT*group, data=our.data, random=~ RT |subject) this yields: Error in MEEM(object, conLin, control$niterEM) : Singularity in backsolve at level 0, block 1 In addition: Warning message: Fewer observations than random effects in all level 1 groups in: lme.formula(betadlpcv ~ RT * group, data = patrizia.data, random = ~RT | what's the problem here? It seems that you only have one observation per subject and you are trying to estimate a model with two random effects per subject plus the per-observation noise term. These terms are completely confounded. thanks for your kind help christoph -- Lust, ein paar Euro nebenbei zu verdienen? Ohne Kosten, ohne Risiko! Satte Provisionen für GMX Partner: http://www.gmx.net/de/go/partner __ 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-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
Re: [R] trouble with wilcox.test
I think your guess about stack overflow is probably correct and I definitely don't think it's worth wasting effort recoding. Peter Ehlers Prof Brian Ripley wrote: If this is stack overflow (and I don't know that yet: when I tried this on Windows the traceback was clearly corrupt, referring to bratio), the issue is that it is impossible to catch such an error, and it is not even AFAIK portably possible to find the stack size limit (or even the current usage) to do some estimates. (The amount of RAM is not relevant.) On Unix-alikes the stack size limit can be controlled from the shell used to launch R so we don't have any a priori knowledge. The underlying code could be rewritten not to use recursion, but that seems not worth the effort involved. All I can see we can do it to put a warning in the help file. [snip] __ 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
Re: [R] help with unknown function
Agnes Gault a écrit : Hello I am working on radio tracking data, with a short programme someone gave me and ... which should, supposedly, work ... In this programme, there is the function : getareahr(kern, levels = 95). But i cannot find any 'getareahr' in R ... could anyone help me? It looks difficult given the (lack of) information you're giving. What are you trying to do ? Best, Renaud thanks! Agnès --- Agnès GAULT graduate student UMR 5173 MNHN-CNRS-UPMC case postale 50 Species Conservation, Restoration and Population Survey (CERSP) 61 rue Buffon, 1er étage 75005 PARIS FRANCE Tel: 33 (0)1 40 79 57 64 Fax: 33 (0)1 40 79 38 35 Email: [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 -- Dr Renaud Lancelot, vétérinaire Projet FSP régional épidémiologie vétérinaire C/0 Ambassade de France - SCAC BP 834 Antananarivo 101 - Madagascar e-mail: [EMAIL PROTECTED] tel.: +261 32 40 165 53 (cell) +261 20 22 665 36 ext. 225 (work) __ 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
Re: [R] help with unknown function
Agnes Gault wrote: Hello I am working on radio tracking data, with a short programme someone gave me and ... which should, supposedly, work ... In this programme, there is the function : getareahr(kern, levels = 95). But i cannot find any 'getareahr' in R ... could anyone help me? Don't know. I'd rather ask someone. ;-) I mean the same someone you got the code from, because getareahr() might be his/her private invention. Uwe Ligges thanks! Agnès --- Agnès GAULT graduate student UMR 5173 MNHN-CNRS-UPMC case postale 50 Species Conservation, Restoration and Population Survey (CERSP) 61 rue Buffon, 1er étage 75005 PARIS FRANCE Tel: 33 (0)1 40 79 57 64 Fax: 33 (0)1 40 79 38 35 Email: [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 __ 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
Re: [R] how to draw an ellipse
Le 18.08.2005 15:45, Guillaume Allain a écrit : That's probably a stupid question, but I'm looking for a low-level command which plots ellipse, specifying only center and deformation axes. The purpose is to illustrate bivariates gaussians with 2D .95 confident levels without using any specific library. Thanxs for your help, Guillaume Hello, Why don't you want to use any specific library ? You can't reinvent the wheel !! There is a package ellipse on CRAN which will do what you are looking for. Have you tried RSiteSearch(ellipse) Cheers, Romain -- visit the R Graph Gallery : http://addictedtor.free.fr/graphiques ~ ~~ Romain FRANCOIS - http://addictedtor.free.fr ~~ Etudiant ISUP - CS3 - Industrie et Services ~~http://www.isup.cicrp.jussieu.fr/ ~~ Stagiaire INRIA Futurs - Equipe SELECT ~~ http://www.inria.fr/recherche/equipes/select.fr.html~~ ~ __ 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] axTicks and window resizing
Dear listers,
I have written a function to facilitate the drawing of altitude profiles
with x (distance), y (altitude) and a z parameter (altitude magnification).
profplot-function(x,y,z=10,...){
op - par()$mai
par(mai=c(0.95625,0.76875,0.76875,0.95625))
plot(x,y*z, type=l,asp=1,las=1,xlab=,ylab=,yaxt=n,...)
axis(2,labels=axTicks(2)/z,las=1)
axis(4,labels=axTicks(2)/z,las=1)
on.exit(par(mai=op))
}
This worked apparently well until I had to resize the graphical window
after plotting. In this case, I get this message:
profplot(prof$dist,prof$alt,col=blue)
Erreur : les longueurs de 'at' et de 'label' diffèrent, 7 != 8
Which means Error: length of 'at' and label' differ, 7!=8 (whish R
2.1.1 could be parametrise 'English' even with a French Windows XP)
At this stage, R crashes (= I cannot get the graphic window
working/resized and must interrupt the process from Windows XP, then
restart R for good work with the graphical window).
The error occur with the difference between the tick number computed
from plot() and the one computed with axTicks(). If still equal (slight
resizing) everything goes smoothly.
Thanks for any comments, even rude... (I am not sure the
problem/programme has been tackled relevantly enough)
Patrick
__
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] How do I make a Sweave + latex table out of this ?
Dear list, I have a table that I would like to convert to latex for inclusion into a Sweave file. round(ftable(prop.table(xtabs(~agemF + votcat + Type , data=work),margin=2))*100,1) Type Voiced Voiceless unaspirated Voiceless aspirated agemF votcat 18 - 24 Prevoiced 2.6 8.7 2.3 Short lag 5.8 6.7 5.1 Long lag 1.0 1.9 2.9 24 - 30 Prevoiced15.1 10.5 1.7 Short lag 9.2 15.3 5.8 Long lag 3.5 8.115.8 30 - 36 Prevoiced12.8 14.0 2.6 Short lag10.2 14.2 3.0 Long lag 2.3 5.522.2 36 - 42 Prevoiced 4.4 6.4 0.6 Short lag 4.0 5.9 1.5 Long lag 1.3 2.9 9.4 42 - 48 Prevoiced 6.4 2.3 0.3 Short lag 3.0 2.8 1.4 Long lag 0.6 7.7 8.8 48 - 54 Prevoiced 4.9 4.1 0.3 Short lag 2.0 2.7 1.3 Long lag 0.3 0.9 4.7 However, I have not been able to use this as a table. The Hmisc latex command only accepts the input if I first convert it to a data.frame format, and that makes the output much more difficult to read as it duplicates the category levels of agemF. Is there a way to do this? /Fredrik Karlsson __ 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
Re: [R] do glm with two data sets
Thanks for your help. # read the two data sets e - as.matrix(read.table(file1.txt, header=TRUE,row.names=1)) g - as.matrix(read.table(file2.txt, header=TRUE,row.names=1)) # solution d1-data.frame(g[1,], e[1,]) fit-glm(e[1,] ~ g[1,], data=d1) summary(fit) I am not sure that is the best solution. Thanks again, Ying -Original Message- From: Gavin Simpson [mailto:[EMAIL PROTECTED] Sent: Wednesday, August 17, 2005 7:01 PM To: Sundar Dorai-Raj Cc: Hu, Ying (NIH/NCI); r-help@stat.math.ethz.ch Subject: Re: [R] do glm with two data sets On Wed, 2005-08-17 at 17:22 -0500, Sundar Dorai-Raj wrote: Hu, Ying (NIH/NCI) wrote: I have two data sets: File1.txt: Name id1 id2 id3 ... N10 1 0 ... N20 1 1 ... N31 1 -1... ... File2.txt: Group id1 id2 id3 ... G1 1.22 1.34 2.44 ... G2 2.33 2.56 2.56 ... G3 1.56 1.99 1.46 ... ... I like to do: x1-c(0,1,0,...) y1-c(1.22,1.34, 2.44, ...) z1-data.frame(x,y) summary(glm(y1~x1,data=z1) But I do the same thing by inputting the data sets from the two files e - read.table(file1.txt, header=TRUE,row.names=1) g - read.table(file2.txt, header=TRUE,row.names=1) e1-exp[1,] g1-geno[1,] d1-data.frame(g, e) summary(glm(e1 ~ g1, data=d1)) the error message is Error in model.frame(formula, rownames, variables, varnames, extras, extranames, : invalid variable type Execution halted Thanks in advance, Ying Hi Ying, That error message is likely caused by having a data.frame on the right hand side (rhs) of the formula. You can't have a data.frame on the rhs of a formula and g1 is still a data frame even if you only choose the first row, e.g.: dat - as.data.frame(matrix(100, 10, 10)) class(dat[1, ]) [1] data.frame You could try: glm(e1 ~ ., data=g1[1, ]) and see if that works, but as Sundar notes, your post is a little difficult to follow, so this may not do what you were trying to achieve. HTH Gav You have several inconsistencies in your example, so it will be difficult to figure out what you are trying to accomplish. e - read.table(file1.txt, header=TRUE,row.names=1) g - read.table(file2.txt, header=TRUE,row.names=1) e1-exp[1,] What's exp? Also it's dangerous to use an R function as a variable name. Most of the time R can tell the difference, but in some cases it cannot. g1-geno[1,] What's geno? d1-data.frame(g, e) d1 is now e and g cbind'ed together? summary(glm(e1 ~ g1, data=d1)) Are e1 and g1 elements of d1? From what you've told us, I don't know where the error is occurring. Also, if you are having errors, you can more easily isolate the problem by doing: fit - glm(e1 ~ g1, data = d1) summary(fit) This will at least tell you the problem is in your call to glm and not summary.glm. --sundar P.S. Please (re-)read the POSTING GUIDE. Most of the time you will figure out problems such as these on your own during the process of creating a reproducible example. __ 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 -- %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~% Gavin Simpson [T] +44 (0)20 7679 5522 ENSIS Research Fellow [F] +44 (0)20 7679 7565 ENSIS Ltd. ECRC [E] gavin.simpsonATNOSPAMucl.ac.uk UCL Department of Geography [W] http://www.ucl.ac.uk/~ucfagls/cv/ 26 Bedford Way[W] http://www.ucl.ac.uk/~ucfagls/ London. WC1H 0AP. %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~% __ 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
Re: [R] lme model: Error in MEEM
Christoph Lehmann christoph.lehmann at gmx.ch writes: Hi, We have data of two groups of subjects: 32 elderly, 14 young adults. for each subject we have 15 observations, each observation consisting of a reaction-time measure (RT) and an activation maesure (betadlpcv). since we want to analyze the influence of (age-)group and RT on the activation, we call: lme(betadlpcv ~ RT*group, data=our.data, random=~ RT |subject) this yields: Error in MEEM(object, conLin, control$niterEM) : Singularity in backsolve at level 0, block 1 If you really have 15 observations (690 lines) it should be enough to estimate the model (see below). Assume you had some degenerate case. From a psychophysical point of view, I am surprised that reaction time is on the right side, but that's off-subject. Dieter sub = data.frame(subject=1:46,group=c(rep(old,32),rep(young,14))) sub$slope = 2.5+as.numeric(sub$group)+rnorm(46,0.5) beta = data.frame( subject=rep(sub$subject,15), group=rep(sub$group,15), slope=rep(sub$slope,15), RT=rnorm(46*15,100,20)) beta$betadlpcv = beta$slope*beta$RT + rnorm(46*15,0.1) library(nlme) beta.lme = lme(betadlpcv ~ RT*group, data=beta, random=~ RT |subject) summary(beta.lme) __ 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
Re: [R] do glm with two data sets
Hu, Ying (NIH/NCI) wrote: Thanks for your help. # read the two data sets e - as.matrix(read.table(file1.txt, header=TRUE,row.names=1)) g - as.matrix(read.table(file2.txt, header=TRUE,row.names=1)) # solution d1-data.frame(g[1,], e[1,]) fit-glm(e[1,] ~ g[1,], data=d1) summary(fit) I am not sure that is the best solution. Thanks again, Ying Hi, Ying, What's wrong with this solution? Do you still get an error? What is your primary goal? A couple of points: 1. It's better to use names in your data.frame: d1 - data.frame(g = g[1,], e = e[1,]) Then in glm: fit - glm(e ~ g, data = d1) 2. Also, you may just be giving us a toy example, but if you don't specify a family argument in glm then you are simply getting the least squares. In that case you should use ?lm instead. HTH, --sundar -Original Message- From: Gavin Simpson [mailto:[EMAIL PROTECTED] Sent: Wednesday, August 17, 2005 7:01 PM To: Sundar Dorai-Raj Cc: Hu, Ying (NIH/NCI); r-help@stat.math.ethz.ch Subject: Re: [R] do glm with two data sets On Wed, 2005-08-17 at 17:22 -0500, Sundar Dorai-Raj wrote: Hu, Ying (NIH/NCI) wrote: I have two data sets: File1.txt: Name id1 id2 id3 ... N10 1 0 ... N20 1 1 ... N31 1 -1... ... File2.txt: Group id1 id2 id3 ... G1 1.22 1.34 2.44 ... G2 2.33 2.56 2.56 ... G3 1.56 1.99 1.46 ... ... I like to do: x1-c(0,1,0,...) y1-c(1.22,1.34, 2.44, ...) z1-data.frame(x,y) summary(glm(y1~x1,data=z1) But I do the same thing by inputting the data sets from the two files e - read.table(file1.txt, header=TRUE,row.names=1) g - read.table(file2.txt, header=TRUE,row.names=1) e1-exp[1,] g1-geno[1,] d1-data.frame(g, e) summary(glm(e1 ~ g1, data=d1)) the error message is Error in model.frame(formula, rownames, variables, varnames, extras, extranames, : invalid variable type Execution halted Thanks in advance, Ying Hi Ying, That error message is likely caused by having a data.frame on the right hand side (rhs) of the formula. You can't have a data.frame on the rhs of a formula and g1 is still a data frame even if you only choose the first row, e.g.: dat - as.data.frame(matrix(100, 10, 10)) class(dat[1, ]) [1] data.frame You could try: glm(e1 ~ ., data=g1[1, ]) and see if that works, but as Sundar notes, your post is a little difficult to follow, so this may not do what you were trying to achieve. HTH Gav You have several inconsistencies in your example, so it will be difficult to figure out what you are trying to accomplish. e - read.table(file1.txt, header=TRUE,row.names=1) g - read.table(file2.txt, header=TRUE,row.names=1) e1-exp[1,] What's exp? Also it's dangerous to use an R function as a variable name. Most of the time R can tell the difference, but in some cases it cannot. g1-geno[1,] What's geno? d1-data.frame(g, e) d1 is now e and g cbind'ed together? summary(glm(e1 ~ g1, data=d1)) Are e1 and g1 elements of d1? From what you've told us, I don't know where the error is occurring. Also, if you are having errors, you can more easily isolate the problem by doing: fit - glm(e1 ~ g1, data = d1) summary(fit) This will at least tell you the problem is in your call to glm and not summary.glm. --sundar P.S. Please (re-)read the POSTING GUIDE. Most of the time you will figure out problems such as these on your own during the process of creating a reproducible example. __ 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-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
Re: [R] do glm with two data sets
On Thu, 2005-08-18 at 10:38 -0400, Hu, Ying (NIH/NCI) wrote: Thanks for your help. # read the two data sets e - as.matrix(read.table(file1.txt, header=TRUE,row.names=1)) g - as.matrix(read.table(file2.txt, header=TRUE,row.names=1)) # solution d1-data.frame(g[1,], e[1,]) This is redundant, as: fit-glm(e[1,] ~ g[1,], data=d1) and: fit - glm(e[1, ] ~ g[1, ]) are equivalent - you don't need data = d1 in this case, e.g: e - matrix(c(0, 1, 0, 0, 1, 1, 1, 1, -1), ncol = 3, byrow = TRUE) e g - matrix(c(1.22, 1.34, 2.44, 2.33, 2.56, 2.56, 1.56, 1.99, 1.46), ncol = 3, byrow = TRUE) g fit - glm(e[1, ] ~ g[1, ]) fit works fine. summary(fit) I am not sure that is the best solution. This seems a strange way of doing this. Why not: pred - g[1, ] resp - e[1, ] fit - glm(resp ~ pred) fit and do your subsetting outside the glm call - makes things clearer no? Unless you plan to do many glm()s one per row of your two matrices. If that is the case, then there are better ways of approaching this. Thanks again, Ying HTH G -Original Message- From: Gavin Simpson [mailto:[EMAIL PROTECTED] Sent: Wednesday, August 17, 2005 7:01 PM To: Sundar Dorai-Raj Cc: Hu, Ying (NIH/NCI); r-help@stat.math.ethz.ch Subject: Re: [R] do glm with two data sets On Wed, 2005-08-17 at 17:22 -0500, Sundar Dorai-Raj wrote: Hu, Ying (NIH/NCI) wrote: I have two data sets: File1.txt: Name id1 id2 id3 ... N10 1 0 ... N20 1 1 ... N31 1 -1... ... File2.txt: Group id1 id2 id3 ... G1 1.22 1.34 2.44 ... G2 2.33 2.56 2.56 ... G3 1.56 1.99 1.46 ... ... I like to do: x1-c(0,1,0,...) y1-c(1.22,1.34, 2.44, ...) z1-data.frame(x,y) summary(glm(y1~x1,data=z1) But I do the same thing by inputting the data sets from the two files e - read.table(file1.txt, header=TRUE,row.names=1) g - read.table(file2.txt, header=TRUE,row.names=1) e1-exp[1,] g1-geno[1,] d1-data.frame(g, e) summary(glm(e1 ~ g1, data=d1)) the error message is Error in model.frame(formula, rownames, variables, varnames, extras, extranames, : invalid variable type Execution halted Thanks in advance, Ying Hi Ying, That error message is likely caused by having a data.frame on the right hand side (rhs) of the formula. You can't have a data.frame on the rhs of a formula and g1 is still a data frame even if you only choose the first row, e.g.: dat - as.data.frame(matrix(100, 10, 10)) class(dat[1, ]) [1] data.frame You could try: glm(e1 ~ ., data=g1[1, ]) and see if that works, but as Sundar notes, your post is a little difficult to follow, so this may not do what you were trying to achieve. HTH Gav You have several inconsistencies in your example, so it will be difficult to figure out what you are trying to accomplish. e - read.table(file1.txt, header=TRUE,row.names=1) g - read.table(file2.txt, header=TRUE,row.names=1) e1-exp[1,] What's exp? Also it's dangerous to use an R function as a variable name. Most of the time R can tell the difference, but in some cases it cannot. g1-geno[1,] What's geno? d1-data.frame(g, e) d1 is now e and g cbind'ed together? summary(glm(e1 ~ g1, data=d1)) Are e1 and g1 elements of d1? From what you've told us, I don't know where the error is occurring. Also, if you are having errors, you can more easily isolate the problem by doing: fit - glm(e1 ~ g1, data = d1) summary(fit) This will at least tell you the problem is in your call to glm and not summary.glm. --sundar P.S. Please (re-)read the POSTING GUIDE. Most of the time you will figure out problems such as these on your own during the process of creating a reproducible example. __ 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 -- %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~% Gavin Simpson [T] +44 (0)20 7679 5522 ENSIS Research Fellow [F] +44 (0)20 7679 7565 ENSIS Ltd. ECRC [E] gavin.simpsonATNOSPAMucl.ac.uk UCL Department of Geography [W] http://www.ucl.ac.uk/~ucfagls/cv/ 26 Bedford Way[W] http://www.ucl.ac.uk/~ucfagls/ London. WC1H 0AP. %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~% __ 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
Re: [R] how to draw an ellipse
Romain Francois wrote: Why don't you want to use any specific library ? You can't reinvent the wheel !! True, but Madonna is always reinventing herself... There is a package ellipse on CRAN which will do what you are looking for. Have you tried RSiteSearch(ellipse) its simple geometry really: theta=seq(0,2*pi,len=100) e=1.5 r=3 x=e*r*cos(theta) y=r*sin(theta) plot(x,y,asp=1) rotate with: phi=pi/4 xr=x*sin(phi)+y*cos(phi) yr=-x*cos(phi)+y*sin(phi) plot(xr,yr,asp=1) or something. Wrap that up into a function and you're done. This is off-the-cuff, I've probably messed something up. So use one prepared earlier from a library... Baz __ 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
Re: [R] how to draw an ellipse
Congratulations guys, you made my day better and funnier! Le 18 août 05, à 17:31, Barry Rowlingson a écrit : Romain Francois wrote: Why don't you want to use any specific library ? You can't reinvent the wheel !! True, but Madonna is always reinventing herself... There is a package ellipse on CRAN which will do what you are looking for. Have you tried RSiteSearch(ellipse) its simple geometry really: theta=seq(0,2*pi,len=100) e=1.5 r=3 x=e*r*cos(theta) y=r*sin(theta) plot(x,y,asp=1) rotate with: phi=pi/4 xr=x*sin(phi)+y*cos(phi) yr=-x*cos(phi)+y*sin(phi) plot(xr,yr,asp=1) or something. Wrap that up into a function and you're done. This is off-the-cuff, I've probably messed something up. So use one prepared earlier from a library... Baz __ Guillaume Allain Carte Blanche Conseil 47 rue de Lancry 75010 Tel : 01 42412121 Mail : [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
[R] How to do log normal regression?
I want to fit a Log-Normal CDF function between two variables, and estimate the parameters. Is there any package/functions designed for this purpose? Basically, I have data for Y and X, and I suspect the relationship between Y and X is Y = CDF Log-Normal (X), and I want to run this regression to verify this and estimate the parameters. Anyone has any thoughts? Any input is valuable to me, so please do not hesitate to share your thoughts. Thank you! Ed. __ 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
Re: [R] display of a loess fitted surface
On Thu, 18 Aug 2005, Marta Colombo wrote: I am Marta Colombo,student at Politecnico,Milan. I am studying local regression models and I am using loess function. My problem is that when I have a loess object I don't know how to display the fitted surface; in fact, while in S when you have a loess object you can see it writing plot(object), in R this dosen't work. Why not use S if that does what you want? There are examples using R in MASS, for example: see the ch04.R and ch15.R scripts. Also I'd like to know if there is something like the S function pointwise that computes upper and lower confidence intervals. (Do you mean upper and lower confidence limits?) Thats a more general question. Many of R's predict() methods can construct confidence (and tolerance) intervals. If not, it is easy to do this yourself, as pred_obj$fit-/+pred_obj$se*qnorm(alpha/2). -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ 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] standard errors for expression intensities
Hello! When I use the functions rma() or justrma() in the Bioconductor affy package, the standard errors for the expression estimates computed by the function se.exprs() is a matrix of NA's. I am wondering if any of you have encountered this problem and if there is a solution. Any help would be appreciated. Thanks, Karthik. [[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
[R] Regular expressions sub
Dear all, I am struggling with the use of regular expression. I got as.character(test$sample.id) [1] 1.11 10.11 11.11 113.31 114.2 114.3 114.8 and need [1] 11 11 11 31 2 3 8 I.e. remove everything before the . . TIA, Bernd __ 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
Re: [R] select previous date
You can try this: dates - df$date[df$temp==15] one.month.before - sapply(strsplit(dates, -), function(x) paste(x[1], sprintf(%02d, as.numeric(x[2])-1), x[3], sep=-)) df[df$date %in% one.month.before,] Ben On 8/18/05, alessandro carletti [EMAIL PROTECTED] wrote: Hi everybody, could anyone help me in finding a way for selecting from a dataframe each row that comes before another that matches a condition? EXAMPLE: df raw.number name date temp 1aaa 2001-04-15 15 2aaa 2001-01-15 12 3aaa 2001-03-15 13 ... i-1 bbb 2002-04-15 15 ibbb 2002-03-15 14 the condition is: df$temp==15 matching raws are 1 and i-1: I need something to select (only) rows where date=one month before the matching raws, so raws 3 and i. (the variable name has more than one level) 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 __ 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
Re: [R] axTicks and window resizing
On Thu, 18 Aug 2005, Patrick Giraudoux wrote:
Dear listers,
I have written a function to facilitate the drawing of altitude profiles
with x (distance), y (altitude) and a z parameter (altitude magnification).
profplot-function(x,y,z=10,...){
op - par()$mai
par(mai=c(0.95625,0.76875,0.76875,0.95625))
plot(x,y*z, type=l,asp=1,las=1,xlab=,ylab=,yaxt=n,...)
axis(2,labels=axTicks(2)/z,las=1)
axis(4,labels=axTicks(2)/z,las=1)
on.exit(par(mai=op))
}
This worked apparently well until I had to resize the graphical window
after plotting. In this case, I get this message:
profplot(prof$dist,prof$alt,col=blue)
Erreur : les longueurs de 'at' et de 'label' diffèrent, 7 != 8
Which means Error: length of 'at' and label' differ, 7!=8 (whish R
2.1.1 could be parametrise 'English' even with a French Windows XP)
If I understand you correctly, it can. Just add LANGUAGE=en to the
shortcut.
At this stage, R crashes (= I cannot get the graphic window
working/resized and must interrupt the process from Windows XP, then
restart R for good work with the graphical window).
The error occur with the difference between the tick number computed
from plot() and the one computed with axTicks(). If still equal (slight
resizing) everything goes smoothly.
The problem is that you need to specify 'at' and 'labels' to axis(): you
cannot safely specify just labels. When re-drawing, 'at' is recomputed,
but your specification of 'labels' is not.
I suspect you can just do dev.off() and open a new graphics window.
Thanks for any comments, even rude... (I am not sure the
problem/programme has been tackled relevantly enough)
Patrick
__
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
--
Brian D. Ripley, [EMAIL PROTECTED]
Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel: +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax: +44 1865 272595__
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
Re: [R] do glm with two data sets
You are right.
# read the two data sets
e - as.matrix(read.table(file1.txt, header=TRUE,row.names=1))
g - as.matrix(read.table(file2.txt, header=TRUE,row.names=1))
# solution 2
summary(glm(e[1,] ~ g[1,]))
summary(glm(e[1,] ~ g[2,]))
...
They work very well.
If I put it in the loop, such as
for (i in 1:50){
for (j in 1:50){
cat(file1 row:, i, file2 row:, j, \n)
print(summary(glm(e[i,] ~ g[j,])))
}
}
Why do I have to use print to print the results? If without print
for (i in 1:50){
for (j in 1:50){
cat(file1 row:, i, file2 row:, j, \n)
summary(glm(e[i,] ~ g[j,]))
}
}
then without the results of glm.
Thanks a lot.
Ying
-Original Message-
From: Gavin Simpson [mailto:[EMAIL PROTECTED]
Sent: Thursday, August 18, 2005 11:00 AM
To: Hu, Ying (NIH/NCI)
Cc: Sundar Dorai-Raj; r-help@stat.math.ethz.ch
Subject: RE: [R] do glm with two data sets
On Thu, 2005-08-18 at 10:38 -0400, Hu, Ying (NIH/NCI) wrote:
Thanks for your help.
# read the two data sets
e - as.matrix(read.table(file1.txt, header=TRUE,row.names=1))
g - as.matrix(read.table(file2.txt, header=TRUE,row.names=1))
# solution
d1-data.frame(g[1,], e[1,])
This is redundant, as:
fit-glm(e[1,] ~ g[1,], data=d1)
and:
fit - glm(e[1, ] ~ g[1, ])
are equivalent - you don't need data = d1 in this case, e.g:
e - matrix(c(0, 1, 0, 0, 1, 1, 1, 1, -1), ncol = 3, byrow = TRUE)
e
g - matrix(c(1.22, 1.34, 2.44, 2.33, 2.56, 2.56, 1.56, 1.99, 1.46),
ncol = 3, byrow = TRUE)
g
fit - glm(e[1, ] ~ g[1, ])
fit
works fine.
summary(fit)
I am not sure that is the best solution.
This seems a strange way of doing this. Why not:
pred - g[1, ]
resp - e[1, ]
fit - glm(resp ~ pred)
fit
and do your subsetting outside the glm call - makes things clearer no?
Unless you plan to do many glm()s one per row of your two matrices. If
that is the case, then there are better ways of approaching this.
Thanks again,
Ying
HTH
G
-Original Message-
From: Gavin Simpson [mailto:[EMAIL PROTECTED]
Sent: Wednesday, August 17, 2005 7:01 PM
To: Sundar Dorai-Raj
Cc: Hu, Ying (NIH/NCI); r-help@stat.math.ethz.ch
Subject: Re: [R] do glm with two data sets
On Wed, 2005-08-17 at 17:22 -0500, Sundar Dorai-Raj wrote:
Hu, Ying (NIH/NCI) wrote:
I have two data sets:
File1.txt:
Name id1 id2 id3 ...
N10 1 0 ...
N20 1 1 ...
N31 1 -1...
...
File2.txt:
Group id1 id2 id3 ...
G1 1.22 1.34 2.44 ...
G2 2.33 2.56 2.56 ...
G3 1.56 1.99 1.46 ...
...
I like to do:
x1-c(0,1,0,...)
y1-c(1.22,1.34, 2.44, ...)
z1-data.frame(x,y)
summary(glm(y1~x1,data=z1)
But I do the same thing by inputting the data sets from the two files
e - read.table(file1.txt, header=TRUE,row.names=1)
g - read.table(file2.txt, header=TRUE,row.names=1)
e1-exp[1,]
g1-geno[1,]
d1-data.frame(g, e)
summary(glm(e1 ~ g1, data=d1))
the error message is
Error in model.frame(formula, rownames, variables, varnames, extras,
extranames, :
invalid variable type
Execution halted
Thanks in advance,
Ying
Hi Ying,
That error message is likely caused by having a data.frame on the right
hand side (rhs) of the formula. You can't have a data.frame on the rhs
of a formula and g1 is still a data frame even if you only choose the
first row, e.g.:
dat - as.data.frame(matrix(100, 10, 10))
class(dat[1, ])
[1] data.frame
You could try:
glm(e1 ~ ., data=g1[1, ])
and see if that works, but as Sundar notes, your post is a little
difficult to follow, so this may not do what you were trying to achieve.
HTH
Gav
You have several inconsistencies in your example, so it will be
difficult to figure out what you are trying to accomplish.
e - read.table(file1.txt, header=TRUE,row.names=1)
g - read.table(file2.txt, header=TRUE,row.names=1)
e1-exp[1,]
What's exp? Also it's dangerous to use an R function as a variable
name. Most of the time R can tell the difference, but in some cases it
cannot.
g1-geno[1,]
What's geno?
d1-data.frame(g, e)
d1 is now e and g cbind'ed together?
summary(glm(e1 ~ g1, data=d1))
Are e1 and g1 elements of d1? From what you've told us, I don't
know where the error is occurring. Also, if you are having errors, you
can more easily isolate the problem by doing:
fit - glm(e1 ~ g1, data = d1)
summary(fit)
This will at least tell you the problem is in your call to glm and not
summary.glm.
--sundar
P.S. Please (re-)read the POSTING GUIDE. Most of the time you will
figure out problems such as these on your own during the process of
creating a reproducible example.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
Re: [R] Regular expressions sub
x - scan(clipboard, what=) Read 7 items x [1] 1.11 10.11 11.11 113.31 114.2 114.3 114.8 gsub([0-9]*\\., , x) [1] 11 11 11 31 2 3 8 Bernd Weiss wrote: Dear all, I am struggling with the use of regular expression. I got as.character(test$sample.id) [1] 1.11 10.11 11.11 113.31 114.2 114.3 114.8 and need [1] 11 11 11 31 2 3 8 I.e. remove everything before the . . TIA, Bernd __ 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-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
Re: [R] Regular expressions sub
One solution is
test - c(1.11,10.11,11.11,113.31,114.2,114.3)
id - unlist(lapply(strsplit(test,[.]),function(x) {x[2]}))
-Original Message-
From: Bernd Weiss [mailto:[EMAIL PROTECTED]
Sent: Thursday, August 18, 2005 12:10 PM
To: r-help@stat.math.ethz.ch
Subject: [R] Regular expressions sub
Dear all,
I am struggling with the use of regular expression. I got
as.character(test$sample.id)
[1] 1.11 10.11 11.11 113.31 114.2 114.3 114.8
and need
[1] 11 11 11 31 2 3 8
I.e. remove everything before the . .
TIA,
Bernd
__
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-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
Re: [R] do glm with two data sets
Hu, Ying (NIH/NCI) wrote:
You are right.
# read the two data sets
e - as.matrix(read.table(file1.txt, header=TRUE,row.names=1))
g - as.matrix(read.table(file2.txt, header=TRUE,row.names=1))
# solution 2
summary(glm(e[1,] ~ g[1,]))
summary(glm(e[1,] ~ g[2,]))
...
They work very well.
If I put it in the loop, such as
for (i in 1:50){
for (j in 1:50){
cat(file1 row:, i, file2 row:, j, \n)
print(summary(glm(e[i,] ~ g[j,])))
}
}
Why do I have to use print to print the results? If without print
for (i in 1:50){
for (j in 1:50){
cat(file1 row:, i, file2 row:, j, \n)
summary(glm(e[i,] ~ g[j,]))
}
}
then without the results of glm.
This is a FAQ 7.16.
See http://cran.r-project.org/doc/FAQ/R-FAQ.html
--sundar
Thanks a lot.
Ying
-Original Message-
From: Gavin Simpson [mailto:[EMAIL PROTECTED]
Sent: Thursday, August 18, 2005 11:00 AM
To: Hu, Ying (NIH/NCI)
Cc: Sundar Dorai-Raj; r-help@stat.math.ethz.ch
Subject: RE: [R] do glm with two data sets
On Thu, 2005-08-18 at 10:38 -0400, Hu, Ying (NIH/NCI) wrote:
Thanks for your help.
# read the two data sets
e - as.matrix(read.table(file1.txt, header=TRUE,row.names=1))
g - as.matrix(read.table(file2.txt, header=TRUE,row.names=1))
# solution
d1-data.frame(g[1,], e[1,])
This is redundant, as:
fit-glm(e[1,] ~ g[1,], data=d1)
and:
fit - glm(e[1, ] ~ g[1, ])
are equivalent - you don't need data = d1 in this case, e.g:
e - matrix(c(0, 1, 0, 0, 1, 1, 1, 1, -1), ncol = 3, byrow = TRUE)
e
g - matrix(c(1.22, 1.34, 2.44, 2.33, 2.56, 2.56, 1.56, 1.99, 1.46),
ncol = 3, byrow = TRUE)
g
fit - glm(e[1, ] ~ g[1, ])
fit
works fine.
summary(fit)
I am not sure that is the best solution.
This seems a strange way of doing this. Why not:
pred - g[1, ]
resp - e[1, ]
fit - glm(resp ~ pred)
fit
and do your subsetting outside the glm call - makes things clearer no?
Unless you plan to do many glm()s one per row of your two matrices. If
that is the case, then there are better ways of approaching this.
Thanks again,
Ying
HTH
G
-Original Message-
From: Gavin Simpson [mailto:[EMAIL PROTECTED]
Sent: Wednesday, August 17, 2005 7:01 PM
To: Sundar Dorai-Raj
Cc: Hu, Ying (NIH/NCI); r-help@stat.math.ethz.ch
Subject: Re: [R] do glm with two data sets
On Wed, 2005-08-17 at 17:22 -0500, Sundar Dorai-Raj wrote:
Hu, Ying (NIH/NCI) wrote:
I have two data sets:
File1.txt:
Name id1 id2 id3 ...
N10 1 0 ...
N20 1 1 ...
N31 1 -1...
...
File2.txt:
Group id1 id2 id3 ...
G1 1.22 1.34 2.44 ...
G2 2.33 2.56 2.56 ...
G3 1.56 1.99 1.46 ...
...
I like to do:
x1-c(0,1,0,...)
y1-c(1.22,1.34, 2.44, ...)
z1-data.frame(x,y)
summary(glm(y1~x1,data=z1)
But I do the same thing by inputting the data sets from the two files
e - read.table(file1.txt, header=TRUE,row.names=1)
g - read.table(file2.txt, header=TRUE,row.names=1)
e1-exp[1,]
g1-geno[1,]
d1-data.frame(g, e)
summary(glm(e1 ~ g1, data=d1))
the error message is
Error in model.frame(formula, rownames, variables, varnames, extras,
extranames, :
invalid variable type
Execution halted
Thanks in advance,
Ying
Hi Ying,
That error message is likely caused by having a data.frame on the right
hand side (rhs) of the formula. You can't have a data.frame on the rhs
of a formula and g1 is still a data frame even if you only choose the
first row, e.g.:
dat - as.data.frame(matrix(100, 10, 10))
class(dat[1, ])
[1] data.frame
You could try:
glm(e1 ~ ., data=g1[1, ])
and see if that works, but as Sundar notes, your post is a little
difficult to follow, so this may not do what you were trying to achieve.
HTH
Gav
You have several inconsistencies in your example, so it will be
difficult to figure out what you are trying to accomplish.
e - read.table(file1.txt, header=TRUE,row.names=1)
g - read.table(file2.txt, header=TRUE,row.names=1)
e1-exp[1,]
What's exp? Also it's dangerous to use an R function as a variable
name. Most of the time R can tell the difference, but in some cases it
cannot.
g1-geno[1,]
What's geno?
d1-data.frame(g, e)
d1 is now e and g cbind'ed together?
summary(glm(e1 ~ g1, data=d1))
Are e1 and g1 elements of d1? From what you've told us, I don't
know where the error is occurring. Also, if you are having errors, you
can more easily isolate the problem by doing:
fit - glm(e1 ~ g1, data = d1)
summary(fit)
This will at least tell you the problem is in your call to glm and not
summary.glm.
--sundar
P.S. Please (re-)read the POSTING GUIDE. Most of the time you will
figure out problems such as these on your own during the process of
creating a reproducible example.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting
Re: [R] line/bar through median in lattice's bwplot?
On 8/18/05, Dieter Menne [EMAIL PROTECTED] wrote: M. K. mk_lists at yahoo.ca writes: Is there a way to render a line through the median point in the boxplot generated by the Lattice command bwplot? The line basically bisects the bar at the median point... bwplot(height~voice.part , pch='|', data=singer, xlab=Height (inches)) How to find this (haven't checked, maybe it's documented) library(lattice) panel.bwplot # no ()! It's not documented (and only available in the very recent versions of lattice, BTW), as I realized last night. It will be documented in a soon-to-come version (and maybe even honour attempts to change its color, line type, etc). Deepayan __ 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
Re: [R] problems when installing R in Fedora core 4
On Thu, 18 Aug 2005, Jonathan Baron wrote: ... Whether you have the needed files already depends on what kind of installation you did. Some of the packages are devel and others are compat. Here is my list of compat rpms that I have installed, and I think I installed all of these just to get R to build: compat-libf2c-32-3.2.3-47.fc4 compat-libstdc++-296-2.96-132.fc4 compat-readline43-4.3-2 compat-gcc-32-3.2.3-47.fc4 My hunch is that I still do not have the optimal installation, but it is possible that the newest versions of gcc have solved some of the problems with the ones that originally came with FC4. I've seen some discussion suggesting that the way to go is to use an older version of gcc, but I did not search for it just now. See https://stat.ethz.ch/pipermail/r-devel/2005-August/034171.html for updates on that advice. It may still be the way to go. -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ 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
Re: [R] Regular expressions sub
Bernd Weiss bernd.weiss at uni-koeln.de writes: I am struggling with the use of regular expression. I got as.character(test$sample.id) [1] 1.11 10.11 11.11 113.31 114.2 114.3 114.8 and need [1] 11 11 11 31 2 3 8 I.e. remove everything before the . . Define the dot as the hard separator, and allow for multiple digits before it: sample.id - c(1.11, 10.11, 11.11, 113.31, 114.2, 114.3, 114.8) gsub(^[0-9]*\., , sample.id) [1] 11 11 11 31 2 3 8 Hope this helps, Dirk __ 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
Re: [R] line/bar through median in lattice's
Dieter Menne dieter.menne at menne-biomed.de writes: Check the code for something like points (which is the default). Find a mysterious '|' and pch Dieter Yes, that was it! Thanks Dieter! The mysterious code is right near the bottom of the function. Too bad this is not documented, at least not in ?panel.bwplot. __ 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] Set R 2.1.1. in English
(whish R 2.1.1 could be parametrise 'English' even with a French Windows XP) If I understand you correctly, it can. Just add LANGUAGE=en to the shortcut. Wonderful hope but not sure to catch what you term shortcut. I tried to add this command in C:\R\rw2011\etc\Rprofile, the .Rprofile in the folder my documents, but this cannot be understood from there... which obviously shows that shortcut is not a general term for profiles! I have also unsuccessfully looked for a file shortcut in the rw2011 folder. I also tried and went through the R-help archive. There were some exchanges on this subject and Asian languages some weeks ago but what I have tried and adapt on this basis did not work. The objective would be to have R in English (thus additionnally allowing MDI mode with R-WinEdt, which is not the case with any other language) and to keep Windows XP and other applications in the foreign (= French, here) language. Thanks for any further hint, Patrick Giraudoux __ 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
Re: [R] Welcome to the R-help mailing list (Digest mode)
Hi, thank you for cooperation. --- Dr. (UA) Tetyana Stepanchuk Department of Electronic Commerce Johann Wolfgang Goethe-University Mertonstrasse 17 D-60054, Frankfurt/Main Germany phone: +49 069 798 22379 http://www.ecommerce.wiwi.uni-frankfurt.de -Ursprüngliche Nachricht- Von: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Im Auftrag von [EMAIL PROTECTED] Gesendet: Donnerstag, 18. August 2005 19:00 An: [EMAIL PROTECTED] Betreff: Welcome to the R-help mailing list (Digest mode) Welcome to the R-help@stat.math.ethz.ch mailing list! To post to this list, send your email to: r-help@stat.math.ethz.ch General information about the mailing list is at: https://stat.ethz.ch/mailman/listinfo/r-help If you ever want to unsubscribe or change your options (eg, switch to or from digest mode, change your password, etc.), visit your subscription page at: https://stat.ethz.ch/mailman/options/r-help/stepanchuk%40wiwi.uni-frankfurt. de You can also make such adjustments via email by sending a message to: [EMAIL PROTECTED] with the word `help' in the subject or body (don't include the quotes), and you will get back a message with instructions. You must know your password to change your options (including changing the password, itself) or to unsubscribe. It is: stepkin78 Normally, Mailman will remind you of your stat.math.ethz.ch mailing list passwords once every month, although you can disable this if you prefer. This reminder will also include instructions on how to unsubscribe or change your account options. There is also a button on your options page that will email your current password to you. __ 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
Re: [R] axTicks and window resizing
OK. Things work now and the window can be resized easy after adding at =
axTicks() in axis() explicitely. This makes:
profplot-function(x,y,z=10,...){
op - par()$mai
par(mai=c(0.95625,0.76875,0.76875,0.95625))
plot(x,y*z, type=l,asp=1,las=1,xlab=,ylab=,yaxt=n,...)
axis(2,at=axTicks(2),labels=axTicks(2)/z,las=1)
axis(4,at=axTicks(2),labels=axTicks(2)/z,las=1)
par(mai=op)
}
Thanks for the hint,
Patrick
Prof Brian Ripley a écrit :
On Thu, 18 Aug 2005, Patrick Giraudoux wrote:
Dear listers,
I have written a function to facilitate the drawing of altitude profiles
with x (distance), y (altitude) and a z parameter (altitude
magnification).
profplot-function(x,y,z=10,...){
op - par()$mai
par(mai=c(0.95625,0.76875,0.76875,0.95625))
plot(x,y*z, type=l,asp=1,las=1,xlab=,ylab=,yaxt=n,...)
axis(2,labels=axTicks(2)/z,las=1)
axis(4,labels=axTicks(2)/z,las=1)
on.exit(par(mai=op))
}
This worked apparently well until I had to resize the graphical window
after plotting. In this case, I get this message:
profplot(prof$dist,prof$alt,col=blue)
Erreur : les longueurs de 'at' et de 'label' diffèrent, 7 != 8
Which means Error: length of 'at' and label' differ, 7!=8 (whish R
2.1.1 could be parametrise 'English' even with a French Windows XP)
If I understand you correctly, it can. Just add LANGUAGE=en to the
shortcut.
At this stage, R crashes (= I cannot get the graphic window
working/resized and must interrupt the process from Windows XP, then
restart R for good work with the graphical window).
The error occur with the difference between the tick number computed
from plot() and the one computed with axTicks(). If still equal (slight
resizing) everything goes smoothly.
The problem is that you need to specify 'at' and 'labels' to axis():
you cannot safely specify just labels. When re-drawing, 'at' is
recomputed, but your specification of 'labels' is not.
I suspect you can just do dev.off() and open a new graphics window.
Thanks for any comments, even rude... (I am not sure the
problem/programme has been tackled relevantly enough)
Patrick
__
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-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
Re: [R] Error messages using LMER
Dear Professor Bates, Here is output R 2.1.1 produced with control = list(EMverbose = TRUE, msVerbose = TRUE). I am getting the new devel version and see what will hapen there: EM iterations 0 85289.766 ( 5407.13: 0.0815) ( 26134.4: 0.00387) 1 84544.322 ( 333.732: 0.137) ( 1462.32: 0.00934) 2 84515.108 ( 129.506: 0.0270) ( 446.306: 0.00481) 3 84514.519 ( 115.592: 0.00355) ( 328.637: 0.00103) 4 84514.505 ( 113.981:0.000505) ( 311.160:0.000165) 5 84514.505 ( 113.755:7.45e-005) ( 308.524:2.50e-005) 6 84514.505 ( 113.722:1.11e-005) ( 308.128:3.77e-006) 7 84514.505 ( 113.717:1.66e-006) ( 308.068:5.66e-007) 8 84514.505 ( 113.716:2.50e-007) ( 308.059:8.49e-008) 9 84514.505 ( 113.716:3.74e-008) ( 308.058:1.27e-008) 10 84514.505 ( 113.716:5.62e-009) ( 308.058:1.91e-009) 11 84514.505 ( 113.716:8.43e-010) ( 308.058:2.87e-010) 12 84514.505 ( 113.716:1.27e-010) ( 308.058:4.31e-011) 13 84514.505 ( 113.716:1.90e-011) ( 308.057:6.47e-012) 14 84514.505 ( 113.716:2.86e-012) ( 308.057:9.73e-013) 15 84514.505 ( 113.716:4.25e-013) ( 308.057:1.44e-013) iter0 value 84514.505044 final value 84514.505044 converged EM iterations 0 83740.342 ( 113.716: -0.0164) ( 308.057:0.000596) 1 83740.273 ( 121.512:-0.00121) ( 298.914:-3.24e-005) 2 83740.272 ( 122.131:-0.000111) ( 299.397:-1.24e-005) iter0 value 83740.272232 final value 83740.272232 converged EM iterations 0 84011.550 ( 122.204:-0.00461) ( 299.576:-0.000256) 1 84011.543 ( 124.624:-0.000459) ( 303.453:-6.41e-005) 2 84011.543 ( 124.870:-5.42e-005) ( 304.440:-1.13e-005) iter0 value 84011.543350 final value 84011.543350 converged EM iterations 0 84018.592 ( 124.915:-6.44e-005) ( 304.548:-1.22e-005) 1 84018.592 ( 124.949:-8.29e-006) ( 304.737:-1.99e-006) 2 84018.592 ( 124.954:-1.15e-006) ( 304.768:-3.08e-007) iter0 value 84018.591624 final value 84018.591624 converged EM iterations 0 84018.612 ( 124.955:3.40e-007) ( 304.770:-1.98e-007) 1 84018.612 ( 124.955:-9.98e-009) ( 304.773:-2.23e-008) 2 84018.612 ( 124.955:-5.47e-009) ( 304.774:-2.86e-009) iter0 value 84018.611512 final value 84018.611512 converged Error in fn(par, ...) : Unable to invert singular factor of downdated X'X In addition: Warning message: Leading minor of size 8 of downdated X'X is indefinite Thanks! Shige On 8/18/05, Douglas Bates [EMAIL PROTECTED] wrote: On 8/18/05, Shige Song [EMAIL PROTECTED] wrote: Dear All, After playing with lmer for couple of days, I have to say that I am amazed! I've been using quite some multilevel/mixed modeling packages, lme4 is a strong candidate for the overall winner, especially for multilevel generzlized linear models. Now go back to my two-level poisson model with cross-classified model. I've been testing various different model specificatios for the past couple of days. Here are the models I tried: 1) Two level random intercept model with level-1 covariates only m1 - lmer(.D ~ offset(log(.Y)) + (1|provn) +(1|bcohort) + x1 + x2 , data, poisson, method=Laplace) 2) Two-level random intercept model with both level-1 and level-2 covariates, but no cross-level interactions: m2 - lmer(.D ~ offset(log(.Y)) + (1|provn) +(1|bcohort) + x1 + x2 + z1 + z2, data, poisson, method=Laplace) 3) Two-level random intercept with cross-level interaction m3 - lmer(.D ~ offset(log(.Y)) + (1|provn) +(1|bcohort) + x1 + x2 + z1 + z2 + x1:z1 + x2:z2, data, poisson, method=Laplace) Both model 1 and 2 run fine. For model 3, I got error message: -- Error in fn(par, ...) : Unable to invert singular factor of downdated X'X In addition: Warning messages: 1: optim or nlminb returned message ERROR: ABNORMAL_TERMINATION_IN_LNSRCH in: LMEopt(x = mer, value = cv) 2: Leading minor of size 1 of downdated X'X is indefinite -- What is going on here? Any workarounds? Thanks! The first thing I would try is set the EMverbose and msVerbose flags in the control list to see what occurs within the optimization. That is append the argument control = list(EMverbose = TRUE, msVerbose = TRUE) to your call to lmer(). You may also want to try the call in a recently compiled R-devel, which will be released as R-2.2.0 in October. You will notice that the first warning message reads optim or nlminb. In R-2.1.1 lmer uses optim for the optimization. Starting with R-2.2.0 the default is to use nlminb. Test compilations of R-devel for Windows are available from CRAN. __ 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
Re: [R] LyX and Sweave
I just wanted to point out that I was there first :) on the Lyx List (Nov 2004): http://www.mail-archive.com/lyx-users@lists.lyx.org/msg36262.html Perhaps somebody who is trying to put all of this together can benefit from both sets of explanations. pj [EMAIL PROTECTED] wrote: On Mon, 25 Jul 2005 14:12:41 +0200, Gorjanc Gregor (GG) wrote: Hello R-users! I have tried to use Sweave within LyX* and found two ways to accomplish this. I have attached LyX source file for both ways as well as generated PDFs. I have copied Gregor's files at http://www.ci.tuwien.ac.at/~leisch/Sweave/LyX for those who didn't get the attachments. LyX looks actually much better and stable then when I last had a look a couple of years ago. -- Paul E. Johnson email: [EMAIL PROTECTED] Dept. of Political Sciencehttp://lark.cc.ku.edu/~pauljohn 1541 Lilac Lane, Rm 504 University of Kansas Office: (785) 864-9086 Lawrence, Kansas 66044-3177 FAX: (785) 864-5700 __ 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
Re: [R] Error messages using LMER
Here is what happened using R-devel: - EM iterations 0 85289.766 ( 5407.13: 0.0815) ( 26134.4: 0.00387) 1 84544.322 ( 333.732: 0.137) ( 1462.32: 0.00934) 2 84515.108 ( 129.506: 0.0270) ( 446.306: 0.00481) 3 84514.519 ( 115.592: 0.00355) ( 328.637: 0.00103) 4 84514.505 ( 113.981:0.000505) ( 311.160:0.000165) 5 84514.505 ( 113.755:7.45e-005) ( 308.524:2.50e-005) 6 84514.505 ( 113.722:1.11e-005) ( 308.128:3.77e-006) 7 84514.505 ( 113.717:1.66e-006) ( 308.068:5.66e-007) 8 84514.505 ( 113.716:2.50e-007) ( 308.059:8.49e-008) 9 84514.505 ( 113.716:3.74e-008) ( 308.058:1.27e-008) 10 84514.505 ( 113.716:5.62e-009) ( 308.058:1.91e-009) 11 84514.505 ( 113.716:8.43e-010) ( 308.058:2.87e-010) 12 84514.505 ( 113.716:1.27e-010) ( 308.058:4.31e-011) 13 84514.505 ( 113.716:1.90e-011) ( 308.057:6.47e-012) 14 84514.505 ( 113.716:2.86e-012) ( 308.057:9.73e-013) 15 84514.505 ( 113.716:4.25e-013) ( 308.057:1.44e-013) EM iterations 0 83740.342 ( 113.716: -0.0164) ( 308.057:0.000596) 1 83740.273 ( 121.512:-0.00121) ( 298.914:-3.24e-005) 2 83740.272 ( 122.131:-0.000111) ( 299.397:-1.24e-005) EM iterations 0 84011.546 ( 122.131:-0.00474) ( 299.397:-0.000265) 1 84011.539 ( 124.616:-0.000472) ( 303.415:-6.62e-005) 2 84011.539 ( 124.869:-5.58e-005) ( 304.433:-1.16e-005) EM iterations 0 84018.589 ( 124.869:-0.000139) ( 304.433:-1.81e-005) 1 84018.589 ( 124.944:-1.62e-005) ( 304.713:-3.26e-006) 2 84018.589 ( 124.953:-2.12e-006) ( 304.764:-5.23e-007) EM iterations 0 84018.611 ( 124.953:-2.38e-006) ( 304.764:-5.44e-007) 1 84018.611 ( 124.954:-3.25e-007) ( 304.772:-8.50e-008) 2 84018.611 ( 124.955:-4.66e-008) ( 304.773:-1.29e-008) EM iterations 0 84018.611 ( 124.955:-4.75e-008) ( 304.773:-1.30e-008) 1 84018.611 ( 124.955:-6.93e-009) ( 304.774:-1.97e-009) 2 84018.611 ( 124.955:-1.03e-009) ( 304.774:-2.96e-010) Warning messages: 1: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 2: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 3: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 4: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 5: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 6: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 7: nlminb failed to converge in: lmer(.D ~ offset(log(.Y)) + (1 | provn) + (1 | bcohort) + agri + - Shige On 8/19/05, Shige Song [EMAIL PROTECTED] wrote: Dear Professor Bates, Here is output R 2.1.1 produced with control = list(EMverbose = TRUE, msVerbose = TRUE). I am getting the new devel version and see what will hapen there: EM iterations 0 85289.766 ( 5407.13: 0.0815) ( 26134.4: 0.00387) 1 84544.322 ( 333.732: 0.137) ( 1462.32: 0.00934) 2 84515.108 ( 129.506: 0.0270) ( 446.306: 0.00481) 3 84514.519 ( 115.592: 0.00355) ( 328.637: 0.00103) 4 84514.505 ( 113.981:0.000505) ( 311.160:0.000165) 5 84514.505 ( 113.755:7.45e-005) ( 308.524:2.50e-005) 6 84514.505 ( 113.722:1.11e-005) ( 308.128:3.77e-006) 7 84514.505 ( 113.717:1.66e-006) ( 308.068:5.66e-007) 8 84514.505 ( 113.716:2.50e-007) ( 308.059:8.49e-008) 9 84514.505 ( 113.716:3.74e-008) ( 308.058:1.27e-008) 10 84514.505 ( 113.716:5.62e-009) ( 308.058:1.91e-009) 11 84514.505 ( 113.716:8.43e-010) ( 308.058:2.87e-010) 12 84514.505 ( 113.716:1.27e-010) ( 308.058:4.31e-011) 13 84514.505 ( 113.716:1.90e-011) ( 308.057:6.47e-012) 14 84514.505 ( 113.716:2.86e-012) ( 308.057:9.73e-013) 15 84514.505 ( 113.716:4.25e-013) ( 308.057:1.44e-013) iter0 value 84514.505044 final value 84514.505044 converged EM iterations 0 83740.342 ( 113.716: -0.0164) ( 308.057:0.000596) 1 83740.273 ( 121.512:-0.00121) ( 298.914:-3.24e-005) 2 83740.272 ( 122.131:-0.000111) ( 299.397:-1.24e-005) iter0 value 83740.272232 final value 83740.272232 converged EM iterations 0 84011.550 ( 122.204:-0.00461) ( 299.576:-0.000256) 1 84011.543 ( 124.624:-0.000459) ( 303.453:-6.41e-005) 2 84011.543 ( 124.870:-5.42e-005) ( 304.440:-1.13e-005) iter0 value 84011.543350 final value 84011.543350 converged EM iterations 0 84018.592 ( 124.915:-6.44e-005) ( 304.548:-1.22e-005) 1 84018.592 ( 124.949:-8.29e-006) ( 304.737:-1.99e-006) 2 84018.592 ( 124.954:-1.15e-006) ( 304.768:-3.08e-007) iter0 value 84018.591624 final value 84018.591624 converged EM iterations 0 84018.612 ( 124.955:3.40e-007) ( 304.770:-1.98e-007) 1 84018.612 ( 124.955:-9.98e-009) ( 304.773:-2.23e-008) 2 84018.612 ( 124.955:-5.47e-009) (
[R] ?
Hello, I have a small problem with developing design matrix X, which I use in estimation the log-likelihood of a multinomial logit model. I have the data: number of observation - 289 number of choice alternative- 3 number of choice specific variables in matrix X -4 matrix X =289x4 I tried to use the function createX, I know that I have to get design matrix 289x12 (am I right?) but it always says bad dim (my code and data in attachment) Where is my fault? Can I use another method in order to create design matrix? Or need it at all here in logmnl (see code in attachment)? Can anyone help me with this issue? Thanks in advance, Tatyana df=read.table(data.dat,header=TRUE) inp=as.matrix(df) Y X1 X2 X3 X4 1 1 1 1 65 20999.89 2 1 1 2 67 2719.60 3 1 1 3 110 3581.09 4 1 1 4 64 1731.63 5 1 1 5 84 4434.97 6 1 6 1 90 691.32 7 1 6 2 31 228.50 8 1 6 3 33 615.12 9 1 6 4 39 910.62 10 1 7 1 169 1246.75 11 1 7 2 183 1183.03 12 1 7 3 203 1345.32 13 1 7 4 177 1088.98 14 1 7 5 169 896.42 15 1 8 1 71 1264.57 16 1 8 2 80 1094.40 17 1 8 3 75 1715.99 18 1 8 4 55 905.37 19 1 8 5 67 1448.17 20 1 10 1 349 1396.77 21 1 10 2 666 2026.89 22 1 10 3 480 774.37 23 1 10 4 456 1972.15 24 1 11 1 500 245.88 25 1 11 2 288 2927.77 26 1 11 3 211 9221.67 27 1 11 4 206 5632.91 28 1 11 5 175 1636.62 29 1 12 1 107 857.06 30 1 12 2 87 789.25 31 1 12 3 103 856.27 32 1 12 4 377 933.74 33 1 12 5 229 1316.31 34 1 13 1 32 149.13 35 1 13 2 19 153.74 36 1 13 3 25 179.60 37 1 13 4 28 252.70 38 1 13 5 22 294.80 39 1 14 1 47 1261.82 40 1 14 2 19 2332.21 41 1 15 1 348 558.91 42 1 15 2 399 550.91 43 1 15 3 388 797.68 44 1 15 4 208 804.76 45 1 15 5 241 673.12 46 1 17 1 70 151.06 47 1 17 2 96 255.22 48 1 17 3 102 1220.30 49 1 17 4 128 793.54 50 1 18 3 10 134.95 51 1 18 4 28 992.30 52 1 21 1 85 1170.71 53 1 21 2 257 464.95 54 1 21 3 353 404.21 55 1 21 4 293 517.64 56 1 21 5 515 1202.68 57 1 22 1 66 372.89 58 1 22 2 79 498.70 59 1 22 3 47 304.83 60 1 22 4 48 430.03 61 1 22 5 52 319.86 62 1 23 1 14 165.28 63 1 23 2 35 2044.52 64 1 23 3 20 499.59 65 1 24 1 94 107.76 66 1 24 2 5961.64 67 1 24 3 47 111.15 68 1 24 4 32 100.75 69 1 25 1 17 142.34 70 1 26 1 144 1105.71 71 1 26 2 196 1445.43 72 1 26 3 328 2297.11 73 1 26 4 517 2143.55 74 1 27 1 85 2457.58 75 1 27 2 99 1921.27 76 1 27 3 65 3380.86 77 1 27 4 88 2218.37 78 1 27 5 100 1881.00 79 1 29 1 107 561.27 80 1 29 2 67 557.43 81 1 29 3 49 387.71 82 1 30 1 77 106.50 83 1 30 2 225 267.87 84 1 30 3 520 502.18 85 1 30 4 552 443.07 86 1 30 5 319 255.50 87 1 31 1 38 6522.32 88 1 31 2 38 632.35 89 1 31 3 50 1615.18 90 1 31 4 53 1657.59 91 1 31 5 25 425.01 92 1 32 1 82 681.77 93 1 32 2 82 605.14 94 1 32 3 117 1068.86 95 1 32 4 90 638.95 96 1 33 1 53 350.89 97 1 33 2 39 378.53 98 1 33 3 44 432.31 99 1 34 1 61 752.13 100 1 34 2 76 1045.36 101 1 34 3 107 1344.42 102 1 34 4 65 1150.82 103 1 34 5 96 973.69 104 1 35 1 132 374.06 105 1 35 2 124 444.83 106 1 35 3 92 142.01 107 1 35 4 69 297.77 108 1 35 5 62 248.21 109 1 36 1 434 374.83 110 1 36 2 183 416.23 111 1 36 3 386 246.27 112 1 36 4 577 527.44 113 1 36 5 457 250.67 114 1 37 1 118 2306.72 115 1 37 2 169 1303.34 116 1 37 3 135 1741.13 117 1 37 4 103 1073.17 118 1 37 5 75 1146.11 119 1 40 1 66 1447.20 120 1 40 2 97 1352.28 121 1 40 3 65 1786.57 122 1 40 4 67 1060.59 123 1 42 1 26 241.23 124 1 42 2 43 334.35 125 1 42 3 65 381.51 126 1 42 4933.14 127 1 43 1 39 1504.44 128 1 43 2 33 1144.56 129 1 43 3 43 870.53 130 1 43 4 43 969.19 131 1 43 5 64 1655.93 132 1 44 12 1555.55 133 1 45 1 2284.39 134 1 46 1 46 996.07 135 1 46 2 33 777.97 136 1 46 3 60 637.64 137 1 46 4 42 1178.10 138 1 46 5 41 1054.84 139 1 47 1 37 1514.12 140 1 47 2 57 2132.21 141 1 47 3 53 2486.14 142 1 47 4 45 1807.57 143 1 47 5 45 1125.80 144 1 48 1 90 449.87 145 1 48 2 1286.38 146 1 48 3 44 159.58 147 1 48 4 42 372.35 148 1 48 5 58 442.60 149 1 49 1 92 645.82 150 1 49 2 82 523.96 151 1 49 3 132 833.91 152 1 49 4 125 490.37 153 1 49 5 89 454.82 154 1 50 1 30 105.94 155 1 50 2 2939.18 156 1 50 3 8016.13 157 1 50 4 185 106.54 158 1 51 1 95 937.76 159 1 51 2 34 1212.81 160 1 51 3 42 1254.46 161 1 51 4 35 644.77 162 1 51 5 36 426.90 163 1 52 1 42 138.73
Re: [R] Error messages using LMER
One thing to be noted: after switching to R-devel, even the simplest model can not converge. I always get this: - Warning messages: 1: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 2: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 3: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 4: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 5: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 6: nlminb failed to converge in: lmer(.D ~ offset(log(.Y)) + (1 | provn) + (1 | bcohort) + educy + - The same model did not have problems converging in R 2.1.1. Shige On 8/19/05, Shige Song [EMAIL PROTECTED] wrote: Here is what happened using R-devel: - EM iterations 0 85289.766 ( 5407.13: 0.0815) ( 26134.4: 0.00387) 1 84544.322 ( 333.732: 0.137) ( 1462.32: 0.00934) 2 84515.108 ( 129.506: 0.0270) ( 446.306: 0.00481) 3 84514.519 ( 115.592: 0.00355) ( 328.637: 0.00103) 4 84514.505 ( 113.981:0.000505) ( 311.160:0.000165) 5 84514.505 ( 113.755:7.45e-005) ( 308.524:2.50e-005) 6 84514.505 ( 113.722:1.11e-005) ( 308.128:3.77e-006) 7 84514.505 ( 113.717:1.66e-006) ( 308.068:5.66e-007) 8 84514.505 ( 113.716:2.50e-007) ( 308.059:8.49e-008) 9 84514.505 ( 113.716:3.74e-008) ( 308.058:1.27e-008) 10 84514.505 ( 113.716:5.62e-009) ( 308.058:1.91e-009) 11 84514.505 ( 113.716:8.43e-010) ( 308.058:2.87e-010) 12 84514.505 ( 113.716:1.27e-010) ( 308.058:4.31e-011) 13 84514.505 ( 113.716:1.90e-011) ( 308.057:6.47e-012) 14 84514.505 ( 113.716:2.86e-012) ( 308.057:9.73e-013) 15 84514.505 ( 113.716:4.25e-013) ( 308.057:1.44e-013) EM iterations 0 83740.342 ( 113.716: -0.0164) ( 308.057:0.000596) 1 83740.273 ( 121.512:-0.00121) ( 298.914:-3.24e-005) 2 83740.272 ( 122.131:-0.000111) ( 299.397:-1.24e-005) EM iterations 0 84011.546 ( 122.131:-0.00474) ( 299.397:-0.000265) 1 84011.539 ( 124.616:-0.000472) ( 303.415:-6.62e-005) 2 84011.539 ( 124.869:-5.58e-005) ( 304.433:-1.16e-005) EM iterations 0 84018.589 ( 124.869:-0.000139) ( 304.433:-1.81e-005) 1 84018.589 ( 124.944:-1.62e-005) ( 304.713:-3.26e-006) 2 84018.589 ( 124.953:-2.12e-006) ( 304.764:-5.23e-007) EM iterations 0 84018.611 ( 124.953:-2.38e-006) ( 304.764:-5.44e-007) 1 84018.611 ( 124.954:-3.25e-007) ( 304.772:-8.50e-008) 2 84018.611 ( 124.955:-4.66e-008) ( 304.773:-1.29e-008) EM iterations 0 84018.611 ( 124.955:-4.75e-008) ( 304.773:-1.30e-008) 1 84018.611 ( 124.955:-6.93e-009) ( 304.774:-1.97e-009) 2 84018.611 ( 124.955:-1.03e-009) ( 304.774:-2.96e-010) Warning messages: 1: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 2: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 3: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 4: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 5: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 6: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 7: nlminb failed to converge in: lmer(.D ~ offset(log(.Y)) + (1 | provn) + (1 | bcohort) + agri + - Shige On 8/19/05, Shige Song [EMAIL PROTECTED] wrote: Dear Professor Bates, Here is output R 2.1.1 produced with control = list(EMverbose = TRUE, msVerbose = TRUE). I am getting the new devel version and see what will hapen there: EM iterations 0 85289.766 ( 5407.13: 0.0815) ( 26134.4: 0.00387) 1 84544.322 ( 333.732: 0.137) ( 1462.32: 0.00934) 2 84515.108 ( 129.506: 0.0270) ( 446.306: 0.00481) 3 84514.519 ( 115.592: 0.00355) ( 328.637: 0.00103) 4 84514.505 ( 113.981:0.000505) ( 311.160:0.000165) 5 84514.505 ( 113.755:7.45e-005) ( 308.524:2.50e-005) 6 84514.505 ( 113.722:1.11e-005) ( 308.128:3.77e-006) 7 84514.505 ( 113.717:1.66e-006) ( 308.068:5.66e-007) 8 84514.505 ( 113.716:2.50e-007) ( 308.059:8.49e-008) 9 84514.505 ( 113.716:3.74e-008) ( 308.058:1.27e-008) 10 84514.505 ( 113.716:5.62e-009) ( 308.058:1.91e-009) 11 84514.505 ( 113.716:8.43e-010) ( 308.058:2.87e-010) 12 84514.505 ( 113.716:1.27e-010) ( 308.058:4.31e-011) 13 84514.505 ( 113.716:1.90e-011) ( 308.057:6.47e-012) 14 84514.505 (
Re: [R] Error messages using LMER
Thanks for including all of that output. I believe that in this version the parameters are the relative variances. This would indicate that somehow you are getting fits with very low residual sums of squares in the weighted least squares problem. It could be that you have too many fixed effects terms in the model and are getting complete separation. On 8/18/05, Shige Song [EMAIL PROTECTED] wrote: Here is what happened using R-devel: - EM iterations 0 85289.766 ( 5407.13: 0.0815) ( 26134.4: 0.00387) 1 84544.322 ( 333.732: 0.137) ( 1462.32: 0.00934) 2 84515.108 ( 129.506: 0.0270) ( 446.306: 0.00481) 3 84514.519 ( 115.592: 0.00355) ( 328.637: 0.00103) 4 84514.505 ( 113.981:0.000505) ( 311.160:0.000165) 5 84514.505 ( 113.755:7.45e-005) ( 308.524:2.50e-005) 6 84514.505 ( 113.722:1.11e-005) ( 308.128:3.77e-006) 7 84514.505 ( 113.717:1.66e-006) ( 308.068:5.66e-007) 8 84514.505 ( 113.716:2.50e-007) ( 308.059:8.49e-008) 9 84514.505 ( 113.716:3.74e-008) ( 308.058:1.27e-008) 10 84514.505 ( 113.716:5.62e-009) ( 308.058:1.91e-009) 11 84514.505 ( 113.716:8.43e-010) ( 308.058:2.87e-010) 12 84514.505 ( 113.716:1.27e-010) ( 308.058:4.31e-011) 13 84514.505 ( 113.716:1.90e-011) ( 308.057:6.47e-012) 14 84514.505 ( 113.716:2.86e-012) ( 308.057:9.73e-013) 15 84514.505 ( 113.716:4.25e-013) ( 308.057:1.44e-013) EM iterations 0 83740.342 ( 113.716: -0.0164) ( 308.057:0.000596) 1 83740.273 ( 121.512:-0.00121) ( 298.914:-3.24e-005) 2 83740.272 ( 122.131:-0.000111) ( 299.397:-1.24e-005) EM iterations 0 84011.546 ( 122.131:-0.00474) ( 299.397:-0.000265) 1 84011.539 ( 124.616:-0.000472) ( 303.415:-6.62e-005) 2 84011.539 ( 124.869:-5.58e-005) ( 304.433:-1.16e-005) EM iterations 0 84018.589 ( 124.869:-0.000139) ( 304.433:-1.81e-005) 1 84018.589 ( 124.944:-1.62e-005) ( 304.713:-3.26e-006) 2 84018.589 ( 124.953:-2.12e-006) ( 304.764:-5.23e-007) EM iterations 0 84018.611 ( 124.953:-2.38e-006) ( 304.764:-5.44e-007) 1 84018.611 ( 124.954:-3.25e-007) ( 304.772:-8.50e-008) 2 84018.611 ( 124.955:-4.66e-008) ( 304.773:-1.29e-008) EM iterations 0 84018.611 ( 124.955:-4.75e-008) ( 304.773:-1.30e-008) 1 84018.611 ( 124.955:-6.93e-009) ( 304.774:-1.97e-009) 2 84018.611 ( 124.955:-1.03e-009) ( 304.774:-2.96e-010) Warning messages: 1: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 2: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 3: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 4: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 5: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 6: optim or nlminb returned message See PORT documentation. Code (27) in: LMEopt(x = mer, value = cv) 7: nlminb failed to converge in: lmer(.D ~ offset(log(.Y)) + (1 | provn) + (1 | bcohort) + agri + - Shige On 8/19/05, Shige Song [EMAIL PROTECTED] wrote: Dear Professor Bates, Here is output R 2.1.1 produced with control = list(EMverbose = TRUE, msVerbose = TRUE). I am getting the new devel version and see what will hapen there: EM iterations 0 85289.766 ( 5407.13: 0.0815) ( 26134.4: 0.00387) 1 84544.322 ( 333.732: 0.137) ( 1462.32: 0.00934) 2 84515.108 ( 129.506: 0.0270) ( 446.306: 0.00481) 3 84514.519 ( 115.592: 0.00355) ( 328.637: 0.00103) 4 84514.505 ( 113.981:0.000505) ( 311.160:0.000165) 5 84514.505 ( 113.755:7.45e-005) ( 308.524:2.50e-005) 6 84514.505 ( 113.722:1.11e-005) ( 308.128:3.77e-006) 7 84514.505 ( 113.717:1.66e-006) ( 308.068:5.66e-007) 8 84514.505 ( 113.716:2.50e-007) ( 308.059:8.49e-008) 9 84514.505 ( 113.716:3.74e-008) ( 308.058:1.27e-008) 10 84514.505 ( 113.716:5.62e-009) ( 308.058:1.91e-009) 11 84514.505 ( 113.716:8.43e-010) ( 308.058:2.87e-010) 12 84514.505 ( 113.716:1.27e-010) ( 308.058:4.31e-011) 13 84514.505 ( 113.716:1.90e-011) ( 308.057:6.47e-012) 14 84514.505 ( 113.716:2.86e-012) ( 308.057:9.73e-013) 15 84514.505 ( 113.716:4.25e-013) ( 308.057:1.44e-013) iter0 value 84514.505044 final value 84514.505044 converged EM iterations 0 83740.342 ( 113.716: -0.0164) ( 308.057:0.000596) 1 83740.273 ( 121.512:-0.00121) ( 298.914:-3.24e-005) 2 83740.272 ( 122.131:-0.000111) ( 299.397:-1.24e-005) iter0 value 83740.272232 final value 83740.272232 converged EM iterations 0 84011.550 ( 122.204:-0.00461) ( 299.576:-0.000256) 1 84011.543 ( 124.624:-0.000459) ( 303.453:-6.41e-005) 2 84011.543 ( 124.870:-5.42e-005) ( 304.440:-1.13e-005) iter0 value
Re: [R] Set R 2.1.1. in English
On Thu, 18 Aug 2005, Patrick Giraudoux wrote: (whish R 2.1.1 could be parametrise 'English' even with a French Windows XP) If I understand you correctly, it can. Just add LANGUAGE=en to the shortcut. Wonderful hope but not sure to catch what you term shortcut. I tried to add this command in C:\R\rw2011\etc\Rprofile, the .Rprofile in the folder my documents, but this cannot be understood from there... which obviously shows that shortcut is not a general term for profiles! I have also See the rw-FAQ Q2.2. `Shortcut' is the standard name for files on Windows with extension .lnk. Another way is to add this to HOME/.Renviron: see the rw-FAQ Q2.17 unsuccessfully looked for a file shortcut in the rw2011 folder. I also tried and went through the R-help archive. There were some exchanges on this subject and Asian languages some weeks ago but what I have tried and adapt on this basis did not work. The objective would be to have R in English (thus additionnally allowing MDI mode with R-WinEdt, which is not the case with any other language) and to keep Windows XP and other applications in the foreign (= French, here) language. Thanks for any further hint, Patrick Giraudoux -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ 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
Re: [R] standard errors for expression intensities
Please ask Bioconductor Qs on their mailing list! PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html This is indeed covered there, as is the (mis)use of HTML mail. On Thu, 18 Aug 2005, Devarajan, Karthik wrote: Hello! When I use the functions rma() or justrma() in the Bioconductor affy package, the standard errors for the expression estimates computed by the function se.exprs() is a matrix of NA's. I am wondering if any of you have encountered this problem and if there is a solution. Any help would be appreciated. Thanks, Karthik. [[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 -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ 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] problem with repeated formal arguments and ...
I want to add an argument if it is not present. Following the Green
Book, p. 337:
test - function(x, ...){
dots - list(...)
if (!hasArg(from))
from - 0
else
from - dots$from
curve(x, from=from, ...)
}
test(sin)
test(sin, from=4)
Error in curve(x, from = from, ...) : formal argument from matched by
multiple actual arguments
The FAQ says, in the section on differences between R and S,
R disallows repeated formal arguments in function calls.
That seems a perfectly reasonable rule, but how do I handle this
situation?
--
Ross Boylan wk: (415) 502-4031
530 Parnassus Avenue (Library) rm 115-4 [EMAIL PROTECTED]
Dept of Epidemiology and Biostatistics fax: (415) 476-9856
University of California, San Francisco
San Francisco, CA 94143-0840 hm: (415) 550-1062
__
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
Re: [R] How do I make a Sweave + latex table out of this ?
Hi Fredrik, You can do this nicely with latex() if you massage this a little so that you end up with a data frame that includes a column for agemF (to indicate which agemF each row belongs to) and then use the rgroup and n.rgroup options of latex(). Dave Fredrik Karlsson wrote: Dear list, I have a table that I would like to convert to latex for inclusion into a Sweave file. round(ftable(prop.table(xtabs(~agemF + votcat + Type , data=work),margin=2))*100,1) Type Voiced Voiceless unaspirated Voiceless aspirated agemF votcat 18 - 24 Prevoiced 2.6 8.7 2.3 Short lag 5.8 6.7 5.1 Long lag 1.0 1.9 2.9 24 - 30 Prevoiced15.1 10.5 1.7 Short lag 9.2 15.3 5.8 Long lag 3.5 8.115.8 30 - 36 Prevoiced12.8 14.0 2.6 Short lag10.2 14.2 3.0 Long lag 2.3 5.522.2 36 - 42 Prevoiced 4.4 6.4 0.6 Short lag 4.0 5.9 1.5 Long lag 1.3 2.9 9.4 42 - 48 Prevoiced 6.4 2.3 0.3 Short lag 3.0 2.8 1.4 Long lag 0.6 7.7 8.8 48 - 54 Prevoiced 4.9 4.1 0.3 Short lag 2.0 2.7 1.3 Long lag 0.3 0.9 4.7 However, I have not been able to use this as a table. The Hmisc latex command only accepts the input if I first convert it to a data.frame format, and that makes the output much more difficult to read as it duplicates the category levels of agemF. Is there a way to do this? /Fredrik Karlsson __ 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 -- David Whiting School of Clinical Medical Sciences, The Medical School University of Newcastle upon Tyne, NE2 4HH, UK. I love deadlines. I love the whooshing noise they make as they go by (Douglas Adams) __ 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
Re: [R] Set R 2.1.1. in English
OK. ¨Perfectly clear. This written for other listers: Shortcut is the shortcuts (alias) everybody use to launch applications (in French: raccourci), and for suckers (as I am): one just have to right-click the shortcut used to launch Rgui.exe, select properties (propriétés in French) and add in the box target (cible in French) LANGUAGE=en to the path. In my case: C:\R\rw2011\bin\Rgui.exe. Thus, it makes C:\R\rw2011\bin\Rgui.exe LANGUAGE=en. Then click OK. It works fantastic Thanks a lot for this. It makes three weeks I was grumling everytime (often) I started R, and frustrated with the SDI mode... Patrick Giraudoux Prof Brian Ripley a écrit : On Thu, 18 Aug 2005, Patrick Giraudoux wrote: (whish R 2.1.1 could be parametrise 'English' even with a French Windows XP) If I understand you correctly, it can. Just add LANGUAGE=en to the shortcut. Wonderful hope but not sure to catch what you term shortcut. I tried to add this command in C:\R\rw2011\etc\Rprofile, the .Rprofile in the folder my documents, but this cannot be understood from there... which obviously shows that shortcut is not a general term for profiles! I have also See the rw-FAQ Q2.2. `Shortcut' is the standard name for files on Windows with extension .lnk. Another way is to add this to HOME/.Renviron: see the rw-FAQ Q2.17 unsuccessfully looked for a file shortcut in the rw2011 folder. I also tried and went through the R-help archive. There were some exchanges on this subject and Asian languages some weeks ago but what I have tried and adapt on this basis did not work. The objective would be to have R in English (thus additionnally allowing MDI mode with R-WinEdt, which is not the case with any other language) and to keep Windows XP and other applications in the foreign (= French, here) language. Thanks for any further hint, Patrick Giraudoux __ 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
Re: [R] problem with repeated formal arguments and ...
Ross Boylan wrote:
I want to add an argument if it is not present. Following the Green
Book, p. 337:
test - function(x, ...){
dots - list(...)
if (!hasArg(from))
from - 0
else
from - dots$from
curve(x, from=from, ...)
}
test(sin)
test(sin, from=4)
Error in curve(x, from = from, ...) : formal argument from matched by
multiple actual arguments
The FAQ says, in the section on differences between R and S,
R disallows repeated formal arguments in function calls.
That seems a perfectly reasonable rule, but how do I handle this
situation?
Hi, Ross,
Add from to your function:
test - function(x, from = 0, ...) {
curve(x, from = from, ...)
}
Or another way:
test - function(x, ...) {
dots - list(...)
if(!hasArg(from)) dots$from - 0
dots$expr - x
do.call(curve, dots)
}
I actually prefer the latter if I'm changing many arguments. I do this
quite often when writing custom lattice plots and I want to override
many of the defaults.
HTH,
--sundar
__
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
Re: [R] trouble with wilcox.test
Ok, I will think more about the appropriateness of the Wilcoxon test here. I was using R version 2.1.1, 2005-06-20 Windows XP SP2 512MB RAM --Greg - Original Message - From: Prof Brian Ripley [EMAIL PROTECTED] To: Greg Hather [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Sent: Wednesday, August 17, 2005 11:45 PM Subject: Re: [R] trouble with wilcox.test On Wed, 17 Aug 2005, Greg Hather wrote: I'm having trouble with the wilcox.test command in R. Are you sure it is not the concepts that are giving 'trouble'? What real problem are you trying to solve here? To demonstrate the anomalous behavior of wilcox.test, consider wilcox.test(c(1.5,5.5), c(1:1), exact = F)$p.value [1] 0.01438390 wilcox.test(c(1.5,5.5), c(1:1), exact = T)$p.value [1] 6.39808e-07 (this calculation takes noticeably longer). wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value (R closes/crashes) I believe that wilcox.test(c(1.5,5.5), c(1:1), exact = F)$p.value yields a bad result because of the normal approximation which R uses when exact = F. Expecting an approximation to be good in the tail for m=2 is pretty unrealistic. But then so is believing the null hypothesis of a common *continuous* distribution. Why worry about the distribution under a hypothesis that is patently false? People often refer to this class of tests as `distribution-free', but they are not. The Wilcoxon test is designed for power against shift alternatives, but here there appears to be a very large difference in spread. So wilcox.test(5000+c(1.5,5.5), c(1:1), exact = T)$p.value [1] 0.9989005 even though the two samples differ in important ways. Any suggestions for how to compute wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value? I get (current R 2.1.1 on Linux) wilcox.test(c(1.5,5.5), c(1:2), exact = T)$p.value [1] 1.59976e-07 and no crash. So the suggestion is to use a machine adequate to the task, and that probably means an OS with adequate stack size. [[alternative HTML version deleted]] PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Please do heed it. What version of R and what machine is this? And do take note of the request about HTML mail. -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ 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
Re: [R] Regular expressions sub
Dirk Eddelbuettel [EMAIL PROTECTED] writes: Bernd Weiss bernd.weiss at uni-koeln.de writes: I am struggling with the use of regular expression. I got as.character(test$sample.id) [1] 1.11 10.11 11.11 113.31 114.2 114.3 114.8 and need [1] 11 11 11 31 2 3 8 I.e. remove everything before the . . Define the dot as the hard separator, and allow for multiple digits before it: sample.id - c(1.11, 10.11, 11.11, 113.31, 114.2, 114.3, 114.8) gsub(^[0-9]*\., , sample.id) [1] 11 11 11 31 2 3 8 Or, more longwinded, but with less assumptions about what goes before the dot: gsub(^.*\\.(.*)$,\\1,sample.id) [1] 11 11 11 31 2 3 8 or, gsub(^.*\\.([^.]*)$,\\1,sample.id) [1] 11 11 11 31 2 3 8 -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ 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] Console
I am at my first steps with R... and I already notice that the console has a quite limited number of lines. Can anyone tell me how to visualise all the information, which is actually present? I only see the last part of the output, which obviosly exceeds the maximum number of rows in the console. Thank you very much for your help! Daniela [[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
[R] How to put factor variables in an nls formula ?
Hello, I want to fit a Gompertz model for tree diameter growth that depends on a 4 levels edaphic factor (Drain) and I dont manage to introduce the factor variable in the formula. Dinc is the annual diameter increment and D is the Diameter. treestab Dinc D Drain [1,] 0.03 26.10 2 [2,] 0.04 13.05 1 [3,] 0.00 24.83 1 [4,] 0.00 15.92 4 [5,] 0.00 12.25 4 [6,] 0.00 11.78 4 [7,] 0.00 16.87 4 [8,] 0.00 15.12 4 [9,] -0.01 13.53 4 [10,] 0.04 16.55 3 [11,] 0.025 16.07 3 [12,] 0.00 30.24 3 [13,] 0.06 15.28 2 etc contrasts(Drain)-contr.sum(4) mymodel-nls(Dinc~r*(1+Drain)*D*log(Asym/D), data=treestab, start(r=0.05, Asym=40)) Error in numericDeriv(form[[3]], names(ind), env) : Missing value or an infinity produced when evaluating the model In addition: Warning messages: 1: + not meaningful for factors in: Ops.factor(1, Drain) 2: + not meaningful for factors in: Ops.factor(1, Drain) Do I need to use another function instead of nls to correctly include the factor Drain ? Thanks, François __ 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
Re: [R] problem with repeated formal arguments and ...
Sundar Dorai-Raj wrote:
Ross Boylan wrote:
I want to add an argument if it is not present. Following the Green
Book, p. 337:
test - function(x, ...){
dots - list(...)
if (!hasArg(from))
from - 0
else
from - dots$from
curve(x, from=from, ...)
}
test(sin)
test(sin, from=4)
Error in curve(x, from = from, ...) : formal argument from matched by
multiple actual arguments
The FAQ says, in the section on differences between R and S,
R disallows repeated formal arguments in function calls.
That seems a perfectly reasonable rule, but how do I handle this
situation?
Hi, Ross,
Add from to your function:
test - function(x, from = 0, ...) {
curve(x, from = from, ...)
}
Or another way:
test - function(x, ...) {
dots - list(...)
if(!hasArg(from)) dots$from - 0
dots$expr - x
do.call(curve, dots)
}
I didn't check this example, but for curve, since `x' is an expression,
we should actually do the following:
test - function(x, ...) {
dots - list(...)
if(!hasArg(from)) dots$from - 0
dots$expr - substitute(x)
invisible(do.call(curve, dots))
}
test(x^3 - 3 * x, to = 5)
test(x^3 - 3 * x, from = -5, to = 5)
HTH,
--sundar
I actually prefer the latter if I'm changing many arguments. I do this
quite often when writing custom lattice plots and I want to override
many of the defaults.
HTH,
--sundar
__
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-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] Set R 2.1.1. in English
Another way is to create a file called Renviron.site with the following single line: LANGUAGE=en and put it in your ...\R\rw2011\etc folder. In your case, echo LANGUAGE=en c:\R\rw2011\etc\Renviron.site would do it. This would result in R being in English not only if you (1) used the shortcut but also if you (2) used the command line to start R or if you (3) used Start Programs R , rw2011 . On 8/18/05, Patrick Giraudoux [EMAIL PROTECTED] wrote: OK. ¨Perfectly clear. This written for other listers: Shortcut is the shortcuts (alias) everybody use to launch applications (in French: raccourci), and for suckers (as I am): one just have to right-click the shortcut used to launch Rgui.exe, select properties (propriétés in French) and add in the box target (cible in French) LANGUAGE=en to the path. In my case: C:\R\rw2011\bin\Rgui.exe. Thus, it makes C:\R\rw2011\bin\Rgui.exe LANGUAGE=en. Then click OK. It works fantastic Thanks a lot for this. It makes three weeks I was grumling everytime (often) I started R, and frustrated with the SDI mode... Patrick Giraudoux Prof Brian Ripley a écrit : On Thu, 18 Aug 2005, Patrick Giraudoux wrote: (whish R 2.1.1 could be parametrise 'English' even with a French Windows XP) If I understand you correctly, it can. Just add LANGUAGE=en to the shortcut. Wonderful hope but not sure to catch what you term shortcut. I tried to add this command in C:\R\rw2011\etc\Rprofile, the .Rprofile in the folder my documents, but this cannot be understood from there... which obviously shows that shortcut is not a general term for profiles! I have also See the rw-FAQ Q2.2. `Shortcut' is the standard name for files on Windows with extension .lnk. Another way is to add this to HOME/.Renviron: see the rw-FAQ Q2.17 unsuccessfully looked for a file shortcut in the rw2011 folder. I also tried and went through the R-help archive. There were some exchanges on this subject and Asian languages some weeks ago but what I have tried and adapt on this basis did not work. The objective would be to have R in English (thus additionnally allowing MDI mode with R-WinEdt, which is not the case with any other language) and to keep Windows XP and other applications in the foreign (= French, here) language. Thanks for any further hint, Patrick Giraudoux __ 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-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] A. Mani : Avoiding loops
Hello, I want to avoid loops in the following situation. There is a 5-col dataframe with col headers alone. two of the columns are non-numeric. The problem is to calculate statistics(scores) for each element of one column. The functions depend on matching in the other non-numeric column. A B C E F 1 2 X Y 1 2 3 G L 1 3 1 G L 5 and so on ...3+ entries. I need scores for col E entries which depend on conditional implications. Thanks, A. Mani Member, Cal. Math. Soc __ 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] Use of contains in S4 classes
setClass(B, representation=representation(B, extra=numeric)) setClass(B, representation=representation(extra=numeric), contains=B) Are these the same? If not, how do they differ? What about setClass(B, representation=representation(B, extra=numeric), contains=B) ? As far as I can tell, the Green Book doesn't talk about a contains argument to setClass. -- Ross Boylan wk: (415) 502-4031 530 Parnassus Avenue (Library) rm 115-4 [EMAIL PROTECTED] Dept of Epidemiology and Biostatistics fax: (415) 476-9856 University of California, San Francisco San Francisco, CA 94143-0840 hm: (415) 550-1062 __ 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] R equivalent to `estimate' in SAS proc mixed
Example: I have the following model model - lmer(response ~ time * trt * bio + (time|id), data = dat) where time = time of observation trt = treatment group (0-no treatment / 1-treated) bio = biological factor (0-absent / 1-present) and I would like to obtain an estimate (with standard error) of the change in response over time for individuals in the treatment group with the biological factor. The estimate is easy, sum(fixef(model)[c(2,5,6,8)]) # ie time + time:trt + time:bio + time:trt:bio but the standard error is a hassle to calculate by hand. Is there some better way to do this? In SAS for example there is an `estimate' option (see sample code below) that will calculate the estimate, SE, df, t statistic, etc... Is there some R equivalent? Thanks, Randy proc mixed data=dat; class id; model response = time + trt + bio + time*trt + time*bio + trt*bio + time*trt*bio; random time; estimate est1 intercept 0 time 1 trt 0 bio 0 time*trt 1 time*bio 1 trt*bio 0 time*trt*bio 1; /* or something like that */ run; ~~ Randy Johnson Laboratory of Genomic Diversity NCI-Frederick Bldg 560, Rm 11-85 Frederick, MD 21702 (301)846-1304 ~~ __ 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
Re: [R] Console
Quoting Daniela Salvini [EMAIL PROTECTED]: I am at my first steps with R... and I already notice that the console has a quite limited number of lines. Can anyone tell me how to visualise all the information, which is actually present? I only see the last part of the output, which obviosly exceeds the maximum number of rows in the console. Thank you very much for your help! Daniela visualize suggests plotting. do you mean how do I look at the whole dataset? you could print a few lines at a time, say X[1:25,]. With bigger datasets I usually look at the file in a text editor like emacs... albyn __ 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
Re: [R] Use of contains in S4 classes
On Thu, 18 Aug 2005, Ross Boylan wrote: setClass(B, representation=representation(B, extra=numeric)) setClass(B, representation=representation(extra=numeric), contains=B) Are these the same? If not, how do they differ? What about setClass(B, representation=representation(B, extra=numeric), contains=B) ? As far as I can tell, the Green Book doesn't talk about a contains argument to setClass. S4 - Composition and Inheritance by Witold Eryk Wolski (a.k.a. Extending.pdf) might be what you're looking for. -- SIGSIG -- signature too long (core dumped) __ 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
Re: [R] Use of contains in S4 classes
Oops, the second class should have been A in the examples. Corrected version: setClass(B, representation=representation(A, extra=numeric)) setClass(B, representation=representation(extra=numeric), contains=A) Are these the same? If not, how do they differ? What about setClass(B, representation=representation(A, extra=numeric), contains=A) ? As far as I can tell, the Green Book doesn't talk about a contains argument to setClass. -- Ross Boylan wk: (415) 502-4031 530 Parnassus Avenue (Library) rm 115-4 [EMAIL PROTECTED] Dept of Epidemiology and Biostatistics fax: (415) 476-9856 University of California, San Francisco San Francisco, CA 94143-0840 hm: (415) 550-1062 __ 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] 0/0, R segfaults
Hi, I noticed that when I was conducting some calculation involving finding correlation coeficients, R stopped abnormally. So I did some research, and find out that 0/0 was the culprit. For sure 0/0 is not a valid expression, but R should give a warning, an error msg or NaN instead of segmentation fault. I am using R 2.1.0 under Gentoo Linux. My GCC version is 3.3.5. Xing __ 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] axis label justified
Hi, I am trying to make my axis labels left justified,
and have used adj=0 in the axis() without success. Can
anyone have a suggestion?
axis(2,at=1:50,labels=paste('a',1:50,sep=''),las=2,cex.axis=0.5,adj=0,tck=0,mgp=c(3,0.5,0))
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
Re: [R] axis label justified
I note that the axis help seems to refer to padj. After playing around it is
obvious that I don't know what is meant by this argument, so maybe I'm doing
something wrong. My practical soultion is
plot(1:50,axes = FALSE,ylab = )
axis(2,at = 1:50,labels = rep(,50),las = 2,padj = 0)
text(rep(-4,5),1:50,labels=paste('a',1:50,sep=''),xpd = TRUE,adj = 0,cex=0.5)
Tom
-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] Behalf Of array chip
Sent: Friday, 19 August 2005 7:09 AM
To: r-help@stat.math.ethz.ch
Subject: [R] axis label justified
Hi, I am trying to make my axis labels left justified,
and have used adj=0 in the axis() without success. Can
anyone have a suggestion?
axis(2,at=1:50,labels=paste('a',1:50,sep=''),las=2,cex.axis=0.
5,adj=0,tck=0,mgp=c(3,0.5,0))
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
__
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
Re: [R] 0/0, R segfaults
On 18 August 2005 at 16:01, Xing Qiu wrote: | Hi, | | I noticed that when I was conducting some calculation involving | finding correlation coeficients, R stopped abnormally. So I did some | research, and find out that 0/0 was the culprit. For sure 0/0 is not | a valid expression, but R should give a warning, an error msg or NaN | instead of segmentation fault. | | I am using R 2.1.0 under Gentoo Linux. My GCC version is 3.3.5. [EMAIL PROTECTED]:~ R R : Copyright 2005, The R Foundation for Statistical Computing Version 2.1.1 (2005-06-20), ISBN 3-900051-07-0 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for a HTML browser interface to help. Type 'q()' to quit R. 0/0 [1] NaN No problem on Debian 'testing' with R 2.1.1. You may want to try a different libc. Dirk -- Statistics: The (futile) attempt to offer certainty about uncertainty. -- Roger Koenker, 'Dictionary of Received Ideas of Statistics' __ 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] Problem with get.hist.quote() in tseries
When using get.hist.quote(), I find the dates are broken. This is with R 2.1.1 on Mac OS X `panther'. library(tseries) Loading required package: quadprog 'tseries' version: 0.9-27 'tseries' is a package for time series analysis and computational finance. See 'library(help=tseries)' for details. x - get.hist.quote(^VIX) trying URL 'http://chart.yahoo.com/table.csv?s=^VIXa=0b=02c=1991d=7e=18f=2005g=dq=qy=0z=^VIXx=.csv' Content type 'text/csv' length unknown opened URL .. .. .. .. .. .. .. .. .. .. .. .. .. .. .. downloaded 150Kb head(x) [1] 26.62 27.93 27.19NANA 28.95 x[1:10,] Open High Low Close [1,] 26.62 26.62 26.62 26.62 [2,] 27.93 27.93 27.93 27.93 [3,] 27.19 27.19 27.19 27.19 [4,]NANANANA [5,]NANANANA [6,] 28.95 28.95 28.95 28.95 [7,] 30.38 30.38 30.38 30.38 [8,] 33.30 33.30 33.30 33.30 [9,] 31.33 31.33 31.33 31.33 [10,] 32.63 32.63 32.63 32.63 plot(x) (dates don't show). str(x) mts [1:5343, 1:4] 26.6 27.9 27.2 NA NA ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr [1:4] Open High Low Close - attr(*, tsp)= num [1:3] 33240 38582 1 - attr(*, class)= chr [1:2] mts ts I wonder what I'm doing wrong. -- Ajay Shah Consultant [EMAIL PROTECTED] Department of Economic Affairs http://www.mayin.org/ajayshah Ministry of Finance, New Delhi __ 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
Re: [R] Regular expressions sub
On 18 Aug 2005 at 21:17, Peter Dalgaard wrote: Dirk Eddelbuettel [EMAIL PROTECTED] writes: Bernd Weiss bernd.weiss at uni-koeln.de writes: I am struggling with the use of regular expression. I got as.character(test$sample.id) [1] 1.11 10.11 11.11 113.31 114.2 114.3 114.8 and need [1] 11 11 11 31 2 3 8 I.e. remove everything before the . . Define the dot as the hard separator, and allow for multiple digits before it: sample.id - c(1.11, 10.11, 11.11, 113.31, 114.2, 114.3, 114.8) gsub(^[0-9]*\., , sample.id) [1] 11 11 11 31 2 3 8 Or, more longwinded, but with less assumptions about what goes before the dot: gsub(^.*\\.(.*)$,\\1,sample.id) [1] 11 11 11 31 2 3 8 Wow, thanks a lot for all the valuable suggestions. Bernd __ 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
