Re: [R] is.category
Applejus gmail.com> writes: > > > Hello, > > Could someone tell me what the SPLUS "is.category" function do and what is > its equivalent in R? > Thank you, I couldn't find any help elsewhere... > > If you Web search for "splus help" or similar, you will find, among others, http://www.uni-muenster.de/ZIV/Mitarbeiter/BennoSueselbeck/s-html/shelp.html and there under "is.category" "Creates or tests for categorical objects (vectors with a "levels" attribute). This function is deprecated; factor or ordered are now preferred." Hans Werner __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Format integer
This is one of those problems where the fine details matter. 1) The version of R. I optimized sprintf() for long inputs and a single format in R 2.7.0 -- the differences are mainly for multiple inputs and where coercion is needed. See also below. 2) The system. My home system with an Intel Core 2 Duo is usually about the same speed as my office desktop with dual Opterons. But not here: Home: system.time(a<-formatC(x,digits=10,flag='0')) user system elapsed 9.705 0.088 9.810 system.time(b<-sprintf("%011d",x)) user system elapsed 0.283 0.000 0.283 Office: system.time(a<-formatC(x,digits=10,flag='0')) user system elapsed 15.851 0.125 16.007 system.time(b<-sprintf("%011d",x)) user system elapsed 0.816 0.001 0.818 and my Windows laptop is similar to the second here. So a speed-up of 95x seems atypical. On Mon, 12 May 2008, Phil Spector wrote: I guess "little" means different things to different people: x = sample(1:100,65,replace=TRUE) system.time(a<-formatC(x,digits=10,flag='0')) user system elapsed 32.854 0.444 34.813 system.time(b<-sprintf("%011d",x)) user system elapsed 0.352 0.012 0.363 If you look at the definitions of the functions, you'll see that formatC is written in R, and sprintf uses a single call to an .Internal function. I Not really: the meat of formatC() is a .C call. In this case it is calling format.default(), also a .Internal. But profiling shows that most of the time here is spent in paste(), another function which was optimized in 2.7.0. (I see 2.7.0 as 1.7x faster than 2.6.2 on formatC here.) But although sprintf is more flexible, on most problems it will be substantially faster. - Phil Spector Statistical Computing Facility Department of Statistics UC Berkeley [EMAIL PROTECTED] On Mon, 12 May 2008, Anh Tran wrote: Yea, thanks all. I checked back and I got a few things mistyped. The array is 650,000 and it took 25 seconds :p. It's acceptable. Just that I had too many variable at the time I ran it. Also, seems like sprintf is a little faster. Thanks all. Anh Tran On Mon, May 12, 2008 at 2:55 PM, Uwe Ligges <[EMAIL PROTECTED]> wrote: Anh Tran wrote: Thanks. formatC(flag) works. But it's awefully slow. I try to do that for 65000 numbers (generating ID for each item) and it seems like forever. On my not that recent laptop: system.time(formatC(1:65000, width=10, flag="0")) user system elapsed 1.920.001.94 I think 2 seconds is less than "forever". Uwe Ligges Is there any faster way? Thank all. Anh Tran On Mon, May 12, 2008 at 2:36 PM, Uwe Ligges < [EMAIL PROTECTED]> wrote: Anh Tran wrote: Hi, What's one way to convert an integer to a string with preceding 0's? such that '13' becomes '013' to be put into a string I've tried formatC, but they removes all the zeros and replace it with blanks Not so for me: formatC(13, digits=10, flag="0") Uwe LIgges Thanks -- Regards, Anh Tran [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- 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@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] R help: problems with step function
Dear List Members, I have encountered two problems when using the step function to select models. To better illustrate the problems, attached is an R image which includes the objects needed to run the code attached. lm.data.frame have factor variables with 3 levels. The following run shows the first problem. AICs (* and **) are different. I noticed that the Df for rs13482096:rs13483699 is 4, while I think Df should be 6, 2 from rs13483699 and 4 from interactions. When I ran add1 directly, I got Df=6 and AIC 848.75. > step2.bic.out <- step(step.bic.out, scope=list(lower=scope.lower2, upper=scope.upper2), + direction="both", k=log(length(step.bic.out$y)), trace=1) Start: AIC=841.18 pheno.dat ~ rs13479085 + rs13480057 + rs13482096 + rs8254221 Df DevianceAIC + rs13482096:rs13483699 4 216.63 840.63 (*) 233.82 841.18 - rs8254221 2 244.08 842.90 - rs13482096 2 245.20 844.31 .. Step: AIC=848.75 (**) pheno.dat ~ rs13479085 + rs13480057 + rs13482096 + rs8254221 + rs13482096:rs13483699 > add1(step.bic.out, scope="rs13482096:rs13483699", k=log(length(step.bic.out$y))) Single term additions Model: pheno.dat ~ rs13479085 + rs13480057 + rs13482096 + rs8254221 Df DevianceAIC 233.82 841.18 rs13482096:rs13483699 6 214.28 848.75 (**) When I used add1 to handle terms to be added together and separately, I got different results. The example is shown below. This might explain the discrepancy shown above. > add1(step.bic.out, scope=int.terms[11:12], k=log(length(step.bic.out$y))) Single term additions Model: pheno.dat ~ rs13479085 + rs13480057 + rs13482096 + rs8254221 Df DevianceAIC 233.82 841.18 rs13479085:rs13475933 6 224.95 863.66 rs13480057:rs13475933 4 226.72 854.62 (***) > add1(step.bic.out, scope=int.terms[11], k=log(length(step.bic.out$y))) Single term additions Model: pheno.dat ~ rs13479085 + rs13480057 + rs13482096 + rs8254221 Df DevianceAIC 233.82 841.18 rs13479085:rs13475933 6 224.95 863.66 > add1(step.bic.out, scope=int.terms[12], k=log(length(step.bic.out$y))) Single term additions Model: pheno.dat ~ rs13479085 + rs13480057 + rs13482096 + rs8254221 Df DevianceAIC 233.82 841.18 rs13480057:rs13475933 6 215.95 851.12 (***) Another problem is that the final model seems to be the first model that satisfies (bAIC >= AIC + 1e-07) if steps haven't used up, rather than the one before that. Please see below. > formula(step2.bic.out) pheno.dat ~ rs13479085 + rs13480057 + rs13482096 + rs8254221 + rs13482096:rs13483699 > step2.bic.out$anova Step Df Deviance Resid. Df Resid. Dev AIC 1 NA NA 298 233.8226 841.1784 Any insights would be greatly appreciated. Thanks much ! Ping Wang __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Legend in scatterplot and adding regression lines
Dear Andrew, The scatterplot() function is in the car package, and xyplot() is in the lattice package. The Rcmdr dialog boxes for these functions don't provide for all possible uses, especially in the case of xyplot(). See the help pages ?scatterplot and ?xyplot for more information. With respect to scatterplot(), if you want a legend in other than the default location in the top margin of the plot, you'll have to add the legend yourself, via the legend() function. As a more general matter, if you want to produce a highly customized graph, it probably makes sense to use R graphics commands to build it up yourself. I hope this helps, John On Tue, 13 May 2008 14:29:07 +1200 "Andrew McFadden" <[EMAIL PROTECTED]> wrote: > Hi R users > > The scatterplot (used package R-Cmdr ) obtained has the legend off > the > plot in the upper left > I would like to decide where to put the legend - is that option still > available and how do I do it? > code: > Dataset <- read.table("J:/Pigmatingsb.txt", header=TRUE, sep="\t", > na.strings="NA", dec=".", strip.white=TRUE) > scatterplot(Farrowing.rate~Time.in.Weeks | Operator, reg.line=lm, > smooth=FALSE, labels=FALSE, boxplots=FALSE, span=0.5, by.groups=TRUE, > ylab="Farrowing Rate", xlab="Time in Weeks", font.lab=2, cex.lab=1.5, > cex.axis=1.2, pch=c(16,15,24,3), data=Dataset) > > The other option was to use the following which allowed manipulation > of > the legend but had its own problems as I wasn't entirely sure how to > produce a simple regression line for the scatter for each group? > > > library(lattice) > setwd("C:\\temp") > Dataset <- read.table("Pigmatingsb.txt", header=TRUE, sep="\t", > na.strings="NA", dec=".", strip.white=TRUE) > head(Dataset) > > xyplot(Farrowing.rate~Time.in.Weeks, groups= Operator, type=c("p"data > = > Dataset,layout = c(1, 1), panel=panel.superpose ,xlab="Time > (weeks)",ylab="Farrowing rate", key = list(points = > Rows(trellis.par.get("superpose.symbol"), 1:4), >text = list(levels(Dataset$Operator)), >columns = 4)) > > Regards > > Andy > > Andrew McFadden MVS BVSc > Incursion Investigator > Investigation & Diagnostic Centres - Wallaceville Biosecurity New > Zealand Ministry of Agriculture and Forestry > > Phone 04 894 5600 Fax 04 894 4973 Mobile 029 894 5611 Postal address: > > Investigation and Diagnostic Centre- Wallaceville Box 40742 Ward St > Upper Hutt > > > > This email message and any attachment(s) is intended solely for the > addressee(s) named above. The information it contains is confidential > and may be legally privileged. Unauthorised use of the message, or > the > information it contains, may be unlawful. If you have received this > message by mistake please call the sender immediately on 64 4 8940100 > or notify us by return email and erase the original message and > attachments. Thank you. > > The Ministry of Agriculture and Forestry accepts no responsibility > for > changes made to this email or to any attachments after transmission > from > the office. > > > [[alternative HTML version deleted]] > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. John Fox, Professor Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Legend in scatterplot and adding regression lines
Hi R users The scatterplot (used package R-Cmdr ) obtained has the legend off the plot in the upper left I would like to decide where to put the legend - is that option still available and how do I do it? code: Dataset <- read.table("J:/Pigmatingsb.txt", header=TRUE, sep="\t", na.strings="NA", dec=".", strip.white=TRUE) scatterplot(Farrowing.rate~Time.in.Weeks | Operator, reg.line=lm, smooth=FALSE, labels=FALSE, boxplots=FALSE, span=0.5, by.groups=TRUE, ylab="Farrowing Rate", xlab="Time in Weeks", font.lab=2, cex.lab=1.5, cex.axis=1.2, pch=c(16,15,24,3), data=Dataset) The other option was to use the following which allowed manipulation of the legend but had its own problems as I wasn't entirely sure how to produce a simple regression line for the scatter for each group? library(lattice) setwd("C:\\temp") Dataset <- read.table("Pigmatingsb.txt", header=TRUE, sep="\t", na.strings="NA", dec=".", strip.white=TRUE) head(Dataset) xyplot(Farrowing.rate~Time.in.Weeks, groups= Operator, type=c("p"data = Dataset,layout = c(1, 1), panel=panel.superpose ,xlab="Time (weeks)",ylab="Farrowing rate", key = list(points = Rows(trellis.par.get("superpose.symbol"), 1:4), text = list(levels(Dataset$Operator)), columns = 4)) Regards Andy Andrew McFadden MVS BVSc Incursion Investigator Investigation & Diagnostic Centres - Wallaceville Biosecurity New Zealand Ministry of Agriculture and Forestry Phone 04 894 5600 Fax 04 894 4973 Mobile 029 894 5611 Postal address: Investigation and Diagnostic Centre- Wallaceville Box 40742 Ward St Upper Hutt This email message and any attachment(s) is intended solely for the addressee(s) named above. The information it contains is confidential and may be legally privileged. Unauthorised use of the message, or the information it contains, may be unlawful. If you have received this message by mistake please call the sender immediately on 64 4 8940100 or notify us by return email and erase the original message and attachments. Thank you. The Ministry of Agriculture and Forestry accepts no responsibility for changes made to this email or to any attachments after transmission from the office. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] help on plot title
Please read the docs (starting with "An Introduction to R") . It is unreasonable to expect to use command line based software without first reading the "how to" manuals. Also, and always, read the man pages: e.g. ?plot, ?plot.default and appropriate links. Finally, as a matter of etiquette, signed emails are more likely to get helpful responses. BTW, is this a homework problem? Cheers, Bert -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Lisa Sent: Monday, May 12, 2008 2:45 PM To: r-help@r-project.org Subject: [R] help on plot title Hi, I am do some plotting and need to add title to each of the plots, Can anyone help me on how to add the variables to the title? here it is the program, the title i want should look like, j=1, j=2.., but if i use title("j=",j), the j will go to the x axis. Can anyone help me on this? thanka s lot for (j in 1:6) { y<-mx.j(x,j) sx<-mx.j(x,j)+sigma*ei plot(x,y, type="l",xlab="x", ylab="") title("j=",j) } [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [OT] xemacs on windows vista
Just in case you think Vincent is pushing his distribution as opposed to some other approach, let me endorse his suggestion as someone who has recently made the change to GNU Emacs. Vincent's distribution is truly very nice. I have installed it on Vista without any problems and for anyone who uses TeX and R it is an excellent option. I haven't checked John Fox's advice concerning XEmacs lately but in the past it has been the best source of information about XEmacs and R if you do still want to stick with XEmacs. David Scott On Mon, 12 May 2008, Vincent Goulet wrote: Le lun. 12 mai à 15:15, Wensui Liu a écrit : Hi, dear all, I just switch to vista (ultimate) and have heard there is some problem for the installation of xemacs on vista. Is there any insight or experience that you could share? Yes: go with GNU Emacs. There doesn't seem to be any compelling reason to prefer XEmacs nowadays. And if you use my distribution http://vgoulet.act.ulaval.ca/en/emacs you get a wizard based installation with up-to-date AUCTeX and ESS built-in. Otherwise, I have never heard of special issues on Vista. HTH Vincent I really appreciate any input. thank you so much! __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. _ David Scott Department of Statistics, Tamaki Campus The University of Auckland, PB 92019 Auckland 1142,NEW ZEALAND Phone: +64 9 373 7599 ext 86830 Fax: +64 9 373 7000 Email: [EMAIL PROTECTED] Graduate Officer, Department of Statistics Director of Consulting, Department of Statistics __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] help on plot title
title(main=paste("j =",j)) On 13/05/2008, at 9:45 AM, Lisa wrote: Hi, I am do some plotting and need to add title to each of the plots, Can anyone help me on how to add the variables to the title? here it is the program, the title i want should look like, j=1, j=2.., but if i use title("j=",j), the j will go to the x axis. Can anyone help me on this? thanka s lot for (j in 1:6) { y<-mx.j(x,j) sx<-mx.j(x,j)+sigma*ei plot(x,y, type="l",xlab="x", ylab="") title("j=",j) } ## Attention:\ This e-mail message is privileged and confid...{{dropped:9}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] combining bar and column graphs?
Michael Kubovy wrote: Perhaps ?mosaic or ?mosaicplot Thanks for these helpful suggestions - I will try them both out - __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] RPM-style install (SLED 10.1)
That's what I ended up doing. I've had complaints about blas, too, but I rpm-ed them out with --nodeps. I've compiled R from the sources (I had to add some -devel- packages, too), and then copied over the stats.*.so library from my fully successful 2.6.2 compilation to semi-successful 2.7.0 installation. So I cannot actually tell if I have 2.7.0 or 2.6.2 on this computer... and I don't know if any other problems may appear later. And yes I have a ThinkPad although it is 32 bit. On 5/12/08, Horace Tso <[EMAIL PROTECTED]> wrote: > Stas, this doesn't solve your problem but may shed some light on what might > have gone wrong. Just recently I installed R 2.7.0 under openSuse 10.3 in a > brand new 64-bit Lenovo ThinkPad. The first install failed because of a > dependency error just like yours. It complained it couldn't find BLAS and > libgfortran. I then got on the varous repositories and installed a couple > flavours (don't remember how many) of fortran and finally it went through. > Right now, in my /usr/lib64, i have libgfortran.so.2.0.0. > > > -Original Message- > From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Stas Kolenikov > Sent: Monday, May 12, 2008 7:33 AM > To: r-help@r-project.org > Subject: [R] RPM-style install (SLED 10.1) > > I am trying to install R on a SLED 10.1 machine. > R-base-2.7.0-7.1-i586.rpm fails with > > [EMAIL PROTECTED]:~/RPMs> rpm -Uvh R-base-2.7.0-7.1.i586.rpm > warning: R-base-2.7.0-7.1.i586.rpm: Header V3 DSA signature: NOKEY, > key ID 14ec5930 > error: Failed dependencies: > libgfortran.so.1 is needed by R-base-2.7.0-7.1.i586 > > I tried to trick it into believing there's the library by setting up the > link: > > [EMAIL PROTECTED]:~/RPMs> ls -l /usr/lib/libgf* > lrwxrwxrwx 1 root root 29 2008-05-08 17:45 > /usr/lib/libgfortran.so.1 -> /usr/lib/libgfortran.so.3.0.0 > lrwxrwxrwx 1 root root 20 2008-05-08 17:37 > /usr/lib/libgfortran.so.3 -> libgfortran.so.3.0.0 > -rwxr-xr-x 1 root root 2251139 2008-05-01 09:43 /usr/lib/libgfortran.so.3.0.0 > > but that did not work, either. I went on and installed GCC's gfortran, > but it did not provide libgfortran.so.1 either: > > [EMAIL PROTECTED]:~/RPMs> ls -l /usr/irun/lib/libgf* > -rw-r--r-- 1 1005 1011 4330512 2008-03-02 02:29 /usr/irun/lib/libgfortran.a > -rwxr-xr-x 1 1005 10111009 2008-03-02 02:29 /usr/irun/lib/libgfortran.la > -rwxr-xr-x 1 1005 1011 2509309 2008-03-02 02:29 /usr/irun/lib/libgfortran.so > -rwxr-xr-x 1 1005 1011 2509309 2008-03-02 02:29 > /usr/irun/lib/libgfortran.so.3 > -rwxr-xr-x 1 1005 1011 2509309 2008-03-02 02:29 > /usr/irun/lib/libgfortran.so.3.0.0 > > Any reason R wants that package instead of the newer one? Of course > rpm --nodeps was an option, but then it fails to load r-stats with the > same message about the missing library: > > Error in dyn.load(file, DLLpath = DLLpath, ...) : > unable to load shared library '/usr/lib/R/library/stats/libs/stats.so': > /usr/lib/R/library/stats/libs/stats.so: undefined symbol: > _gfortran_pow_r8_i4 > During startup - Warning message: > package stats in options("defaultPackages") was not found > > I am stuck... I am trying to compile R from 2.6.2 sources in the > meantime, but that's painful for a Linux newbie like me on a brand-new > computer that does not necessarily have all the packages. I almost > want to switch to Ubuntu as R installation back there was a single > click :)) > > -- > Stas Kolenikov, also found at http://stas.kolenikov.name > Small print: Please do not reply to my Gmail address as I don't check > it regularly. > > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > -- Stas Kolenikov, also found at http://stas.kolenikov.name Small print: Please do not reply to my Gmail address as I don't check it regularly. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [R-sig-ME] lme nesting/interaction advice
On 13/05/2008, at 4:09 AM, Douglas Bates wrote: I'm entering this discussion late so I may be discussing issues that have already been addressed. As I understand it, Federico, you began by describing a model for data in which two factors have a fixed set of levels and one factor has an extensible, or "random", set of levels and you wanted to fit a model that you described as y ~ effect1 * effect2 * effect3 The problem is that this specification is not complete. At *last* (as Owl said to Rabbit) we're getting somewhere!!! I always knew that there was some basic fundamental point about this business that I (and I believe many others) were simply missing. But I could not for the life of me get anyone to explain to me what that point was. Or to put it another way, I was never able to frame a question that would illuminate just what it was that I wasn't getting. I now may be at a stage where I can start asking the right questions. An interaction of factors with fixed levels and a factor with random levels can mean, in the lmer specification, lmer(y ~ effect1 * effect2 + (1| effect3) + (1| effect1:effect2:effect3), ...) or lmer(y ~ effect1 * effect2 + (effect1*effect2 | effect3), ...) or other variations. When you specify a random effect or an random interaction term you must, either explicitly or implicitly, specify the form of the variance-covariance matrix associated with those random effects. The "advantage" that other software may provide for you is that it chooses the model for you but that, of course, means that you only have the one choice. Now may I start asking what I hope are questions that will lift the fog a bit? Let us for specificity consider a three-way model with two fixed effects and one random effect from the good old Rothamstead style agricultural experiment context: Suppose we have a number of species/breeds of wheat (say) and a number of fertilizers. These are fixed effects. And we have a number of fields (blocks) --- a random effect. Each breed-fertilizer combination is applied a number of times in each field. We ***assume*** that that the field or block effect is homogeneous throughout. This may or may not be a ``good'' assumption, but it's not completely ridiculous and would often be made in practice. And probably *was* made at Rothamstead. The response would be something like yield in bushels per acre. The way that I would write the ``full'' model for this setting, in mathematical style is: Y_ijkl = mu + alpha_i + beta_j + (alpha.beta)_ij + C_k + (alpha.C)_ik + (beta.C)_jk + (alpha.beta.C)_ijk + E_ijkl The alpha_i and beta_j are parameters corresponding to breed and fertilizer respectively; the C_k are random effects corresponding to fields or blocks. Any effect ``involving'' C is also random. The assumptions made by the Package-Which-Must-Not-Be-Named are (I think) that C_k ~ N(0,sigma_C^2) (alpha.C)_ik ~ N(0,sigma_aC^2) (beta.C)jk ~ N(0,sigma_bC^2) (alpha.beta.C)_ijk ~ N(0,sigma_abC^2) E_ijkl ~ N(0,sigma^2) and these random variables are *all independent*. A ... perhaps I'm on the way to answering my own question. Is it this assumption of ``all independent'' which is questionable? It seemed innocent enough when I first learned about this stuff, lo these many years ago. But maybe not! To start with: What would be the lmer syntax to fit the foregoing (possibly naive) model? I am sorry, but I really cannot get my head around the syntax of lmer model specification, and I've tried. I really have. Hard. I know I must be starting from the wrong place, but I haven't a clue as to what the *right* place to start from is. And if I'm in that boat, I will wager Euros to pretzels that there are others in it. I know that I'm not the brightest bulb in the chandelier, but I'm not the dullest either. Having got there: Presuming that I'm more-or-less on the right track in my foregoing conjecture that it's the over-simple dependence structure that is the problem with what's delivered by the Package-Which-Must- Not-Be-Named, how might one go about being less simple-minded? I.e. what might be some more realistic dependence structures, and how would one specify these in lmer? And how would one assess whether the assumed dependence structure gives a reasonable fit to the data? If you can describe how many variance components you think should be estimated in your model and what they would represent then I think it will be easier to describe how to fit the model.
Re: [R] Format integer
I guess "little" means different things to different people: x = sample(1:100,65,replace=TRUE) system.time(a<-formatC(x,digits=10,flag='0')) user system elapsed 32.854 0.444 34.813 system.time(b<-sprintf("%011d",x)) user system elapsed 0.352 0.012 0.363 If you look at the definitions of the functions, you'll see that formatC is written in R, and sprintf uses a single call to an .Internal function. I - Phil Spector Statistical Computing Facility Department of Statistics UC Berkeley [EMAIL PROTECTED] On Mon, 12 May 2008, Anh Tran wrote: Yea, thanks all. I checked back and I got a few things mistyped. The array is 650,000 and it took 25 seconds :p. It's acceptable. Just that I had too many variable at the time I ran it. Also, seems like sprintf is a little faster. Thanks all. Anh Tran On Mon, May 12, 2008 at 2:55 PM, Uwe Ligges <[EMAIL PROTECTED]> wrote: Anh Tran wrote: Thanks. formatC(flag) works. But it's awefully slow. I try to do that for 65000 numbers (generating ID for each item) and it seems like forever. On my not that recent laptop: system.time(formatC(1:65000, width=10, flag="0")) user system elapsed 1.920.001.94 I think 2 seconds is less than "forever". Uwe Ligges Is there any faster way? Thank all. Anh Tran On Mon, May 12, 2008 at 2:36 PM, Uwe Ligges < [EMAIL PROTECTED]> wrote: Anh Tran wrote: Hi, What's one way to convert an integer to a string with preceding 0's? such that '13' becomes '013' to be put into a string I've tried formatC, but they removes all the zeros and replace it with blanks Not so for me: formatC(13, digits=10, flag="0") Uwe LIgges Thanks -- Regards, Anh Tran [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Format integer
Yea, thanks all. I checked back and I got a few things mistyped. The array is 650,000 and it took 25 seconds :p. It's acceptable. Just that I had too many variable at the time I ran it. Also, seems like sprintf is a little faster. Thanks all. Anh Tran On Mon, May 12, 2008 at 2:55 PM, Uwe Ligges <[EMAIL PROTECTED]> wrote: > > > Anh Tran wrote: > > > Thanks. formatC(flag) works. > > > > But it's awefully slow. I try to do that for 65000 numbers (generating > > ID > > for each item) and it seems like forever. > > > > On my not that recent laptop: > > > system.time(formatC(1:65000, width=10, flag="0")) > user system elapsed > 1.920.001.94 > > > I think 2 seconds is less than "forever". > > Uwe Ligges > > > > > > > Is there any faster way? > > > > Thank all. > > > > Anh Tran > > > > On Mon, May 12, 2008 at 2:36 PM, Uwe Ligges < > > [EMAIL PROTECTED]> wrote: > > > > > > > Anh Tran wrote: > > > > > > Hi, > > > > What's one way to convert an integer to a string with preceding 0's? > > > > such that > > > > '13' becomes '013' > > > > to be put into a string > > > > > > > > I've tried formatC, but they removes all the zeros and replace it > > > > with > > > > blanks > > > > > > > > Not so for me: > > > > > > formatC(13, digits=10, flag="0") > > > > > > Uwe LIgges > > > > > > > > > > > > Thanks > > > > > > > > > > > > > > > > -- Regards, Anh Tran [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Format integer
Anh Tran wrote: Hi, What's one way to convert an integer to a string with preceding 0's? such that '13' becomes '013' to be put into a string I've tried formatC, but they removes all the zeros and replace it with blanks Not so for me: formatC(13, digits=10, flag="0") Uwe LIgges Thanks __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Format integer
Anh Tran wrote: Thanks. formatC(flag) works. But it's awefully slow. I try to do that for 65000 numbers (generating ID for each item) and it seems like forever. On my not that recent laptop: > system.time(formatC(1:65000, width=10, flag="0")) user system elapsed 1.920.001.94 I think 2 seconds is less than "forever". Uwe Ligges Is there any faster way? Thank all. Anh Tran On Mon, May 12, 2008 at 2:36 PM, Uwe Ligges < [EMAIL PROTECTED]> wrote: Anh Tran wrote: Hi, What's one way to convert an integer to a string with preceding 0's? such that '13' becomes '013' to be put into a string I've tried formatC, but they removes all the zeros and replace it with blanks Not so for me: formatC(13, digits=10, flag="0") Uwe LIgges Thanks __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Format integer
On May 12, 2008, at 5:22 PM, Anh Tran wrote: Hi, What's one way to convert an integer to a string with preceding 0's? such that '13' becomes '013' to be put into a string I've tried formatC, but they removes all the zeros and replace it with blanks formatC(13, width=10, format="d", flag="0") Thanks -- Regards, Anh Tran Haris Skiadas Department of Mathematics and Computer Science Hanover College __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] test
You probably should check this section in your R-help subscription options (via https://stat.ethz.ch/mailman/options/r-help/, I think): Receive your own posts to the list? Ordinarily, you will get a copy of every message you post to the list. If you don't want to receive this copy, set this option to No. I see 5 identical posts with the subject "Several questions about MCMClogit" on R-help recently. -- Tony Plate j t wrote: Sorry to bother your. I am trying to post my question for more than 10 times, but I still didn't see it. It drives my crazy!!! It is a test for posting some simple pure text. Chao __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] help on plot title
Hi, I am do some plotting and need to add title to each of the plots, Can anyone help me on how to add the variables to the title? here it is the program, the title i want should look like, j=1, j=2.., but if i use title("j=",j), the j will go to the x axis. Can anyone help me on this? thanka s lot for (j in 1:6) { y<-mx.j(x,j) sx<-mx.j(x,j)+sigma*ei plot(x,y, type="l",xlab="x", ylab="") title("j=",j) } [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Format integer
Thanks. formatC(flag) works. But it's awefully slow. I try to do that for 65000 numbers (generating ID for each item) and it seems like forever. Is there any faster way? Thank all. Anh Tran On Mon, May 12, 2008 at 2:36 PM, Uwe Ligges < [EMAIL PROTECTED]> wrote: > > > Anh Tran wrote: > > > Hi, > > What's one way to convert an integer to a string with preceding 0's? > > such that > > '13' becomes '013' > > to be put into a string > > > > I've tried formatC, but they removes all the zeros and replace it with > > blanks > > > > Not so for me: > > formatC(13, digits=10, flag="0") > > Uwe LIgges > > > > > Thanks > > > > -- Regards, Anh Tran [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] is.category
Hello, Could someone tell me what the SPLUS "is.category" function do and what is its equivalent in R? Thank you, I couldn't find any help elsewhere... -- View this message in context: http://www.nabble.com/is.category-tp1719p1719.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Format integer
Try something like one of these (as documented in ?formatC) > formatC(13, flag="0", width=10) [1] "13" > sprintf("%010g", 13) [1] "13" > Anh Tran wrote: Hi, What's one way to convert an integer to a string with preceding 0's? such that '13' becomes '013' to be put into a string I've tried formatC, but they removes all the zeros and replace it with blanks Thanks __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Converting qqplot2 qplot() to grammar?
Hello all, I've been using the following qplot command: qplot(pixX,pixY, data=som, geom="tile", fill=rgb) + scale_fill_identity() + opts(aspect.ratio = .75) + facet_grid(unitX ~ unitY) Now I would like to convert it into the explicit ggplot grammar, so I can remove the extras: axes, labels, background, borders, facet labels, and extra white-space around the plot. (If anyone has suggestions on removing these please let me know) I've come up with the following but it is not behaving the same way as the qplot above: ggplot(data = somdf, mapping = aes(x = pixX, y = pixY)) + layer(data = somdf, geom = "tile", fill=rgb) + scale_y_continuous(name=" ",breaks=" ") + scale_x_continuous(name=" ",breaks=" ") + scale_fill_identity() + coord_cartesian() + opts(aspect.ratio = .75) + facet_grid(unitX ~ unitY,margins=FALSE) The result is a plot where all tiles are filled with grey50, and not the data values. I've also tried this variation with the same results: ggplot(data = somdf, mapping = aes(x = pixX, y = pixY)) + geom_tile(data = somdf, fill=rgb) + scale_y_continuous(name=" ",breaks=" ") + scale_x_continuous(name=" ",breaks=" ") + scale_fill_identity() + coord_cartesian() + opts(aspect.ratio = .75) + facet_grid(unitX ~ unitY,margins=FALSE) The other issue I'm having is that the above seems to create multiple plots, another, apparently identical, plot gets drawn after I call dev.off(). I have no idea why this would be so. Setting name and breaks to " " seems like a hack to get rid of the axis stuff, is there a better way? Oh, and I can't find documentation for opts() on the ggplot2 website, where is it available? Thanks all, Hadley in particular, B. Bogart __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Format integer
Hi, What's one way to convert an integer to a string with preceding 0's? such that '13' becomes '013' to be put into a string I've tried formatC, but they removes all the zeros and replace it with blanks Thanks -- Regards, Anh Tran [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] help with rpart
Some answers are on the help pages for plot.rpart and text.rpart. On Mon, 12 May 2008, Linus An wrote: Hi, I am using rpart as a part of my masters' project. I am trying to print out the resulting model using plot() function along with text() function. I am having difficulties with labels being cut-off. In text() function, I am using use.n=T option to get the number of people in each nodes but the on the lower and left part of the plot, the numbers get cut off. Thanks! Linus [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. PLEASE do. -- 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@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] test
Sorry to bother your. I am trying to post my question for more than 10 times, but I still didn't see it. It drives my crazy!!! It is a test for posting some simple pure text. Chao __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [OT] xemacs on windows vista
Wensui Liu wrote: Hi, dear all, I just switch to vista (ultimate) and have heard there is some problem for the installation of xemacs on vista. Is there any insight or experience that you could share? I really appreciate any input. thank you so much! Hi, I don't know about XEmacs, but I am using a specially bundled version of GNU Emacs by Vincet Goulet which is great. I was using a plain-vanilla version of GNU Emacs before, but this one has ESS and spell check support built in :-) I am using it under XP .. I believe this also works under Vista Check out: http://vgoulet.act.ulaval.ca/en/ressources/emacs/ HTH, Esmail __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Several questions about MCMClogit
Hello everybody, I'm new to MCMClogit. I'm trying to use MCMClogit to fit a logistic regression model but I got some warnings I can't understand. My input data X is 32(tissue sample)*20(genes) matrix, each element in this matrix corresponds to the expression value of one particular gene in one of 32 samples. And the Y presents the corresponding classes (0-non cancer, 1-cancer) for those 32 samples. The formula is Y~. All other parameters are default. This is the output message: @@ The Metropolis acceptance rate for beta was 0.0 @@ Warning messages: 1: In glm.fit(x = X, y = Y, weights = weights, start = start, etastart = etastart, : fitted probabilities numerically 0 or 1 occurred 2: In glm.fit(x = X, y = Y, weights = weights, start = start, etastart = etastart, : fitted probabilities numerically 0 or 1 occurred Question 1: What are the meanings of two warning messages? Question 2: How could I use the regression coefficients to predict for other data, in other words, how could I extract those regression coefficients from the result of MCMClogit. I know maybe my questions are some basic but it already bothered me for several days. I hope somebody can give me some hint about them. BTW, is there any specific manual or examples for using MCMClogit, the manual in package is not enough to figure out my problem. Thanks and regards, Chao __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Several questions about MCMClogit
Hello everybody, I'm new to MCMClogit. I'm trying to use MCMClogit to fit a logistic regression model but I got some warnings I can't understand. My input data X is 32(tissue sample)*20(genes) matrix, each element in this matrix corresponds to the expression value of one particular gene in one of 32 samples. And the Y presents the corresponding classes (0-non cancer, 1-cancer) for those 32 samples. The formula is Y~. All other parameters are default. This is the output message: @@ The Metropolis acceptance rate for beta was 0.0 @@ Warning messages: 1: In glm.fit(x = X, y = Y, weights = weights, start = start, etastart = etastart, : fitted probabilities numerically 0 or 1 occurred 2: In glm.fit(x = X, y = Y, weights = weights, start = start, etastart = etastart, : fitted probabilities numerically 0 or 1 occurred Question 1: What are the meanings of two warning messages? Question 2: How could I use the regression coefficients to predict for other data, in other words, how could I extract those regression coefficients from the result of MCMClogit. I know maybe my questions are some basic but it already bothered me for several days. I hope somebody can give me some hint about them. BTW, is there any specific manual or examples for using MCMClogit, the manual in package is not enough to figure out my problem. Thanks and regards, Chao __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Several questions about MCMClogit
Hello everybody, I'm new to MCMClogit. I'm trying to use MCMClogit to fit a logistic regression model but I got some warnings I can't understand. My input data X is 32(tissue sample)*20(genes) matrix, each element in this matrix corresponds to the expression value of one particular gene in one of 32 samples. And the Y presents the corresponding classes (0-non cancer, 1-cancer) for those 32 samples. The formula is Y~. All other parameters are default. This is the output message: @@ The Metropolis acceptance rate for beta was 0.0 @@ Warning messages: 1: In glm.fit(x = X, y = Y, weights = weights, start = start, etastart = etastart, : fitted probabilities numerically 0 or 1 occurred 2: In glm.fit(x = X, y = Y, weights = weights, start = start, etastart = etastart, : fitted probabilities numerically 0 or 1 occurred Question 1: What are the meanings of two warning messages? Question 2: How could I use the regression coefficients to predict for other data, in other words, how could I extract those regression coefficients from the result of MCMClogit. I know maybe my questions are some basic but it already bothered me for several days. I hope somebody can give me some hint about them. BTW, is there any specific manual or examples for using MCMClogit, the manual in package is not enough to figure out my problem. Thanks and regards, Chao [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Converting qqplot2 qplot() to grammar?
Hello all, I've been using the following qplot command: qplot(pixX,pixY, data=som, geom="tile", fill=rgb) + scale_fill_identity() + opts(aspect.ratio = .75) + facet_grid(unitX ~ unitY) Now I would like to convert it into the explicit ggplot grammar, so I can remove the extras: axes, labels, background, borders, facet labels, and extra white-space around the plot. (If anyone has suggestions on removing these please let me know) I've come up with the following but it is not behaving the same way as the qplot above: ggplot(data = somdf, mapping = aes(x = pixX, y = pixY)) + layer(data = somdf, geom = "tile", fill=rgb) + scale_y_continuous(name=" ",breaks=" ") + scale_x_continuous(name=" ",breaks=" ") + scale_fill_identity() + coord_cartesian() + opts(aspect.ratio = .75) + facet_grid(unitX ~ unitY,margins=FALSE) The result is a plot where all tiles are filled with grey50, and not the data values. I've also tried this variation with the same results: ggplot(data = somdf, mapping = aes(x = pixX, y = pixY)) + geom_tile(data = somdf, fill=rgb) + scale_y_continuous(name=" ",breaks=" ") + scale_x_continuous(name=" ",breaks=" ") + scale_fill_identity() + coord_cartesian() + opts(aspect.ratio = .75) + facet_grid(unitX ~ unitY,margins=FALSE) The other issue I'm having is that the above seems to create multiple plots, another, apparently identical, plot gets drawn after I call dev.off(). I have no idea why this would be so. Setting name and breaks to " " seems like a hack to get rid of the axis stuff, is there a better way? Oh, and I can't find documentation for opts() on the ggplot2 website, where is it available? Thanks all, Hadley in particular, B. Bogart __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Several questions about MCMClogit
Hello everybody, I'm new to MCMClogit. I'm trying to use MCMClogit to fit a logistic regression model but I got some warnings I can't understand. My input data X is 32(tissue sample)*20(genes) matrix, each element in this matrix corresponds to the expression value of one particular gene in one of 32 samples. And the Y presents the corresponding classes (0-non cancer, 1-cancer) for those 32 samples. The formula is Y~. All other parameters are default. This is the output message: @@ The Metropolis acceptance rate for beta was 0.0 @@ Warning messages: 1: In glm.fit(x = X, y = Y, weights = weights, start = start, etastart = etastart, : fitted probabilities numerically 0 or 1 occurred 2: In glm.fit(x = X, y = Y, weights = weights, start = start, etastart = etastart, : fitted probabilities numerically 0 or 1 occurred Question 1: What are the meanings of two warning messages? Question 2: How could I use the regression coefficients to predict for other data, in other words, how could I extract those regression coefficients from the result of MCMClogit. I know maybe my questions are some basic but it already bothered me for several days. I hope somebody can give me some hint about them. BTW, is there any specific manual or examples for using MCMClogit, the manual in package is not enough to figure out my problem. Thanks and regards, Chao [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] RPM-style install (SLED 10.1)
Stas, this doesn't solve your problem but may shed some light on what might have gone wrong. Just recently I installed R 2.7.0 under openSuse 10.3 in a brand new 64-bit Lenovo ThinkPad. The first install failed because of a dependency error just like yours. It complained it couldn't find BLAS and libgfortran. I then got on the varous repositories and installed a couple flavours (don't remember how many) of fortran and finally it went through. Right now, in my /usr/lib64, i have libgfortran.so.2.0.0. HTH. Horace -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Stas Kolenikov Sent: Monday, May 12, 2008 7:33 AM To: r-help@r-project.org Subject: [R] RPM-style install (SLED 10.1) I am trying to install R on a SLED 10.1 machine. R-base-2.7.0-7.1-i586.rpm fails with [EMAIL PROTECTED]:~/RPMs> rpm -Uvh R-base-2.7.0-7.1.i586.rpm warning: R-base-2.7.0-7.1.i586.rpm: Header V3 DSA signature: NOKEY, key ID 14ec5930 error: Failed dependencies: libgfortran.so.1 is needed by R-base-2.7.0-7.1.i586 I tried to trick it into believing there's the library by setting up the link: [EMAIL PROTECTED]:~/RPMs> ls -l /usr/lib/libgf* lrwxrwxrwx 1 root root 29 2008-05-08 17:45 /usr/lib/libgfortran.so.1 -> /usr/lib/libgfortran.so.3.0.0 lrwxrwxrwx 1 root root 20 2008-05-08 17:37 /usr/lib/libgfortran.so.3 -> libgfortran.so.3.0.0 -rwxr-xr-x 1 root root 2251139 2008-05-01 09:43 /usr/lib/libgfortran.so.3.0.0 but that did not work, either. I went on and installed GCC's gfortran, but it did not provide libgfortran.so.1 either: [EMAIL PROTECTED]:~/RPMs> ls -l /usr/irun/lib/libgf* -rw-r--r-- 1 1005 1011 4330512 2008-03-02 02:29 /usr/irun/lib/libgfortran.a -rwxr-xr-x 1 1005 10111009 2008-03-02 02:29 /usr/irun/lib/libgfortran.la -rwxr-xr-x 1 1005 1011 2509309 2008-03-02 02:29 /usr/irun/lib/libgfortran.so -rwxr-xr-x 1 1005 1011 2509309 2008-03-02 02:29 /usr/irun/lib/libgfortran.so.3 -rwxr-xr-x 1 1005 1011 2509309 2008-03-02 02:29 /usr/irun/lib/libgfortran.so.3.0.0 Any reason R wants that package instead of the newer one? Of course rpm --nodeps was an option, but then it fails to load r-stats with the same message about the missing library: Error in dyn.load(file, DLLpath = DLLpath, ...) : unable to load shared library '/usr/lib/R/library/stats/libs/stats.so': /usr/lib/R/library/stats/libs/stats.so: undefined symbol: _gfortran_pow_r8_i4 During startup - Warning message: package stats in options("defaultPackages") was not found I am stuck... I am trying to compile R from 2.6.2 sources in the meantime, but that's painful for a Linux newbie like me on a brand-new computer that does not necessarily have all the packages. I almost want to switch to Ubuntu as R installation back there was a single click :)) -- Stas Kolenikov, also found at http://stas.kolenikov.name Small print: Please do not reply to my Gmail address as I don't check it regularly. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] writing a table in the device (pdf in this case)
You probably want to look at "Sweave" in the "utils" package or the "odfWeave" package. Both let you set up a planned set of commands interspersed with text (notes, explanations, full report, etc.) and then you process the file and get the output (and commands) in either a LaTeX file or an OpenOffice writer document, either of which can then be processed to a pdf file. Both give several options for what input/output to return, inclusion of tables and graphs, and much more. Hope this helps, -- Gregory (Greg) L. Snow Ph.D. Statistical Data Center Intermountain Healthcare [EMAIL PROTECTED] (801) 408-8111 > -Original Message- > From: [EMAIL PROTECTED] > [mailto:[EMAIL PROTECTED] On Behalf Of Julien Colomb > Sent: Saturday, May 10, 2008 9:47 AM > To: r-help@r-project.org > Subject: [R] writing a table in the device (pdf in this case) > > hello all, > I would like to introduce a summary table into the pdf along > with the plots (in order to archive my data into single files > automatically). > Similarly, It would be great to have the result of the > statistical analysis (for instance anova) in the same file. > > Is there a way to do that? > > > example: > pdf("example.pdf") > layout (matrix(1:2,1,2)) > plot (groups, scores) > > Result <- cbind (samplesize, mean_score, sem_score) > #calculated before ??? > > dev.off () > -- > ___ > > Dr. Julien Colomb > > Genes et Dynamique des Systemes de Memoire CNRS UMR 7637 > ESPCI 10 rue Vauquelin > 75005 Paris > France > www.gdsm.espci.fr > > tel +33 (0)1 40 79 51 23 > fax +33 (0)1 40 79 52 29 > > AIM (ichat) or skype: thrawny1 > [[alternative HTML version deleted]] > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] help with rpart
Hi, I am using rpart as a part of my masters' project. I am trying to print out the resulting model using plot() function along with text() function. I am having difficulties with labels being cut-off. In text() function, I am using use.n=T option to get the number of people in each nodes but the on the lower and left part of the plot, the numbers get cut off. Thanks! Linus [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Several questions about MCMClogit
Hello everybody, I'm new to MCMClogit. I'm trying to use MCMClogit to fit a logistic regression model but I got some warnings I can't understand. My input data X is 32(tissue sample)*20(genes) matrix, each element in this matrix corresponds to the expression value of one particular gene in one of 32 samples. And the Y presents the corresponding classes (0-non cancer, 1-cancer) for those 32 samples. The formula is Y~. All other parameters are default. This is the output message: @@ The Metropolis acceptance rate for beta was 0.0 @@ Warning messages: 1: In glm.fit(x = X, y = Y, weights = weights, start = start, etastart = etastart, : fitted probabilities numerically 0 or 1 occurred 2: In glm.fit(x = X, y = Y, weights = weights, start = start, etastart = etastart, : fitted probabilities numerically 0 or 1 occurred Question 1: What are the meanings of two warning messages? Question 2: How could I use the regression coefficients to predict for other data, in other words, how could I extract those regression coefficients from the result of MCMClogit. I know maybe my questions are some basic but it already bothered me for several days. I hope somebody can give me some hint about them. BTW, is there any specific manual or examples for using MCMClogit, the manual in package is not enough to figure out my problem. Thanks and regards, Chao [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [OT] xemacs on windows vista
Le lun. 12 mai à 15:15, Wensui Liu a écrit : Hi, dear all, I just switch to vista (ultimate) and have heard there is some problem for the installation of xemacs on vista. Is there any insight or experience that you could share? Yes: go with GNU Emacs. There doesn't seem to be any compelling reason to prefer XEmacs nowadays. And if you use my distribution http://vgoulet.act.ulaval.ca/en/emacs you get a wizard based installation with up-to-date AUCTeX and ESS built-in. Otherwise, I have never heard of special issues on Vista. HTH Vincent I really appreciate any input. thank you so much! __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Left censored responses in mixed effects models
Bert Gunter gene.com> writes: > Dear R Fellow-Travellers: > > What is your recommended way of dealing with a left-censored response > (non-detects) in (linear Gaussian) mixed effects models? Your description of the data calls for a tobit model http://en.wikipedia.org/wiki/Tobit_model I think you need to take a look in the survival package. Regards, Gregor __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] lexicographic comparison of two vectors
On 5/12/2008 2:58 PM, Gabriel Valiente wrote: Is there any built-in way to lexicographically compare two vectors of the same length in R? The textbook algorithm could be coded as follows: lex.cmp <- function (vec1,vec2) { for (j in 1:length(vec1)) { if (vec1[j] < vec2[j]) { return(-1) } if (vec1[j] > vec2[j]) { return(1) } } return(0) } I don't think there's any standard function for this. You could write one as above, or slightly faster as lex.cmp <- function(vec1, vec2) { index <- which.min(vec1 == vec2) # find the first diff sign(vec1[index] - vec2[index]) # assumes numeric } If you don't want to assume numeric data, you may need to expand that last line to a series of comparisons like yours, but with index in place of j, e.g. if (vec1[index] < vec2[index]) { return(-1) } if (vec1[index] > vec2[index]) { return(1) } return(0) (unless there's a compare function in some package or other.) Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] [OT] xemacs on windows vista
Hi, dear all, I just switch to vista (ultimate) and have heard there is some problem for the installation of xemacs on vista. Is there any insight or experience that you could share? I really appreciate any input. thank you so much! __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] lexicographic comparison of two vectors
Is there any built-in way to lexicographically compare two vectors of the same length in R? The textbook algorithm could be coded as follows: lex.cmp <- function (vec1,vec2) { for (j in 1:length(vec1)) { if (vec1[j] < vec2[j]) { return(-1) } if (vec1[j] > vec2[j]) { return(1) } } return(0) } Thanks, Gabriel __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Left censored responses in mixed effects models
I assume you've looked at the NADA package(?) While I don't believe it goes as far as dealing the mixed effects models, it might give you a starting point, and possibly some additional references. -Don At 9:08 AM -0700 5/12/08, Bert Gunter wrote: Dear R Fellow-Travellers: What is your recommended way of dealing with a left-censored response (non-detects) in (linear Gaussian) mixed effects models? Specifics: Response is a numeric positive measurement (of volume, actually); but when it falls below some unknown and slightly random value (depending on how the sample is prepared and measured), it cannot be measured and is recorded as 0. There is some statistical literature on this, but I was unable to find anything that appeared to me to implement a strategy in any R package. If it matters, I am less interested in inference than in removing possible bias in estimation. Feel free to respond off-list if you feel that this would not be of general interest. Cheers, Bert Gunter Genentech __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- -- Don MacQueen Environmental Protection Department Lawrence Livermore National Laboratory Livermore, CA, USA 925-423-1062 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to not sort factors when plotting
Thanks for all of the help! Everything's working beautifully now, and I've accomplished in a few hours what it takes most of my colleagues weeks to do, so I think I'll stick with R after all! Lydia On Mon, May 12, 2008 at 4:41 PM, Jorge Ivan Velez <[EMAIL PROTECTED]> wrote: > > > Hi Lydia, > > I compared my ratio function with Dimitris and Phil's suggestions. Please do > NOT use my approach because it's painfully slow for a large vector (as Phil > told me). Here is why (using Win XP SP2, Intel Core- 2 Duo 2.4 GHz, R 2.7.0 > Patched): > > > # Vector > x=rnorm(10,0,1) > > # Suggestion > new.ratio=function(x) x[2:length(x)]/x[1:(length(x)-1)] > > # My horrible function > my.ratio=function(x){ > > temp=NULL > for (n in 1:length(x)) temp=c(temp,x[n]/x[n-1]) > temp > } > > # System time > t=system.time(my.ratio(x)) > tnr=system.time(new.ratio(x)) > t >user system elapsed > 38.790.06 39.31 > tnr >user system elapsed > 0 0 0 > > > Thanks to all, > > Jorge > > > > On Mon, May 12, 2008 at 11:15 AM, Phil Spector <[EMAIL PROTECTED]> > wrote: > > Another alternative would be to take advantage of R's vectorization: > > > > > > > > > > > > > x=c(1,2,3,2,1,2,3) > > > x[2:length(x)]/x[1:(length(x)-1)] > > > > > > > > > [1] 2.000 1.500 0.667 0.500 2.000 1.500 > > > > The solution using your ratio function will be painfully slow > > for a large vector. > > > > - Phil Spector > > Statistical Computing Facility > > Department of Statistics > > UC Berkeley > > [EMAIL PROTECTED] > > > > > __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Compact Indicator Matrices
Thanks. It works! I think I found another solution, working straight with the indicator matrix. > count <- factor(table(apply(ind, 1, paste, collapse=""))) However, that way I can't store the indices of the collapsed rows. -Angelos Markos __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Limiting size of pairs plots
Dear R-users, This is a follow-up on a quite old post of mine which dealt with margins in the pairs function. I thought my problem was solved, but it doesn't seem so (see the code below). I use the pairs function to produce matrix plots, where distinct groups are represented by dots of different colors within each panel. My troubles start when I add a legend at the bottom of the plot. The steps are: 1- I build my dataframe 2- I plot the graph and define the margins within pairs with oma (I want to add 5 + my number of group lines at the bottom) 3- I use gridBase functions to add the legend (symbol + text) The result is an overlay of the legend on the plot. It looks like the lines do not have the same size when using pairs or the gridBase functions. Could anyone explain this discrepancy to me? Thank you in advance for your help. Sebastien ### START library(gridBase) my.Group.number=12 my.df<-data.frame(rnorm(my.Group.number*20, mean=5, sd=2),rnorm(my.Group.number*20, mean=4, sd=1.6),rnorm(my.Group.number*20, mean=8, sd=3),rnorm(my.Group.number*20, mean=6, sd=1),rep(letters[1:my.Group.number],each=20)) names(my.df)<-c("Var1","Var2","Var2","Var2","Group") pairs(my.df[1:4], main = "My Title", pch = 21, col =rainbow(n=my.Group.number)[unclass(my.df$Group)], bg = rainbow(n=my.Group.number)[unclass(my.df$Group)], oma = c(5 + my.Group.number, 3, 5, 3) ) mylegend = paste("Group ",letters[1:my.Group.number],sep="") mylegend.width = strwidth(mylegend[which.max(nchar(mylegend))], "figure", mylegend.Cex(mylegend,myTarget = 0.90)) vps <- baseViewports() pushViewport(vps$inner) for (j in 1:my.Group.number) { grid.text(mylegend[j], x = unit((1-mylegend.width)/2,"npc"), y = unit(1 + my.Group.number - j,"lines"), gp = gpar(cex = 1), just = "left", draw=TRUE) grid.points(x = unit((1-mylegend.width)/2,"npc")-convertUnit(unit(0.5,"lines"),"npc"), y = unit(1 + my.Group.number - j,"lines"), size = unit(0.5, "lines"), default.units = "lines", pch = 21, gp = gpar(col = rainbow(n=my.Group.number)[j],fill = rainbow(n=my.Group.number)[j]), draw=TRUE) } END Prof Brian Ripley a écrit : On Tue, 28 Aug 2007, Sébastien wrote: Thanks you very much. I did not thought about calling oma inside the pairs function... Do you know by any chance why 'pairs' behaves differently from 'plot', with regards to the plot dimension ? Well, it is an array of figures so why should it be the same? See chapter 12 of 'An Introduction to R'. Prof Brian Ripley a écrit : From ?pairs The graphical parameter 'oma' will be set by 'pairs.default' unless supplied as an argument. so try pairs(iris[1:4], main = "Anderson's Iris Data -- 3 species", pch = 21, bg = c("red", "green3", "blue")[unclass(iris$Species)], oma = c(8,3,5,3)) On Tue, 28 Aug 2007, Sébastien wrote: Dear R-users, I would like to add a legend at the bottom of pairs plots (it's my first use of this function). With the plot function, I usually add some additional space at the bottom when I define the size of the graphical device (using mar); grid functions then allows me to draw my legend as I want. Unfortunatley, this technique does not seem to work with the pairs function as the generated plots use all the available space on the device (see below). I guess I am missing a key argument... my attempts to modify the oma, mar, usr arguments were unsuccesfull, and I could not find any helpful threads on the archives. As usual, any advice would be greatly appreciated Sebastien pdf(file="C:/test.pdf", width=6, height= 6 + 0.2*6) par(mar=c(5 + 6,4,4,2)+0.1) pairs(iris[1:4], main = "Anderson's Iris Data -- 3 species", pch = 21, bg = c("red", "green3", "blue")[unclass(iris$Species)]) dev.off() __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] what kind of residuals are the ones calculated in coxph?
Hi Gurus: In the coxph() objects in Survival package, there is an attribute called residuals. Usually, there are several kinds for censored survival data. I can't seem to find in the documentation as to which one this is calculating. Anyone knows? Karen _ [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] predicting from coxph with pspline
-== begin included message - Hello. I get a bit confused by the output from the predict function when used on an object from coxph in combination with p-spline, e.g. fit <- coxph(Surv(time1, time2, status)~pspline(x), Data) predict(fit, newdata=data.frame(x=1:2)) - end included -- Yes, you should be confused. Coxph still retains a well known problem with S models, namely that the prediction is incorrect when there are data-dependent transformations in the formula such as ns(), poly() or pspline(). That is, the set of basis functions chosen by pspline(x) depends on the range of x; for a new data prediction the basis functions are re-calculated, giving results that are wrong (unless the new x happens to have the exact same lower and upper limits). This is on my to-be-fixed list. Once a few other things are cleared away (actually a long list of other things). These fixes have been applied to lm and others for long enough now, coxph needs to catch up. Terry Therneau __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to check linearity in Cox regression
Use pspline within a Cox model. It includes a fairly general test for nonlinearity, that is similar to GAM models. Terry Therneau > coxph(Surv(time, status) ~ ph.ecog + pspline(age), lung) Call: coxph(formula = Surv(time, status) ~ ph.ecog + pspline(age), data = lung) coef se(coef) se2 Chisq DF p ph.ecog 0.4505 0.11766 0.11723 14.66 1.00 0.00013 pspline(age), linear 0.0112 0.00927 0.00927 1.45 1.00 0.23000 pspline(age), nonlin 2.96 3.08 0.41000 Iterations: 4 outer, 10 Newton-Raphson Theta= 0.797 Degrees of freedom for terms= 1.0 4.1 Likelihood ratio test=22.7 on 5.07 df, p=0.000412 n=227 (1 observation deleted due to missingness) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [R-sig-ME] lme nesting/interaction advice
On Mon, May 12, 2008 at 11:22 AM, Federico Calboli <[EMAIL PROTECTED]> wrote: > On 12 May 2008, at 17:09, Douglas Bates wrote: > >> I'm entering this discussion late so I may be discussing issues that >> have already been addressed. >> >> As I understand it, Federico, you began by describing a model for data >> in which two factors have a fixed set of levels and one factor has an >> extensible, or "random", set of levels and you wanted to fit a model >> that you described as >> >> y ~ effect1 * effect2 * effect3 >> >> The problem is that this specification is not complete. > > My apologies for that, I thought that the above formula was the shorthand > for what I would call the 'full' model, i.e. the single factors and the 2 > and 3 ways interactions. As I indicated, the trick is that the interaction of a fixed factor and a random factor can be defined in more than one way. It sounds as if what you want is lmer(y ~ factor1 * factor2 + (1|factor3) + (1|factor1:factor3) + (1|factor2:factor3) + (1|factor1:factor2:factor3), ...) but I'm not sure. >> An interaction of factors with fixed levels and a factor with random >> levels can mean, in the lmer specification, >> >> lmer(y ~ effect1 * effect2 + (1| effect3) + (1|effect1:effect2:effect3), >> ...) >> >> or >> >> lmer(y ~ effect1 * effect2 + (effect1*effect2 | effect3), ...) >> >> or other variations. When you specify a random effect or an random >> interaction term you must, either explicitly or implicitly, specify >> the form of the variance-covariance matrix associated with those >> random effects. > > I'll play around with this and see what I can get. >> >> The "advantage" that other software may provide for you is that it >> chooses the model for you but that, of course, means that you only >> have the one choice. > > I'm more than happy to stick to R, and to put more legwork into my models >> >> If you can describe how many variance components you think should be >> estimated in your model and what they would represent then I think it >> will be easier to describe how to fit the model. > > I'll work on that. Incidentally, what/where is the most comprehensive and up > to date documentation for lme4? the pdfs coming with the package? I suspect > knowing which are the right docs will help a lot in keeping me within the > boundaries of civility and prevent me from annoying anyone (which is not > something I sent forth to do on purpose). Documentation for lme4 is pretty sketchy at present. I hope to remedy that during our summer break. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Monty Hall simulation
The common "Monty Hall" problem (where the MC helps out the contestant) is not random! See: http://en.wikipedia.org/wiki/Monty_Hall_problem My edits were removed as I had no references. :-( Maybe you can verify my statements, included below. :-) Another analysis considers three types of hosts and three prize levels. The Benevolent Host always shows the worst remaining prize after you choose, the Random Host randomly picks a remaining door to show, and the Malevolent Host always shows the best remaining prize. The prizes are bad, middle, best; e.g. Goat, Luggage, Car. The Player is unaware of which prize is which. He may expect to be choosing among Pigs, Goats, Blenders, Luggage, Cars, and Houses. If you always switch, the results for each host are: Results after switching, expanded behaviors Host V /Prize ->Bad Middle Good Benevolent Host 0% 33% 67% Random Host 33% 33% 33% Malevolent Host 67% 33% 0% If each host were equally likely, the total probability for each prize would be 33%-the same as not switching. Without knowing the type of Host and the prize mix, you can make no meaningful statement about the success of a switching strategy. The only way for a Player to "improve the odds" is if he or she can get some meaningful information from the prize shown. I.e., if you know what the three prizes are and what type of Host you have, then you can develop a winning strategy. If you were wrong about either the Host or the prize mix, that strategy may be harmful. Robert Farley Metro www.Metro.net -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Johannes Huesing Sent: Sunday, May 11, 2008 01:26 To: r-help@r-project.org Subject: Re: [R] Monty Hall simulation cirrus74 <[EMAIL PROTECTED]> [Sun, May 11, 2008 at 03:44:46AM CEST]: > > Is it possible to simulate the Monty Hall problem using R? If so, could > someone please show me how? Thanks for any help rendered. The kind of simulation, as any thinking about this seeingly paradoxical situation, depends on your mindset. To my mind, niter <- 999 prize <- sample(c("car", "car", "goat"), niter, replace=TRUE) would be a perfect simulation. -- Johannes Hüsing There is something fascinating about science. One gets such wholesale returns of conjecture mailto:[EMAIL PROTECTED] from such a trifling investment of fact. http://derwisch.wikidot.com (Mark Twain, "Life on the Mississippi") __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] hessian in constrained optimization (constrOptim)
Dear helpers, I am using the function "constrOptim" to estimate a model with ML with an inequality constraint using the option method='Nelder-Mead'. When I specify the option: hessian = TRUE I obtain the response: Error in f(theta, ...) : unused argument(s) (hessian = TRUE) I guess the function "constrOptim" does not allow this argument which, on the other hand, is allowed in "optim". I would be extremely grateful if anybody could suggest a way I could use to I obtain the values of the hessian matrix... Many thanks, Carlo ** Carlo Fezzi Senior Research Associate Centre for Social Research on the Global Environment (CSERGE) Department of Environmental Sciences University of East Anglia Norwich (UK) NR2 7TJ Telephone: +44(0)1603 591408 Fax: +44(0)1603 593739 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Compact Indicator Matrices
On Mon, May 12, 2008 at 11:27 AM, amarkos <[EMAIL PROTECTED]> wrote: > Thanks, it works! > Could you please provide the direct method you mentioned for the > multivariate case? I'm not sure what you mean. I looked at what I wrote and I don't see anything that would fit that description. May I suggest that you continue to cc: the R-help list on the discussion. I can't always respond rapidly to requests and there are many who read the list that can. > On May 12, 4:30 pm, "Douglas Bates" <[EMAIL PROTECTED]> wrote: >> On Sun, May 11, 2008 at 9:49 AM, amarkos <[EMAIL PROTECTED]> wrote: >> > On May 11, 4:47 pm, "Douglas Bates" <[EMAIL PROTECTED]> wrote: >> >> >> Do you mean that you want to collapse similar rows into a single row >> >> and perhaps a count of the number of times that this row occurs? >> >> > Let me rephrase the problem by providing an example. >> >> > Input: >> >> > A = >> > [,1] [,2] >> > [1,]11 >> > [2,]13 >> > [3,]21 >> > [4,]12 >> > [5,]21 >> > [6,]12 >> > [7,]11 >> > [8,]12 >> > [9,]13 >> > [10,]21 >> >> An important question here is do you start with two or more variables >> like the columns of your matrix A? If so, there is a more direct >> method of getting the answers that you want. The natural way to store >> such variables in R is as factors. I prefer to use letters instead of >> numbers to represent the levels of a factor (that way I don't confuse >> a factor with a numeric variable when I look at rows) so I would >> create a data frame with two factors instead of a matrix. >> >> > V1 <- factor(c(1,1,2,1,2,1,1,1,1,2), labels = LETTERS[1:2]) >> > V2 <- factor(c(1,3,1,2,1,2,1,2,3,1), labels = letters[1:3]) >> > df <- data.frame(f1 = V1, f2 = V2) >> > df >> >>f1 f2 >> 1 A a >> 2 A c >> 3 B a >> 4 A b >> 5 B a >> 6 A b >> 7 A a >> 8 A b >> 9 A c >> 10 B a >> >> You could produce the indicator matrix and check for unique rows, etc. >> - I will show that below - but all you need is the interaction of the >> two factors >> >> > df$f12 <- with(df, f1:f2)[drop = TRUE] >> > df >> >>f1 f2 f12 >> 1 A a A:a >> 2 A c A:c >> 3 B a B:a >> 4 A b A:b >> 5 B a B:a >> 6 A b A:b >> 7 A a A:a >> 8 A b A:b >> 9 A c A:c >> 10 B a B:a> str(df) >> >> 'data.frame': 10 obs. of 3 variables: >> $ f1 : Factor w/ 2 levels "A","B": 1 1 2 1 2 1 1 1 1 2 >> $ f2 : Factor w/ 3 levels "a","b","c": 1 3 1 2 1 2 1 2 3 1 >> $ f12: Factor w/ 4 levels "A:a","A:b","A:c",..: 1 3 4 2 4 2 1 2 3 4 >> >> > table(df$f12) >> >> A:a A:b A:c B:a >> 2 3 2 3> as.numeric(df$f12) >> >> [1] 1 3 4 2 4 2 1 2 3 4 >> >> Notice that this shows you that there are four distinct combinations >> that occur 2, 3, 2 and 3 times respectively; the first combination >> occurs in rows 1 and 7, it consists of the first level of f1 and the >> first level of f2, etc. >> >> If you really do want the indicator matrix you could generate it as >> >> > (ind <- cbind(model.matrix(~ 0 + f1, df), model.matrix(~ 0 + f2, df))) >> >>f1A f1B f2a f2b f2c >> 11 0 1 0 0 >> 21 0 0 0 1 >> 30 1 1 0 0 >> 41 0 0 1 0 >> 50 1 1 0 0 >> 61 0 0 1 0 >> 71 0 1 0 0 >> 81 0 0 1 0 >> 91 0 0 0 1 >> 10 0 1 1 0 0> unique(ind) >> >> f1A f1B f2a f2b f2c >> 1 1 0 1 0 0 >> 2 1 0 0 0 1 >> 3 0 1 1 0 0 >> 4 1 0 0 1 0 >> >> but working with the factors is generally much simpler than working >> with the indicators. >> >> >> >> > # Indicator matrix >> > A <- data.frame(lapply(data.frame(obj), as.factor)) >> >> > nocases <- dim(obj)[1] >> > novars <- dim(obj)[2] >> >> > # variable levels >> > levels.n <- sapply(obj, nlevels) >> > n<- cumsum(levels.n) >> >> > # Indicator matrix calculations >> > Z<- matrix(0, nrow = nocases, ncol = n[length(n)]) >> > newdat <- lapply(obj, as.numeric) >> > offset <- (c(0, n[-length(n)])) >> > for (i in 1:novars) >> > Z[1:nocases + (nocases * (offset[i] + newdat[[i]] - 1))] <- 1 >> >> > ### >> >> > Output: >> >> > Z = >> >> >[,1] [,2] [,3] [,4] [,5] >> > [1,]10100 >> > [2,]10001 >> > [3,]01100 >> > [4,]10010 >> > [5,]01100 >> > [6,]10010 >> > [7,]10100 >> > [8,]10010 >> > [9,]10001 >> > [10,]01100 >> >> > Z is an indicator matrix in the Multiple Correspondence Analysis >> > framework. >> > My problem is to collapse identical rows (e.g. 2 and 9) into a single >> > row and >> > store the row ids. >> >> > __ >> > [EMAIL PROTECTED] mailing list >> >https://stat.ethz.ch/mailman/listinfo/r-help >> > PLEASE do read the posting guidehttp://www.R-project.org/pos
Re: [R] how to not sort factors when plotting
On Mon, 2008-05-12 at 11:41 -0400, Jorge Ivan Velez wrote: > Hi Lydia, [I'm struggling to see what this has to do with the subject line?] > > I compared my ratio function with Dimitris and Phil's suggestions. Please do > NOT use my approach because it's painfully slow for a large vector (as Phil > told me). Here is why (using Win XP SP2, Intel Core- 2 Duo 2.4 GHz, R 2.7.0 > Patched): Jorge, If you pre-allocate storage for temp rather than the concatenation approach you use, this is reasonably speedy for the example you quote: Note: you have 1:length(x) which throws and error, the loop condition needs to be 2:length(x) (or more robustly: seq_along(x[-1]) + 1) > x <- rnorm(10,0,1) > my.ratio <- function(x){ + temp <- numeric(length(x)) + for(n in 2:length(x)) + temp[n] <- x[n] / x[n-1] + temp + } > system.time(my.ratio(x)) user system elapsed 0.757 0.000 0.758 > new.ratio <- function(x) x[2:length(x)]/x[1:(length(x)-1)] > system.time(new.ratio(x)) user system elapsed 0.011 0.003 0.013 OK, so it isn't faster than the vectorised approach, but it isn't bad. For those more familiar with C-type programming than the R vector approach, you can do reasonably well with a for loop as long as you do proper allocation of the result vector/object first. Note that I'm not advocating that people shouldn't bother to learn to use R to its advantages and code for R rather than as you might have learnt from other languages (I'm not!), but for some problems a loop is just fine unless you are the sort of person who needs that extra 0.7 of a second to do something else with... ;-) HTH G > > > # Vector > x=rnorm(10,0,1) > > # Suggestion > new.ratio=function(x) x[2:length(x)]/x[1:(length(x)-1)] > > # My horrible function > my.ratio=function(x){ > temp=NULL > for (n in 1:length(x)) temp=c(temp,x[n]/x[n-1]) > temp > } > > # System time > t=system.time(my.ratio(x)) > tnr=system.time(new.ratio(x)) > t >user system elapsed > 38.790.06 39.31 > tnr >user system elapsed > 0 0 0 > > > Thanks to all, > > Jorge > > > > On Mon, May 12, 2008 at 11:15 AM, Phil Spector <[EMAIL PROTECTED]> > wrote: > > > Another alternative would be to take advantage of R's vectorization: > > > > x=c(1,2,3,2,1,2,3) > > > x[2:length(x)]/x[1:(length(x)-1)] > > > > > [1] 2.000 1.500 0.667 0.500 2.000 1.500 > > > > The solution using your ratio function will be painfully slow > > for a large vector. > > > > - Phil Spector > > Statistical Computing Facility > > Department of Statistics > > UC Berkeley > > [EMAIL PROTECTED] > > > > > > [[alternative HTML version deleted]] > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. -- %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~% Dr. Gavin Simpson [t] +44 (0)20 7679 0522 ECRC, UCL Geography, [f] +44 (0)20 7679 0565 Pearson Building, [e] gavin.simpsonATNOSPAMucl.ac.uk Gower Street, London [w] http://www.ucl.ac.uk/~ucfagls/ UK. WC1E 6BT. [w] http://www.freshwaters.org.uk %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~% __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Fonts in Quartz Devices
El 12/05/2008, a las 18:18, Prof Brian Ripley escribió: On Mon, 12 May 2008, Arcadio Rubio García wrote: Hi, I'm new to R. I'm using a Mac OS X 10.5 and R 2.7.0. I'm trying to change the font family for a plot in a quartz device. Simply passing the desired font by using the family argument works with other devices, but not with quartz. Am I missing anything? I've already checked the docs for quartz and quartzFonts, but didn't find a solution. Can you tell us what you have tried? Some families worked not long before release, at least. I tried quartz(family="Monaco") and quartz(family="Courier") with no success. I've just downgraded back to R 2.6.2 and both commands work fine. Such MacOS-specific questions are more likely to get an informed answer on R-sig-mac. Any help much appreciated. Regards, A. Rubio -- 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@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [R-sig-ME] lme nesting/interaction advice
On 12 May 2008, at 17:09, Douglas Bates wrote: I'm entering this discussion late so I may be discussing issues that have already been addressed. As I understand it, Federico, you began by describing a model for data in which two factors have a fixed set of levels and one factor has an extensible, or "random", set of levels and you wanted to fit a model that you described as y ~ effect1 * effect2 * effect3 The problem is that this specification is not complete. My apologies for that, I thought that the above formula was the shorthand for what I would call the 'full' model, i.e. the single factors and the 2 and 3 ways interactions. An interaction of factors with fixed levels and a factor with random levels can mean, in the lmer specification, lmer(y ~ effect1 * effect2 + (1| effect3) + (1| effect1:effect2:effect3), ...) or lmer(y ~ effect1 * effect2 + (effect1*effect2 | effect3), ...) or other variations. When you specify a random effect or an random interaction term you must, either explicitly or implicitly, specify the form of the variance-covariance matrix associated with those random effects. I'll play around with this and see what I can get. The "advantage" that other software may provide for you is that it chooses the model for you but that, of course, means that you only have the one choice. I'm more than happy to stick to R, and to put more legwork into my models If you can describe how many variance components you think should be estimated in your model and what they would represent then I think it will be easier to describe how to fit the model. I'll work on that. Incidentally, what/where is the most comprehensive and up to date documentation for lme4? the pdfs coming with the package? I suspect knowing which are the right docs will help a lot in keeping me within the boundaries of civility and prevent me from annoying anyone (which is not something I sent forth to do on purpose). Best regards, Federico -- Federico C. F. Calboli Department of Epidemiology and Public Health Imperial College, St. Mary's Campus Norfolk Place, London W2 1PG Tel +44 (0)20 75941602 Fax +44 (0)20 75943193 f.calboli [.a.t] imperial.ac.uk f.calboli [.a.t] gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Fonts in Quartz Devices
On Mon, 12 May 2008, Arcadio Rubio García wrote: Hi, I'm new to R. I'm using a Mac OS X 10.5 and R 2.7.0. I'm trying to change the font family for a plot in a quartz device. Simply passing the desired font by using the family argument works with other devices, but not with quartz. Am I missing anything? I've already checked the docs for quartz and quartzFonts, but didn't find a solution. Can you tell us what you have tried? Some families worked not long before release, at least. Such MacOS-specific questions are more likely to get an informed answer on R-sig-mac. Any help much appreciated. Regards, A. Rubio -- 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@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Mathematical annotation in lattice strip: Is it possible?
On 5/12/08, Andrewjohnclose <[EMAIL PROTECTED]> wrote: > > I have tried without success to find a way including the square root symbol > in lattice strips as part of my conditioning labels. I have tried > supplementing by creating a list of vectors using the var.name function > coupled with the expression function used in xlab/ylab. > > xyplot(adjusted_Rand_index~cluster|distance_measure, main="Level of > agreement between partitions: Wards Method", ylab="Coefficient value > (adjusted rand index)", xlab="number of clusters", type="l", data=randA1, > strip=strip.custom(varnames=c(expression(sqrt(Bray-Curtis))) You say 'var.names', but use 'varnames'. ?strip.custom also says: var.name: vector of character strings or expressions as long as the number of conditioning variables. The contents are interpreted as names for the conditioning variables. Whether they are shown on the strip depends on the values of 'strip.names' and 'style' (see below). By default, the names are shown for shingles, but not for factors. This seems to work as expected: xyplot(1 ~ 1 | gl(1, 1), strip = strip.custom(var.name = expression(sqrt(Bray-Curtis)), strip.names = TRUE, strip.levels = FALSE)) -Deepayan __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [R-sig-ME] lme nesting/interaction advice
I'm entering this discussion late so I may be discussing issues that have already been addressed. As I understand it, Federico, you began by describing a model for data in which two factors have a fixed set of levels and one factor has an extensible, or "random", set of levels and you wanted to fit a model that you described as y ~ effect1 * effect2 * effect3 The problem is that this specification is not complete. An interaction of factors with fixed levels and a factor with random levels can mean, in the lmer specification, lmer(y ~ effect1 * effect2 + (1| effect3) + (1|effect1:effect2:effect3), ...) or lmer(y ~ effect1 * effect2 + (effect1*effect2 | effect3), ...) or other variations. When you specify a random effect or an random interaction term you must, either explicitly or implicitly, specify the form of the variance-covariance matrix associated with those random effects. The "advantage" that other software may provide for you is that it chooses the model for you but that, of course, means that you only have the one choice. If you can describe how many variance components you think should be estimated in your model and what they would represent then I think it will be easier to describe how to fit the model. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Fundamental formula and dataframe question.
I would have thought that: > lm( C1 ~ M^2, data=DF ) Would give the main effects and 2 way interaction(s) (but a quick test did not match my expectation). Possibly a feature request is in order if people plan to use this a lot. -- Gregory (Greg) L. Snow Ph.D. Statistical Data Center Intermountain Healthcare [EMAIL PROTECTED] (801) 408-8111 > -Original Message- > From: [EMAIL PROTECTED] > [mailto:[EMAIL PROTECTED] On Behalf Of Ted Harding > Sent: Sunday, May 11, 2008 2:07 PM > To: Myers, Brent > Cc: r-help@r-project.org > Subject: Re: [R] Fundamental formula and dataframe question. > > On 11-May-08 18:58:45, Myers, Brent wrote: > > There is a very useful and apparently fundamental feature > of R (or of > > the package pls) which I don't understand. > > > > For datasets with many independent (X) variables such as > chemometric > > datasets there is a convenient formula and dataframe > construction that > > allows one to access the entire X matrix with a single term. > > > > Consider the gasoline dataset available in the pls package. For the > > model statement in the plsr function one can write: Octane ~ NIR > > > > NIR refers to a (wide) matrix which is a portion of a > dataframe. The > > naming of the columns is of the form: 'NIR. nm' > > > > names(gasoline) returns... > > > > $names > > [1] "octane" "NIR" > > > > instead of... > > > > $names > > [1] "octane" "NIR.1000 nm" "NIR.1001 nm" ... > > > > How do I construct and manipulate such dataframes and the > column names > > that go with? > > > > Does the use of these types of formulas and dataframes > generalize to > > other modeling functions? > > > > Some specific clues on a help search might be enough, I've > tried many. > > > > Regards, > > Brent > > I don't have the 'gasoline' dataset to hand, but I can > produce something to which your descrption applies as follows: > > C1 <- c(1.1,1.2,1.3,1.4) > C2 <- c(2.1,2.2,2.3,2.4) >M <- cbind(M1=c(11.1,11.2,11.3,11.4), > M2=c(12.1,12.2,12.3,12.4)) > DF <- data.frame(C1=C1,C2=C2,M=M) > DF > #C1 C2 M.M1 M.M2 > # 1 1.1 2.1 11.1 12.1 > # 2 1.2 2.2 11.2 12.2 > # 3 1.3 2.3 11.3 12.3 > # 4 1.4 2.4 11.4 12.4 > > so the two columns C1 and C2 have gone in as named, and the > matrix M (with named columns M1 and M2) has gone in with > columns M.M1, M.M2 > > Now let's fuzz the numbers a bit, so that the lm() fit makes sense: > > C1 <- C1 + round(0.1*runif(4),2) > C1 <- C1 + round(0.1*runif(4),2) >M <- cbind(M1=c(11.1,11.2,11.3,11.4), > M2=c(12.1,12.2,12.3,12.4)) + > round(0.1*runif(8),2) > DF <- data.frame(C1=C1,C2=C2,M=M) > DF > # C1 C2 M.M1 M.M2 > # 1 1.21 2.1 11.19 12.13 > # 2 1.34 2.2 11.23 12.23 > # 3 1.38 2.3 11.36 12.30 > # 4 1.50 2.4 11.43 12.48 > > summary(lm(C1 ~ M),data=DF) > # Call: > # lm(formula = C1 ~ M) > # Residuals: > #1234 > # -0.02422 0.02448 0.01309 -0.01335 > # Coefficients: > # Estimate Std. Error t value Pr(>|t|) > # (Intercept) -8.284352.48952 -3.3280.186 > # MM1 -0.054110.66909 -0.0810.949 > # MM2 0.834630.50687 1.6470.347 > # Residual standard error: 0.03919 on 1 degrees of freedom > # Multiple R-Squared: 0.9642, Adjusted R-squared: 0.8925 > # F-statistic: 13.46 on 2 and 1 DF, p-value: 0.1893 > > In other words, a perfectly standard LM fit, equivalent to > > summary(lm(C1 ~ M[,1]+M[,2])) > > (as you can check). So all that looks straightforward. > > One thing, however, is not clear to me in this scenario. > Suppose, for example, that the columns M1 and M2 of M were > factors (and that you had more rows than I've used above, so > that the fit is non-trivial). > > Then, in the standard specification of an LM, you could write > > summary(lm(C1 ~ M[,1]*M[,2])) > > and get the main effects and interactions. But how would you > do that in the other type of specification: > > Where you used > summary(lm(C1 ~ M, data=DF)) > to get the equivalent of > summary(lm(C1 ~ M[,1]+M[,2])) > what would you use to get the equivalent of > summary(lm(C1 ~ M[,1]*M[,2]))?? > > Would you have to "spell out" the interaction term[s] in > additional columns of M? > > Hmmm, interesting! I hadn't been aware of this aspect of > formula and dataframe construction for modellinng, until you > pointed it out! > > Best wishes, > Ted. > > > E-Mail: (Ted Harding) <[EMAIL PROTECTED]> > Fax-to-email: +44 (0)870 094 0861 > Date: 11-May-08 Time: 21:06:49 > -- XFMail -- > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. >
[R] Left censored responses in mixed effects models
Dear R Fellow-Travellers: What is your recommended way of dealing with a left-censored response (non-detects) in (linear Gaussian) mixed effects models? Specifics: Response is a numeric positive measurement (of volume, actually); but when it falls below some unknown and slightly random value (depending on how the sample is prepared and measured), it cannot be measured and is recorded as 0. There is some statistical literature on this, but I was unable to find anything that appeared to me to implement a strategy in any R package. If it matters, I am less interested in inference than in removing possible bias in estimation. Feel free to respond off-list if you feel that this would not be of general interest. Cheers, Bert Gunter Genentech __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Random number generation
> -Original Message- > From: [EMAIL PROTECTED] > [mailto:[EMAIL PROTECTED] On Behalf Of Esmail Bonakdarian > Sent: Sunday, May 11, 2008 7:25 AM > To: Prof Brian Ripley > Cc: [EMAIL PROTECTED] > Subject: Re: [R] Random number generation [snip] > What I read doesn't seem to be incorrect however (it may even > have been an archived message here), the *language* itself > does not seem to support block *comments*. Using conditional > constructs, or an IDE/editor to achieve similar results is a > work around - but not the same. I don't mean to nitpick, but > as a computer scientist I see this as different :-) I am not a computer scientist, so correct me if I am wrong, but from what I remember (and a quick glance at my copy of Kernighan Ritchie), the C *language* itself does not support block *comments*, rather the preproccessor replaces the comments with a single space and the compilor does not even see them. Since R is optimized for interactive use rather than compilation, running everything through a preproccessor is not really an option. However as an additional work around you could always run your R scripts through the C preproccessor and have it strip the block comments for you. Given the complexity of implementing block commenting (even deciding on the syntax) and the ease of current work arounds, the cost benefit ratio probably puts this very near the bottom of the priority list. -- Gregory (Greg) L. Snow Ph.D. Statistical Data Center Intermountain Healthcare [EMAIL PROTECTED] (801) 408-8111 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] inserting mathematical symbols in lattice strip
> I have tried without success to find a way including the square root symbol > in lattice strips as part of my conditioning labels. I have tried > supplementing by creating a list of vectors using the var.name function > coupled with the expression function used in xlab/ylab. > > xyplot(adjusted_Rand_index~cluster|distance_measure, main="Level of > agreement between partitions: Wards Method", ylab="Coefficient value > (adjusted rand index)", xlab="number of clusters", type="l", data=randA1, > strip=strip.custom(varnames=c(expression(sqrt(Bray-Curtis))) > > Is there a way of generating the square root symbol inside the strip or am I > wasting my time. i.e. converting Bray-Curtis to ... [sqrt symbol] > Bray-Curtis. You want factor.level, not var.name in strip.custom. Here's an example with the iris data. Apologies for the second line; there's almost certainly an easier way to do it, but you (hopefully) get the idea. flevels <- levels(iris$Species) foo <- paste("c(", paste("expression(paste(", flevels, ", sqrt(", 1:3, ")))", sep="", collapse=","), ")") xyplot(Sepal.Length ~ Petal.Length | Species, data = iris, strip=strip.custom(factor.levels=eval(parse(text=foo Regards, Richie. Mathematical Sciences Unit HSL ATTENTION: This message contains privileged and confidential inform...{{dropped:20}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to not sort factors when plotting
Hi Lydia, I compared my ratio function with Dimitris and Phil's suggestions. Please do NOT use my approach because it's painfully slow for a large vector (as Phil told me). Here is why (using Win XP SP2, Intel Core- 2 Duo 2.4 GHz, R 2.7.0 Patched): # Vector x=rnorm(10,0,1) # Suggestion new.ratio=function(x) x[2:length(x)]/x[1:(length(x)-1)] # My horrible function my.ratio=function(x){ temp=NULL for (n in 1:length(x)) temp=c(temp,x[n]/x[n-1]) temp } # System time t=system.time(my.ratio(x)) tnr=system.time(new.ratio(x)) t user system elapsed 38.790.06 39.31 tnr user system elapsed 0 0 0 Thanks to all, Jorge On Mon, May 12, 2008 at 11:15 AM, Phil Spector <[EMAIL PROTECTED]> wrote: > Another alternative would be to take advantage of R's vectorization: > > x=c(1,2,3,2,1,2,3) > > x[2:length(x)]/x[1:(length(x)-1)] > > > [1] 2.000 1.500 0.667 0.500 2.000 1.500 > > The solution using your ratio function will be painfully slow > for a large vector. > > - Phil Spector > Statistical Computing Facility > Department of Statistics > UC Berkeley > [EMAIL PROTECTED] > > [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] ggplot2: font size mismatch for pdf output
> More generally, how do I control the size of fonts used in legends > and axis labels? There is no general way (yet) - it is on my customisation to do list, which I hope to make progress on over summer. Hadley -- http://had.co.nz/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to not sort factors when plotting
this can be more efficiently coded using indexing, e.g., x <- c(1,2,3,2,1,2,3) n <- length(x) x[2:n] / x[1:(n-1)] Best, Dimitris Dimitris Rizopoulos Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven, Belgium Tel: +32/(0)16/336899 Fax: +32/(0)16/337015 Web: http://med.kuleuven.be/biostat/ http://www.student.kuleuven.be/~m0390867/dimitris.htm Quoting Jorge Ivan Velez <[EMAIL PROTECTED]>: Hi Lydia, Try this: # Function ratio=function(x){ temp=NULL for (n in 1:length(x)) temp=c(temp,x[n]/x[n-1]) temp } # Example x=c(1,2,3,2,1,2,3) ratio(x) [1] 2.000 1.500 0.667 0.500 2.000 1.500 HTH, Jorge On Mon, May 12, 2008 at 9:52 AM, Lydia N. Slobodian <[EMAIL PROTECTED]> wrote: Hello. I'm trying find the ratios between each of the integers in a vector. I have: for (n in x) { ratio <- (x[n]/x[n-1]) ratio.all <- c(ratio.all, ratio) } Of course this doesn't work, nor does diff(n)/diff(n-1). Is there a way to specify a pair of integers in a vector? Thank you, Lydia __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Fonts in Quartz Devices
In relation to my previous email, I've discovered that quartz(family = "Monaco") works fine with 2.6.2. However, it doesn't with 2.7.0. Moreover, with 2.6.2 I get lot's of warnings, that previously I didn't. Maybe an update of the OS has broken some code leading to a bug? El 12/05/2008, a las 16:43, Arcadio Rubio García escribió: Hi, I'm new to R. I'm using a Mac OS X 10.5 and R 2.7.0. I'm trying to change the font family for a plot in a quartz device. Simply passing the desired font by using the family argument works with other devices, but not with quartz. Am I missing anything? I've already checked the docs for quartz and quartzFonts, but didn't find a solution. Any help much appreciated. Regards, A. Rubio __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to not sort factors when plotting
Hi Lydia, Try this: # Function ratio=function(x){ temp=NULL for (n in 1:length(x)) temp=c(temp,x[n]/x[n-1]) temp } # Example x=c(1,2,3,2,1,2,3) ratio(x) [1] 2.000 1.500 0.667 0.500 2.000 1.500 HTH, Jorge On Mon, May 12, 2008 at 9:52 AM, Lydia N. Slobodian <[EMAIL PROTECTED]> wrote: > Hello. I'm trying find the ratios between each of the integers in a > vector. I have: > > for (n in x) { > ratio <- (x[n]/x[n-1]) > ratio.all <- c(ratio.all, ratio) > } > > Of course this doesn't work, nor does diff(n)/diff(n-1). Is there a > way to specify a pair of integers in a vector? > > Thank you, > > Lydia > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] ggplot2: font size mismatch for pdf output
Sorry for not providing this in my initial posting. I'm using R version 2.7.0 (2008-04-22), on Win XP Pro. As I said, the .png output matched what I saw on screen. It was the .pdf output for which the font size was noticeably larger, enough to make the legend run off the screen. -Michael Prof Brian Ripley wrote: What version of R and OS is this? Prior to R 2.7.0 there was little attempt to match output dimensions from various devices, and one of the png devices in 2.7.0 has an error in doing so, fixed in R-patched (see NEWS). On Mon, 12 May 2008, Michael Friendly wrote: Hi In the following, the graph I see on the screen and the .png output coincide. However, in the .pdf file, the fonts seem to be scaled fairly larger, resulting in the label for the top legend disappearing. Is this an infelicity or bug, or is there something I've missed? More generally, how do I control the size of fonts used in legends and axis labels? library(car) library(ggplot2) qp <-qplot (education , income , shape=type , size=women , colour=prestige , xlab="Education" , ylab="Income", data=Prestige) + scale_y_continuous(limits=c(NA, 2)) Hmm, you can't break the line before '+'. qp + scale_size(to=c(1,8)) ggsave(file="prestige-ggplot.png", width=6, height=5) # OK ggsave(file="prestige-ggplot.pdf", width=6, height=5) # fonts too large I would not expect you to be able to specify a smaller size without also reducing 'pointsize'. E.g. dev.print() does so, but ggsave seems not to. -Michael __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Michael Friendly Email: friendly AT yorku DOT ca Professor, Psychology Dept. York University Voice: 416 736-5115 x66249 Fax: 416 736-5814 4700 Keele Streethttp://www.math.yorku.ca/SCS/friendly.html Toronto, ONT M3J 1P3 CANADA __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] help with calculating the differences between dates
Tom Cohen wrote: Dear list, How can I calculate the difference in days between the eventdate and basedate in the below dataset? id basedate outcome.3 eventdate daydiff 1 1001 1999-09-28 2 1999-10-013 2 1002 1999-09-22 1 3 1003 2000-01-19 1 4 1004 2004-01-25 2 2004-02-039 5 1005 2005-08-11 1 6 1006 2000-07-04 1 2001-05-29 7 1007 2004-02-12 1 2004-11-18 8 1008 2006-01-18 2 2006-02-02 9 1009 2005-04-29 2 2005-06-14 10 1010 2006-03-17 2 2006-03-31 11 1011 2000-03-21 2 2000-03-28 12 1012 2004-07-12 1 2006-11-28 13 1013 2000-02-24 1 14 1014 2003-04-17 1 15 1015 2000-04-05 1 Given they are in some appropriate time format (e.g. when reading specify colClasses = "POSIXct" for those columns or use strptime() on the colummns later on): difftime(x$eventdate, x$basedate, "days") Uwe Ligges Thanks for any help, Tom - Går det långsamt? Skaffa dig en snabbare bredbandsuppkoppling. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Fonts in Quartz Devices
Hi, I'm new to R. I'm using a Mac OS X 10.5 and R 2.7.0. I'm trying to change the font family for a plot in a quartz device. Simply passing the desired font by using the family argument works with other devices, but not with quartz. Am I missing anything? I've already checked the docs for quartz and quartzFonts, but didn't find a solution. Any help much appreciated. Regards, A. Rubio __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] ggplot2: font size mismatch for pdf output
What version of R and OS is this? Prior to R 2.7.0 there was little attempt to match output dimensions from various devices, and one of the png devices in 2.7.0 has an error in doing so, fixed in R-patched (see NEWS). On Mon, 12 May 2008, Michael Friendly wrote: Hi In the following, the graph I see on the screen and the .png output coincide. However, in the .pdf file, the fonts seem to be scaled fairly larger, resulting in the label for the top legend disappearing. Is this an infelicity or bug, or is there something I've missed? More generally, how do I control the size of fonts used in legends and axis labels? library(car) library(ggplot2) qp <-qplot (education , income , shape=type , size=women , colour=prestige , xlab="Education" , ylab="Income", data=Prestige) + scale_y_continuous(limits=c(NA, 2)) Hmm, you can't break the line before '+'. qp + scale_size(to=c(1,8)) ggsave(file="prestige-ggplot.png", width=6, height=5) # OK ggsave(file="prestige-ggplot.pdf", width=6, height=5) # fonts too large I would not expect you to be able to specify a smaller size without also reducing 'pointsize'. E.g. dev.print() does so, but ggsave seems not to. -Michael __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- 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@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] inserting mathematical symbols in lattice strip
I have tried without success to find a way including the square root symbol in lattice strips as part of my conditioning labels. I have tried supplementing by creating a list of vectors using the var.name function coupled with the expression function used in xlab/ylab. xyplot(adjusted_Rand_index~cluster|distance_measure, main="Level of agreement between partitions: Wards Method", ylab="Coefficient value (adjusted rand index)", xlab="number of clusters", type="l", data=randA1, strip=strip.custom(varnames=c(expression(sqrt(Bray-Curtis))) Is there a way of generating the square root symbol inside the strip or am I wasting my time. i.e. converting Bray-Curtis to ... [sqrt symbol] Bray-Curtis. Thank you very much Regards Andrew http://www.nabble.com/file/p17188291/randA1.csv randA1.csv -- View this message in context: http://www.nabble.com/inserting-mathematical-symbols-in-lattice-strip-tp17188291p17188291.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] help with calculating the differences between dates
Dear list, How can I calculate the difference in days between the eventdate and basedate in the below dataset? id basedate outcome.3 eventdate daydiff 1 1001 1999-09-28 2 1999-10-013 2 1002 1999-09-22 1 3 1003 2000-01-19 1 4 1004 2004-01-25 2 2004-02-039 5 1005 2005-08-11 1 6 1006 2000-07-04 1 2001-05-29 7 1007 2004-02-12 1 2004-11-18 8 1008 2006-01-18 2 2006-02-02 9 1009 2005-04-29 2 2005-06-14 10 1010 2006-03-17 2 2006-03-31 11 1011 2000-03-21 2 2000-03-28 12 1012 2004-07-12 1 2006-11-28 13 1013 2000-02-24 1 14 1014 2003-04-17 1 15 1015 2000-04-05 1 Thanks for any help, Tom - Går det långsamt? Skaffa dig en snabbare bredbandsuppkoppling. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] RPM-style install (SLED 10.1)
I am trying to install R on a SLED 10.1 machine. R-base-2.7.0-7.1-i586.rpm fails with [EMAIL PROTECTED]:~/RPMs> rpm -Uvh R-base-2.7.0-7.1.i586.rpm warning: R-base-2.7.0-7.1.i586.rpm: Header V3 DSA signature: NOKEY, key ID 14ec5930 error: Failed dependencies: libgfortran.so.1 is needed by R-base-2.7.0-7.1.i586 I tried to trick it into believing there's the library by setting up the link: [EMAIL PROTECTED]:~/RPMs> ls -l /usr/lib/libgf* lrwxrwxrwx 1 root root 29 2008-05-08 17:45 /usr/lib/libgfortran.so.1 -> /usr/lib/libgfortran.so.3.0.0 lrwxrwxrwx 1 root root 20 2008-05-08 17:37 /usr/lib/libgfortran.so.3 -> libgfortran.so.3.0.0 -rwxr-xr-x 1 root root 2251139 2008-05-01 09:43 /usr/lib/libgfortran.so.3.0.0 but that did not work, either. I went on and installed GCC's gfortran, but it did not provide libgfortran.so.1 either: [EMAIL PROTECTED]:~/RPMs> ls -l /usr/irun/lib/libgf* -rw-r--r-- 1 1005 1011 4330512 2008-03-02 02:29 /usr/irun/lib/libgfortran.a -rwxr-xr-x 1 1005 10111009 2008-03-02 02:29 /usr/irun/lib/libgfortran.la -rwxr-xr-x 1 1005 1011 2509309 2008-03-02 02:29 /usr/irun/lib/libgfortran.so -rwxr-xr-x 1 1005 1011 2509309 2008-03-02 02:29 /usr/irun/lib/libgfortran.so.3 -rwxr-xr-x 1 1005 1011 2509309 2008-03-02 02:29 /usr/irun/lib/libgfortran.so.3.0.0 Any reason R wants that package instead of the newer one? Of course rpm --nodeps was an option, but then it fails to load r-stats with the same message about the missing library: Error in dyn.load(file, DLLpath = DLLpath, ...) : unable to load shared library '/usr/lib/R/library/stats/libs/stats.so': /usr/lib/R/library/stats/libs/stats.so: undefined symbol: _gfortran_pow_r8_i4 During startup - Warning message: package stats in options("defaultPackages") was not found I am stuck... I am trying to compile R from 2.6.2 sources in the meantime, but that's painful for a Linux newbie like me on a brand-new computer that does not necessarily have all the packages. I almost want to switch to Ubuntu as R installation back there was a single click :)) -- Stas Kolenikov, also found at http://stas.kolenikov.name Small print: Please do not reply to my Gmail address as I don't check it regularly. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] rpart with circular data?
On Mon, May 12, 2008 at 6:29 AM, Bálint Czúcz <[EMAIL PROTECTED]> wrote: > Hi Rainer, > > In a similar situation I used the two components of the normal vector > of the surface ("northing" & "easting"). I.e. for a horizontal plane > both are 0, for a vertical slope facing south northing=-1and > easting=0, etc. This descartian decomposition of the slope vector > avoids the problem of circularity present in the widespreadly used > polar (aspect, slope) decomposition, and thus seems to suit ecological > problems much better to me. However I have not looked into this much, > I am also very interested in the opinion of others. > > Bálint > > > On Mon, May 12, 2008 at 11:35 AM, Rainer M Krug <[EMAIL PROTECTED]> wrote: >> Hi >> >> I am planning to use a CART analysis with rpart() to analyse the >> impact of, among others, slope, altitude and aspect on mortality >> rates. >> My question is: >> Is there a p[roblem with using aspect as a predictor as it is circular? >> And if it is a problem (which I suspect), is there a transformation I >> could use to transform aspect? >> >> Thanks >> >> Rainer >> >> -- >> Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation >> Biology, UCT), Dipl. Phys. (Germany) >> >> Plant Conservation Unit >> Department of Botany >> University of Cape Town >> Rondebosch 7701 >> South Africa >> >> __ >> R-help@r-project.org mailing list >> https://stat.ethz.ch/mailman/listinfo/r-help >> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html >> and provide commented, minimal, self-contained, reproducible code. >> > > I would recommend going straight to the source-- modeling solar radiation. This gives a purely physical-based generalization of terrain-induced variation in microclimate. Keep an eye out for this paper, soon to be published in Soil. Sci, Soc, America: D.E. Beaudette, and A.T. O'Geen. Quantifying the Aspect Effect. Soil Sci. Soc. Am. J. March 3, 2008. (In Review) This paper uses a GRASS-based solution to the problem. Cheers Dylan __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [R-sig-ME] lme nesting/interaction advice
On 12 May 2008, at 14:37, Doran, Harold wrote: I haven't followed this thread carefully, so apologies if I'm too off base. But, in response to Rolf's questions/issues. First, SAS cannot handle models with crossed random effects (at least well at all). SAS is horribly incapable of handling even the simplest of models (especially generalized linear mixed models). I can cite numerous (recent) examples of SAS coming to a complete halt (proc nlmixed) for an analyses we were recently working on. R (and Ubuntu) was the only solution to our problem First off, let's keep SAS out of this. I never used it, never wanted to use it and did not mention anywhere I wanted to get SAS-like results! Although, seeing how easily it creeps up, I can sympathise with those who have strog feelings about it! [for those with strong feelings about me, this is meant to be something joke-like] Now, lme is not optimized for crossed random effects, but lmer is. That is why lmer is supported and lme is not really supported much. lmer is optimized for models with nested random effects and crossed random effects. When working with models with nested random effects, and software optimized for those problems (e.g., HLM, SAS, mlWin) the variance/covariance matrix forms a special, and simple structure that can be easily worked with. This is not the case for models with crossed random effects. Software packages designed for nested random effects can be tricked into handling models with crossed random effects, but this kludge is slow and really inefficient. If you want complete transparency into the why and how, here is a citation for your review. Thank you very much. I'll read the paper and hopefully get the answers I was looking for. Best, Federico Best Harold @article{Doran:Bates:Bliese:Dowling:2007:JSSOBK:v20i02, author = "Harold Doran and Douglas Bates and Paul Bliese and Maritza Dowling", title = "Estimating the Multilevel Rasch Model: With the lme4 Package", journal = "Journal of Statistical Software", volume = "20", number = "2", pages = "1--18", day = "22", month = "2", year ="2007", CODEN = "JSSOBK", ISSN ="1548-7660", bibdate = "2007-02-22", URL = "http://www.jstatsoft.org/v20/i02";, accepted ="2007-02-22", acknowledgement = "", keywords ="", submitted = "2006-10-01", } ___ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models -- Federico C. F. Calboli Department of Epidemiology and Public Health Imperial College, St. Mary's Campus Norfolk Place, London W2 1PG Tel +44 (0)20 75941602 Fax +44 (0)20 75943193 f.calboli [.a.t] imperial.ac.uk f.calboli [.a.t] gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [R-sig-ME] lme nesting/interaction advice
On 12 May 2008, at 12:21, Nick Isaac wrote: I *think* the syntax for the model Federico wants is this: lmer(y~selection*males+ (selection|month) + (males|month)) I'll try and check against some back of the envelope calculations -- as I said, the model is, per se, nothing really new, and my data is fully balanced. My lme syntax is a bit rusty, so I'm not confident how to recode with nested random effects, as in P&B p24. Two quick points: 1. I think Federico has caused some confusion on account of the way he used the term 'crossing' Sorry about that, I'll try to avoid causing such confusion in the future. Cheers, /F -- Federico C. F. Calboli Department of Epidemiology and Public Health Imperial College, St. Mary's Campus Norfolk Place, London W2 1PG Tel +44 (0)20 75941602 Fax +44 (0)20 75943193 f.calboli [.a.t] imperial.ac.uk f.calboli [.a.t] gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] ggplot2: font size mismatch for pdf output
Hi In the following, the graph I see on the screen and the .png output coincide. However, in the .pdf file, the fonts seem to be scaled fairly larger, resulting in the label for the top legend disappearing. Is this an infelicity or bug, or is there something I've missed? More generally, how do I control the size of fonts used in legends and axis labels? library(car) library(ggplot2) qp <-qplot (education , income , shape=type , size=women , colour=prestige , xlab="Education" , ylab="Income", data=Prestige) + scale_y_continuous(limits=c(NA, 2)) qp + scale_size(to=c(1,8)) ggsave(file="prestige-ggplot.png", width=6, height=5) # OK ggsave(file="prestige-ggplot.pdf", width=6, height=5) # fonts too large -Michael __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Collection of lm()s
try something like the following: x <- runif(1000, -3, 3) y <- 2 + 3 * x + rnorm(1000) f <- gl(10, 100) lm.lis <- vector("list", 10) for (i in 1:10) { lm.lis[[i]] <- lm(y ~ x, subset = f == as.character(i)) } sapply(lm.lis, coef) sapply(lm.lis, fitted) I hope it helps. Best, Dimitris Dimitris Rizopoulos Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven, Belgium Tel: +32/(0)16/336899 Fax: +32/(0)16/337015 Web: http://med.kuleuven.be/biostat/ http://www.student.kuleuven.be/~m0390867/dimitris.htm Quoting Chip Barnaby <[EMAIL PROTECTED]>: Hello, I would like to create a subscriptable collection (presumably a list) of lm() models. I have a data frame DX containing 6 groups of data. The general idea is (NOT RUN) ... for (i in 1:6) { DXS = subset( DX, ); LMX[ i] = lm( , data = DXS); } Now access model results by subscript ... e.g. coefficients( LMX[ 2]). Or would it be [[ 2]]? I have experimented with various schemes, attempting to "pre-allocate" a list etc. without success. Also, I assume that lm() does not make a copy of input data? That is, if I want "predict( LMX[ i])", I would have to retain a copy the associated DXS subsets? TIA, Chip Barnaby - Chip Barnaby [EMAIL PROTECTED] Vice President of Research Wrightsoft Corp. 781-862-8719 x118 voice 131 Hartwell Ave 781-861-2058 fax Lexington, MA 02421 www.wrightsoft.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Collection of lm()s
On Mon, 12 May 2008, Chip Barnaby wrote: Hello, I would like to create a subscriptable collection (presumably a list) of lm() models. I have a data frame DX containing 6 groups of data. The general idea is (NOT RUN) ... for (i in 1:6) { DXS = subset( DX, ); LMX[ i] = lm( , data = DXS); } Now access model results by subscript ... e.g. coefficients( LMX[ 2]). Or would it be [[ 2]]? I have experimented with various schemes, attempting to "pre-allocate" a list etc. without success. You want LMX <- list(6) ... LMX[[1]] <- lm( , data = DX, subset=) (No trailing ';' needed, [[]] not [] for list elements.) Also, I assume that lm() does not make a copy of input data? It may. This is controlled by the 'model' argument which defaults to TRUE (in recent versions of R: I think it once defaulted to FALSE). That is, if I want "predict( LMX[ i])", I would have to retain a copy the associated DXS subsets? It does retain the fitted values: see the 'Value' section of ?lm. TIA, Chip Barnaby - Chip Barnaby [EMAIL PROTECTED] Vice President of Research Wrightsoft Corp. 781-862-8719 x118 voice 131 Hartwell Ave 781-861-2058 fax Lexington, MA 02421 www.wrightsoft.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- 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@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to not sort factors when plotting
Hello. I'm trying find the ratios between each of the integers in a vector. I have: for (n in x) { ratio <- (x[n]/x[n-1]) ratio.all <- c(ratio.all, ratio) } Of course this doesn't work, nor does diff(n)/diff(n-1). Is there a way to specify a pair of integers in a vector? Thank you, Lydia __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Mathematical annotation in lattice strip: Is it possible?
I have tried without success to find a way including the square root symbol in lattice strips as part of my conditioning labels. I have tried supplementing by creating a list of vectors using the var.name function coupled with the expression function used in xlab/ylab. xyplot(adjusted_Rand_index~cluster|distance_measure, main="Level of agreement between partitions: Wards Method", ylab="Coefficient value (adjusted rand index)", xlab="number of clusters", type="l", data=randA1, strip=strip.custom(varnames=c(expression(sqrt(Bray-Curtis))) Is there a way of generating the square root symbol inside the strip or am I wasting my time. Thank you very much Regards Andrew http://www.nabble.com/file/p17187888/randA1.csv randA1.csv -- View this message in context: http://www.nabble.com/Mathematical-annotation-in-lattice-strip%3A-Is-it-possible--tp17187888p17187888.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] combining bar and column graphs?
Perhaps ?mosaic or ?mosaicplot _ Professor Michael Kubovy University of Virginia Department of Psychology USPS: P.O.Box 400400Charlottesville, VA 22904-4400 Parcels:Room 102Gilmer Hall McCormick RoadCharlottesville, VA 22903 Office:B011+1-434-982-4729 Lab:B019+1-434-982-4751 Fax:+1-434-982-4766 WWW:http://www.people.virginia.edu/~mk9y/ On May 11, 2008, at 11:56 AM, Me wrote: > Hi, I'm hoping to find out whether R, or an R add-on, can generate a > particular type of graph. And, more basically, whether such a type > of graph even makes sense. > > I'm looking for something resembling both a column chart and a bar > chart, where the basic visual "unit" is a solid rectangle of color > that can be extended either horizontally, vertically, or both. The > data that needs to be graphed consists of the relative contributions > of a number (6 or 8) of companies (entities, whatever) to a common > fund, over the course of a number of years (say, 1990-2008). I'm > picturing years on the X-axis, and dollar amounts on the Y-axis > (say, $0-$100,000). From a "temporal" perspective, every year will > have at least one contributor, starting with dollar zero, but some > years will have multiple contributors. From a "company" > perspective, some companies will contribute, e.g., dollars $1,001- > $5,000 for several years running, visually forming a horizontal > "block" riding on top of whatever happens to be below. > > So as a simple example, between the years 2000 and 2001, Company A > might inhabit a solid block extending from dollar zero to dollar > 1000, two years wide. In year 2000, Company B might contribute > dollars $1,001-$2,000, while right next door in year 2001, a > different Company C might contribute dollars $1,001-$10,000. > > Is it possible to have this sort of horizontal/vertical chart > generated automatically, or is this impossible? Do I need to > generate a year-on-year column graph and manually elide the > boundaries between companies' contributions in successive years, > thus forming the horizontal "blocks" I have in mind manually? Is > there perhaps another software tool that would be good for this? > > Thanks very much - this is a long question... > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Collection of lm()s
Hello, I would like to create a subscriptable collection (presumably a list) of lm() models. I have a data frame DX containing 6 groups of data. The general idea is (NOT RUN) ... for (i in 1:6) { DXS = subset( DX, ); LMX[ i] = lm( , data = DXS); } Now access model results by subscript ... e.g. coefficients( LMX[ 2]). Or would it be [[ 2]]? I have experimented with various schemes, attempting to "pre-allocate" a list etc. without success. Also, I assume that lm() does not make a copy of input data? That is, if I want "predict( LMX[ i])", I would have to retain a copy the associated DXS subsets? TIA, Chip Barnaby - Chip Barnaby [EMAIL PROTECTED] Vice President of Research Wrightsoft Corp. 781-862-8719 x118 voice 131 Hartwell Ave 781-861-2058 fax Lexington, MA 02421 www.wrightsoft.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [R-sig-ME] lme nesting/interaction advice
> But that avoids the question as to *why* it isn't very well > set up for crossed random effects? What's the problem? > What are the issues? The model is indeed bog-standard. > It would seem not unreasonable to expect that it could be > fitted in a straightforward manner, and it is irritating to > find that it cannot be. If SAS and Minitab can do it at > the touch of a button, why can't R do it? I haven't followed this thread carefully, so apologies if I'm too off base. But, in response to Rolf's questions/issues. First, SAS cannot handle models with crossed random effects (at least well at all). SAS is horribly incapable of handling even the simplest of models (especially generalized linear mixed models). I can cite numerous (recent) examples of SAS coming to a complete halt (proc nlmixed) for an analyses we were recently working on. R (and Ubuntu) was the only solution to our problem. Now, lme is not optimized for crossed random effects, but lmer is. That is why lmer is supported and lme is not really supported much. lmer is optimized for models with nested random effects and crossed random effects. When working with models with nested random effects, and software optimized for those problems (e.g., HLM, SAS, mlWin) the variance/covariance matrix forms a special, and simple structure that can be easily worked with. This is not the case for models with crossed random effects. Software packages designed for nested random effects can be tricked into handling models with crossed random effects, but this kludge is slow and really inefficient. If you want complete transparency into the why and how, here is a citation for your review. Best Harold @article{Doran:Bates:Bliese:Dowling:2007:JSSOBK:v20i02, author = "Harold Doran and Douglas Bates and Paul Bliese and Maritza Dowling", title = "Estimating the Multilevel Rasch Model: With the lme4 Package", journal = "Journal of Statistical Software", volume = "20", number = "2", pages = "1--18", day = "22", month = "2", year ="2007", CODEN = "JSSOBK", ISSN ="1548-7660", bibdate = "2007-02-22", URL = "http://www.jstatsoft.org/v20/i02";, accepted ="2007-02-22", acknowledgement = "", keywords ="", submitted = "2006-10-01", } __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Compact Indicator Matrices
On Sun, May 11, 2008 at 9:49 AM, amarkos <[EMAIL PROTECTED]> wrote: > On May 11, 4:47 pm, "Douglas Bates" <[EMAIL PROTECTED]> wrote: > >> Do you mean that you want to collapse similar rows into a single row >> and perhaps a count of the number of times that this row occurs? > > Let me rephrase the problem by providing an example. > > Input: > > A = > [,1] [,2] > [1,]11 > [2,]13 > [3,]21 > [4,]12 > [5,]21 > [6,]12 > [7,]11 > [8,]12 > [9,]13 > [10,]21 An important question here is do you start with two or more variables like the columns of your matrix A? If so, there is a more direct method of getting the answers that you want. The natural way to store such variables in R is as factors. I prefer to use letters instead of numbers to represent the levels of a factor (that way I don't confuse a factor with a numeric variable when I look at rows) so I would create a data frame with two factors instead of a matrix. > V1 <- factor(c(1,1,2,1,2,1,1,1,1,2), labels = LETTERS[1:2]) > V2 <- factor(c(1,3,1,2,1,2,1,2,3,1), labels = letters[1:3]) > df <- data.frame(f1 = V1, f2 = V2) > df f1 f2 1 A a 2 A c 3 B a 4 A b 5 B a 6 A b 7 A a 8 A b 9 A c 10 B a You could produce the indicator matrix and check for unique rows, etc. - I will show that below - but all you need is the interaction of the two factors > df$f12 <- with(df, f1:f2)[drop = TRUE] > df f1 f2 f12 1 A a A:a 2 A c A:c 3 B a B:a 4 A b A:b 5 B a B:a 6 A b A:b 7 A a A:a 8 A b A:b 9 A c A:c 10 B a B:a > str(df) 'data.frame': 10 obs. of 3 variables: $ f1 : Factor w/ 2 levels "A","B": 1 1 2 1 2 1 1 1 1 2 $ f2 : Factor w/ 3 levels "a","b","c": 1 3 1 2 1 2 1 2 3 1 $ f12: Factor w/ 4 levels "A:a","A:b","A:c",..: 1 3 4 2 4 2 1 2 3 4 > table(df$f12) A:a A:b A:c B:a 2 3 2 3 > as.numeric(df$f12) [1] 1 3 4 2 4 2 1 2 3 4 Notice that this shows you that there are four distinct combinations that occur 2, 3, 2 and 3 times respectively; the first combination occurs in rows 1 and 7, it consists of the first level of f1 and the first level of f2, etc. If you really do want the indicator matrix you could generate it as > (ind <- cbind(model.matrix(~ 0 + f1, df), model.matrix(~ 0 + f2, df))) f1A f1B f2a f2b f2c 11 0 1 0 0 21 0 0 0 1 30 1 1 0 0 41 0 0 1 0 50 1 1 0 0 61 0 0 1 0 71 0 1 0 0 81 0 0 1 0 91 0 0 0 1 10 0 1 1 0 0 > unique(ind) f1A f1B f2a f2b f2c 1 1 0 1 0 0 2 1 0 0 0 1 3 0 1 1 0 0 4 1 0 0 1 0 but working with the factors is generally much simpler than working with the indicators. > # Indicator matrix > A <- data.frame(lapply(data.frame(obj), as.factor)) > > nocases <- dim(obj)[1] > novars <- dim(obj)[2] > > # variable levels > levels.n <- sapply(obj, nlevels) > n<- cumsum(levels.n) > > # Indicator matrix calculations > Z<- matrix(0, nrow = nocases, ncol = n[length(n)]) > newdat <- lapply(obj, as.numeric) > offset <- (c(0, n[-length(n)])) > for (i in 1:novars) > Z[1:nocases + (nocases * (offset[i] + newdat[[i]] - 1))] <- 1 > > ### > > Output: > > Z = > >[,1] [,2] [,3] [,4] [,5] > [1,]10100 > [2,]10001 > [3,]01100 > [4,]10010 > [5,]01100 > [6,]10010 > [7,]10100 > [8,]10010 > [9,]10001 > [10,]01100 > > > Z is an indicator matrix in the Multiple Correspondence Analysis > framework. > My problem is to collapse identical rows (e.g. 2 and 9) into a single > row and > store the row ids. > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] rpart with circular data?
Hi Rainer, In a similar situation I used the two components of the normal vector of the surface ("northing" & "easting"). I.e. for a horizontal plane both are 0, for a vertical slope facing south northing=-1and easting=0, etc. This descartian decomposition of the slope vector avoids the problem of circularity present in the widespreadly used polar (aspect, slope) decomposition, and thus seems to suit ecological problems much better to me. However I have not looked into this much, I am also very interested in the opinion of others. Bálint On Mon, May 12, 2008 at 11:35 AM, Rainer M Krug <[EMAIL PROTECTED]> wrote: > Hi > > I am planning to use a CART analysis with rpart() to analyse the > impact of, among others, slope, altitude and aspect on mortality > rates. > My question is: > Is there a p[roblem with using aspect as a predictor as it is circular? > And if it is a problem (which I suspect), is there a transformation I > could use to transform aspect? > > Thanks > > Rainer > > -- > Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation > Biology, UCT), Dipl. Phys. (Germany) > > Plant Conservation Unit > Department of Botany > University of Cape Town > Rondebosch 7701 > South Africa > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > -- Bálint Czúcz Institute of Ecology and Botany of the Hungarian Academy of Sciences H-2163 Vácrátót, Alkotmány u. 2-4. HUNGARY Tel: +36 28 360122/157 +36 70 7034692 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Optimization problem, to minimize the length(rle(B)$lengths) for all the rows and columns
Hi, how can I order the rows and columns of a matrix A to generate B, in order to minimize the length(rle(B)$lengths) for all the rows and columns ? > set.seed(5) > a <- matrix(rnorm(200), nrow=20) > a[a<=0] <- 0 > a[a>0] <- 1 > a [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,]011010100 1 [2,]110101101 0 [3,]010011011 1 [4,]111101011 1 [5,]110010100 0 [6,]001001100 0 [7,]010010010 0 [8,]010000100 0 [9,]000000000 1 [10,]100001000 0 [11,]111010110 1 [12,]011001011 1 [13,]011111000 0 [14,]010111010 1 [15,]010110100 0 [16,]010011111 1 [17,]001100101 1 [18,]000101010 1 [19,]100101101 0 [20,]001000110 1 [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [R-sig-ME] lme nesting/interaction advice
I *think* the syntax for the model Federico wants is this: lmer(y~selection*males+ (selection|month) + (males|month)) My lme syntax is a bit rusty, so I'm not confident how to recode with nested random effects, as in P&B p24. Two quick points: 1. I think Federico has caused some confusion on account of the way he used the term 'crossing'2. Whilst interesting reading, some of the content of this thread does not conform to the list's rules and regs. Best wishes, Nick [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Insert a recorde into a table using SQL
It works. Thanks very much. Best On Mon, May 12, 2008 at 7:02 PM, Gabor Grothendieck <[EMAIL PROTECTED]> wrote: > Create an encoding function which replaces single quotes with > two single quotes: > > # first string is a single character consisting of single quote > # second string is two characters consisting of two single quotes > enc <- function(x) gsub("'", "''", x) > dbGetQuery(con,sprintf("insert into dd (txt) values ('%s')", enc(dd[2,1]))) > > > On Mon, May 12, 2008 at 1:54 AM, ronggui <[EMAIL PROTECTED]> wrote: >> Dear list, >> >> I want to insert a recorde into a SQLite table by dbGetQuery(), but >> there is a problem when the value contains quotation mark. >> >> > dd<-data.frame(txt=c("having both ' and \" in character.","OK")) >> > library(RSQLite) >> Loading required package: DBI >> > con<-dbConnect(dbDriver("SQLite"),":memory:") >> > dbWriteTable(con,"dd",dd,over=T) >> [1] TRUE >> > dbGetQuery(con,sprintf("insert into dd (txt) values (\"%s\")",dd[2,1])) >> NULL >> > dbGetQuery(con,sprintf("insert into dd (txt) values (\"%s\")",dd[1,1])) >> Error in sqliteExecStatement(con, statement, bind.data) : >> RS-DBI driver: (error in statement: unrecognized token: "")") >> >> How can I insert a (key, value) pair into a table by dbGetQuery? Thanks. >> >> -- >> HUANG Ronggui, Wincent >> Bachelor of Social Work, Fudan University, China >> Master of sociology, Fudan University, China >> Ph.D. Candidate, CityU of HK. >> >> __ >> R-help@r-project.org mailing list >> https://stat.ethz.ch/mailman/listinfo/r-help >> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html >> and provide commented, minimal, self-contained, reproducible code. >> > -- HUANG Ronggui, Wincent Bachelor of Social Work, Fudan University, China Master of sociology, Fudan University, China Ph.D. Candidate, CityU of HK, http://www.cityu.edu.hk/sa/psa_web2006/students/rdegree/huangronggui.html __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Quadratic Constraints
Hello R, By any chance, Rdonlp2 package of R defined in http://arumat.net/Rdonlp2/tutorial.html#SECTION0002 serves the below purpose? BR, Shubha -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Shubha Vishwanath Karanth Sent: Monday, May 12, 2008 3:27 PM To: [EMAIL PROTECTED] Subject: [R] Quadratic Constraints Hi R, A quick question How can I optimize the objective function constrained to quadratic constraints? Which function of R is useful for quadratic constraints? Many Thanks, Shubha This e-mail may contain confidential and/or privileged i...{{dropped:13}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. This e-mail may contain confidential and/or privileged i...{{dropped:10}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Insert a recorde into a table using SQL
Create an encoding function which replaces single quotes with two single quotes: # first string is a single character consisting of single quote # second string is two characters consisting of two single quotes enc <- function(x) gsub("'", "''", x) dbGetQuery(con,sprintf("insert into dd (txt) values ('%s')", enc(dd[2,1]))) On Mon, May 12, 2008 at 1:54 AM, ronggui <[EMAIL PROTECTED]> wrote: > Dear list, > > I want to insert a recorde into a SQLite table by dbGetQuery(), but > there is a problem when the value contains quotation mark. > > > dd<-data.frame(txt=c("having both ' and \" in character.","OK")) > > library(RSQLite) > Loading required package: DBI > > con<-dbConnect(dbDriver("SQLite"),":memory:") > > dbWriteTable(con,"dd",dd,over=T) > [1] TRUE > > dbGetQuery(con,sprintf("insert into dd (txt) values (\"%s\")",dd[2,1])) > NULL > > dbGetQuery(con,sprintf("insert into dd (txt) values (\"%s\")",dd[1,1])) > Error in sqliteExecStatement(con, statement, bind.data) : > RS-DBI driver: (error in statement: unrecognized token: "")") > > How can I insert a (key, value) pair into a table by dbGetQuery? Thanks. > > -- > HUANG Ronggui, Wincent > Bachelor of Social Work, Fudan University, China > Master of sociology, Fudan University, China > Ph.D. Candidate, CityU of HK. > > __ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] lme nesting/interaction advice
On 12 May 2008, at 11:16, Andrew Robinson wrote: Well. I have documentation relevant to nlme that goes back about 10 years. I don't know when it was first added to S-plus, but I assume that it was about then. Now, do you think that if the thing that you want to do was really bog standard, that noone would have raised a fuss or solved it within 10 years? I'm pretty unpleasant, more so in person, so I'll tell you this. If people raised the issue and got the answer I got, I would not be surprised if they'd migrated to 'any other stats software' in droves. I have no doubt that, given the cryptic and sparse nature of the documentation of the issue, most people migrated well before asking -- on the grounds most people have a job to do, papers to publish, grants to write, kids to pick up from school and pretty little time for RTFM and all that sanctimonious attitude. Once people stop nagging about 'whatever', the reason is because they finally got the message things ain't gonna improve, so cut your losses and look elsewhere. Being unpleasant, thick skinned and cheap I will keep nagging and use R (the fact I do like it very much might be a factor). But given the selection users go through, it will be Vogon time sooner or later ;). /F -- Federico C. F. Calboli Department of Epidemiology and Public Health Imperial College, St. Mary's Campus Norfolk Place, London W2 1PG Tel +44 (0)20 75941602 Fax +44 (0)20 75943193 f.calboli [.a.t] imperial.ac.uk f.calboli [.a.t] gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] lme nesting/interaction advice
On Mon, May 12, 2008 at 10:50:03AM +0100, Federico Calboli wrote: > > On 12 May 2008, at 01:05, Andrew Robinson wrote: > > >On Mon, May 12, 2008 at 10:34:40AM +1200, Rolf Turner wrote: > >> > >>On 12/05/2008, at 9:45 AM, Andrew Robinson wrote: > >> > >>>On Sun, May 11, 2008 at 07:52:50PM +0100, Federico Calboli wrote: > > The main point of my question is, having a 3 way anova (or > ancova, if > you prefer), with *no* nesting, 2 fixed effects and 1 random > effect, > why is it so boneheaded difficult to specify a bog standard fully > crossed model? I'm not talking about some rarified esoteric model > here, we're talking about stuff tought in a first year Biology > Stats > course here[1]. > >>> > >>>That may be so, but I've never needed to use one. > >> > >>So what? This is still a standard, common, garden-variety > >>model that you will encounter in exercises in many (if not > >>all!) textbooks on experimental design and anova. > > > >To reply in similar vein, so what? Why should R-core or the R > >community feel it necessary to reproduce every textbook example? How > >many times have *you* used such a model in real statistical work, > >Rolf? > > There is a very important reason why R (or any other stats package) > should *easily* face the challenge of bog standard models: because it > is a *tool* for an end (i.e. the analysis of data to figure out what > the heck they tell us) rather than a end in itself. But a tool that mostly (entirely?) appears in textbooks. > Bog standard models are *likely* to be used over and over again > because they are *bog standard*, and they became such by being used > *lots*. Well. I have documentation relevant to nlme that goes back about 10 years. I don't know when it was first added to S-plus, but I assume that it was about then. Now, do you think that if the thing that you want to do was really bog standard, that noone would have raised a fuss or solved it within 10 years? > If someone with a relatively easy model cannot use R for his job s/he > will use something else, and the R community will *not* increase in > numbers. Since R is a *community driven project*, you do the math on > what that would mean in the long run. Fewer pestering questions? ;) Andrew -- Andrew Robinson Department of Mathematics and StatisticsTel: +61-3-8344-6410 University of Melbourne, VIC 3010 Australia Fax: +61-3-8344-4599 http://www.ms.unimelb.edu.au/~andrewpr http://blogs.mbs.edu/fishing-in-the-bay/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [R-sig-ME] lme nesting/interaction advice
On 12 May 2008, at 10:05, Ken Beath wrote: There is only one random effect, so where does the crossing come from ? The fixed effects vary across blocks, but they are fixed so are just covariates. For this type of data the usual model in lme4 is y~fixed1+fixed2+1|group and for lme split into fixed and random parts. First off, whoa, an helpful reply! thanks for that, I hope I won't sound sarcastic or aggressive because I do not mean to be either. Regarding your comment, the experiment was replicated three times, in 3 different months. I would argue that for the fixed effects to be meaningful, they must have an effect over an above the effect:month interaction (because each fixed effect, and their interaction, might vary between each replicate). I would then argue I need to calculate 1) fixed.effect1:random.effect 2) fixed.effect2:random.effect 3) fixed.effect1:fixed.effect2:random.effect to test if fixed.effect1 is meaningful (and use 1) as the error); if fixed.effect2 has is meaningful (and use 2) as the error); fixed.effect1:fixed.effect2 is meaningful (and use 3) as the error). I'm happy to be correct if I am wrong here. The problems seems to be that you want lme to work in the same way as an ANOVA table and it doesn't. The secret with lme and lme4 is to think about the structure of the data and describe with an equation. Then each term in the equation corresponds to part of the model definition in R. I'll try to do that. Once I have sorted how to specify such trivial model I'll face the horror of the nesting, in any case I attach a toy dataset I created especially to test how to specify the correct model (silly me). I'm a bit lost with your data file, it has 4 covariates, which is more than enough for 2 fixed effects, assuming block is the grouping and y the outcome. In the data file, 'selection' and 'males' are fixed effects, and 'month' is the effect I am using for the model we are discussing here. The y was generatde with runif() just to have something, I'm not expecting any intersting result, just to understand how to fit the right model. In the dataset 'line' is nested within 'selection' and 'block' is nested within 'month'. That's the nesting I will have to take into account once I get the more straightforward (sic!) model we're discussing right. Best, Federico -- Federico C. F. Calboli Department of Epidemiology and Public Health Imperial College, St. Mary's Campus Norfolk Place, London W2 1PG Tel +44 (0)20 75941602 Fax +44 (0)20 75943193 f.calboli [.a.t] imperial.ac.uk f.calboli [.a.t] gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Quadratic Constraints
Hi R, A quick question How can I optimize the objective function constrained to quadratic constraints? Which function of R is useful for quadratic constraints? Many Thanks, Shubha This e-mail may contain confidential and/or privileged i...{{dropped:13}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] mgcv::gam shrinkage of smooths
On Tuesday 06 May 2008 23:34, David Katz wrote: > In Dr. Wood's book on GAM, he suggests in section 4.1.6 that it might be > useful to shrink a single smooth by adding S=S+epsilon*I to the penalty > matrix S. The context was the need to be able to shrink the term to zero if > appropriate. I'd like to do this in order to shrink the coefficients > towards zero (irrespective of the penalty for "wiggliness") - but not > necessarily all the way to zero. IE, my informal prior is to keep the > contribution of a specific term small. > > 1) Is adding eps*I to the penalty matrix an effective way to achieve this > goal? > > 2) How do I accomplish this in practice using mgcv::gam? Are you saying that you would like to specify the amount of shrinkage directly, with the degree of shrinkage set `by hand' or do you want mgcv::gam to estimate the degree of shrinkage? If you want to supply a fixed ridge penalty then use argument `H' of mgcv::gam, which is designed for just this sort of purpose. If you want what's suggested in section 4.1.6 then try s(...,bs="ts") or s(...,bs="cs"). If you want to add an extra ridge penalty to a smooth and have mgcv::gam estimate its smoothing parameter, then you would need to add write a smoother class to do this (by modifying one of the existing ones). See ?p.spline, or exercise 8, chapter 5 of Wood (2006) "GAMs:An Intro with R". best, Simon -- > Simon Wood, Mathematical Sciences, University of Bath, Bath, BA2 7AY UK > +44 1225 386603 www.maths.bath.ac.uk/~sw283 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] lme nesting/interaction advice
On 12 May 2008, at 01:05, Andrew Robinson wrote: On Mon, May 12, 2008 at 10:34:40AM +1200, Rolf Turner wrote: On 12/05/2008, at 9:45 AM, Andrew Robinson wrote: On Sun, May 11, 2008 at 07:52:50PM +0100, Federico Calboli wrote: The main point of my question is, having a 3 way anova (or ancova, if you prefer), with *no* nesting, 2 fixed effects and 1 random effect, why is it so boneheaded difficult to specify a bog standard fully crossed model? I'm not talking about some rarified esoteric model here, we're talking about stuff tought in a first year Biology Stats course here[1]. That may be so, but I've never needed to use one. So what? This is still a standard, common, garden-variety model that you will encounter in exercises in many (if not all!) textbooks on experimental design and anova. To reply in similar vein, so what? Why should R-core or the R community feel it necessary to reproduce every textbook example? How many times have *you* used such a model in real statistical work, Rolf? There is a very important reason why R (or any other stats package) should *easily* face the challenge of bog standard models: because it is a *tool* for an end (i.e. the analysis of data to figure out what the heck they tell us) rather than a end in itself. Bog standard models are *likely* to be used over and over again because they are *bog standard*, and they became such by being used *lots*. If someone with a relatively easy model cannot use R for his job s/he will use something else, and the R community will *not* increase in numbers. Since R is a *community driven project*, you do the math on what that would mean in the long run. Regards, Federico -- Federico C. F. Calboli Department of Epidemiology and Public Health Imperial College, St. Mary's Campus Norfolk Place, London W2 1PG Tel +44 (0)20 75941602 Fax +44 (0)20 75943193 f.calboli [.a.t] imperial.ac.uk f.calboli [.a.t] gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] lme nesting/interaction advice
On 11 May 2008, at 23:34, Rolf Turner wrote: It doesn't seem to me to be a complaint as such. It is a request for insight. I too would like some insight as to what on earth is going on. And why do you say Federico shows no evidence of having searched the archives? One can search till one is blue in the face and come away no wiser on this issue. cheers, Rolf Turner Cheers for the support Rolf. I have searched the archives, and have the Pinheiro and Bates book in front of my nose (plus MASS4 and many others). The bottom line here is, I'm pretty cool with RTFM and all that, my problem is that I do not have a clear FM to read about my issue, and hence I have to pester (because that how people seem to feel) the list. I apologise for asking an inconvenient question. Fede -- Federico C. F. Calboli Department of Epidemiology and Public Health Imperial College, St. Mary's Campus Norfolk Place, London W2 1PG Tel +44 (0)20 75941602 Fax +44 (0)20 75943193 f.calboli [.a.t] imperial.ac.uk f.calboli [.a.t] gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] lme nesting/interaction advice
On 12 May 2008, at 09:29, Dieter Menne wrote: Federico: First, mixed models are different from "standard 101 Anova", and quite a lot of the nesting stuff I used to ponder about 30 year ago when I started teaching this is no longer relevant and works implicitely when you code the parameters correctly. with effect3 being random (all all the jazz that comes from this fact). I fully apprecciate that the only reasonable F-tests would be for effect1, effect2 and effect1:effect2, but there is no way I can use lme to specify such simple thing without getting the *wrong* denDF. >> Good to know that you are sure what is "right"; probably == SAS. Since most people active in the lme-business have read http://wiki.r-project.org/rwiki/doku.php?id=guides:lmer-tests http://finzi.psych.upenn.edu/R/Rhelp02a/archive/76742.html carefully, you might be rather lonely. I will. While I do, feel free to have a look at Appendix A.3 (page App6, at the end of the book) of the Zar 'Biostatistical Analysis', IV ed, second table from the top. That's where I get the feeling for what's right or wrong. I surely cannot get it from SAS because I never had it. I never had the budget for it, so much so I had to lear how to use R from the start because it was free and that was the budget of my department had for stats software All in all, if you feel statistical analysis has moved forth from such humble beginnings (the book I mean, not SAS), and you can convince of that every ref for every paper you submit, please do tell me how you do it, it would be more valuable than knowing how to fit my model. Cheers, Federico -- Federico C. F. Calboli Department of Epidemiology and Public Health Imperial College, St. Mary's Campus Norfolk Place, London W2 1PG Tel +44 (0)20 75941602 Fax +44 (0)20 75943193 f.calboli [.a.t] imperial.ac.uk f.calboli [.a.t] gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Cumulative lattice histograms
Ola Caster gmail.com> writes: > It's fairly straightforward to plot cumulative histograms using the hist() > function. You do something like: > > h <- hist(rnorm(100), plot=FALSE) > h$counts<- cumsum(h$counts) > plot(h) > > However, I have failed to find any example where this is done using the > lattice histogram() function. This is not a full solution to your problem with the red line, but at least comes close (hope so) to your original hist with lattic library(lattice) h = sort(rnorm(100),decreasing=TRUE) df = data.frame(h=h,cum=cumsum(h)) histogram(h~cum,data=df) Dieter __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.