[R] using mclapply (multi core apply) to do matrix multiplication
Dear all, I am trying to multiply three different matrices and each matrice is of size 16384,16384 the normal %*% multiplciation operator has not finished one day now. As I am running a system with many cores (and it seems that R is using only one of those) I would like to write fast a brief function that converts the typical for loops of a matrix multiplication to a set of lapply sets (mclapply uses the lapply syntax but it applies the work to many cores). If my thinking is correct , in the sense that this will speed up things a lot, I want you to help me covert the first matrix in rows the second in columns and convert those in a format that lapply would like to work with. I would like to thank you in advance for your help Regards Alex [[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] A question on p-value of Unit Root Tests using “urca” toolbox
A question on p-value of Unit Root Tests using urca toolbox Theres a function punitroot on unit root probability value punitroot(q, N = Inf, trend = c(c, nc, ct, ctt), statistic = c(t, n), na.rm = FALSE) To my knowledge, c means const or intercept, nc means no const or intercept, and ct means time trend with an intercept. What does ctt mean? Isnt punitroot a function of lag period? Why isnt lag an input of the function? BTW, in the output of the unit root test, does tau2 means t statistic? Thank you for your attention! [[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] Writing to a file
Thank you a lot for answering so fast! but..what do you mean by example? I 've mentioned above the loop I used and I also show how the file looks like 'cause its huge. the way i read the file is x=read.table(filename.txt,header=FALSE,sep=\t,fill=TRUE) y=x[1:45,] (i use only some rows in order to test if it works ) -- View this message in context: http://r.789695.n4.nabble.com/Writing-to-a-file-tp3070617p4364034.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] Flexmix new data classification
Hi I built a flexmix GLM binomial model with 200 observations and the model gave me 2 clusters, so if the model is named as newModel then i get the cluster index for each row using newModel@clusters. Now is there any way to predict which cluster the new observation or 201 observation belongs to using the above built model (newModel) ie so 201 observation can either belong to cluster 1 or cluster 2. Thanks -- View this message in context: http://r.789695.n4.nabble.com/Flexmix-new-data-classification-tp4363996p4363996.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] Flexmix new data classification
?clusters tells you how to get class membership for new data (and the appropriate way of getting cluster membership is to use clusters(newModel) rather newModel@clusters) hth, Ingmar On Tue, Feb 7, 2012 at 8:48 AM, loyolite270 loyolite...@gmail.com wrote: Hi I built a flexmix GLM binomial model with 200 observations and the model gave me 2 clusters, so if the model is named as newModel then i get the cluster index for each row using newModel@clusters. Now is there any way to predict which cluster the new observation or 201 observation belongs to using the above built model (newModel) ie so 201 observation can either belong to cluster 1 or cluster 2. Thanks -- View this message in context: http://r.789695.n4.nabble.com/Flexmix-new-data-classification-tp4363996p4363996.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. [[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] Flexmix new data classification
Hi, I built a flexmix GLM binomial model with 200 observations and the model gave me 2 clusters, so if the model is named as newModel then i get the cluster index for each row using newModel@clusters. Now is there any way to predict which cluster the new observation or 201 observation belongs to using the above built model (newModel) ie so 201 observation can either belong to cluster 1 or cluster 2. You can obtain the predicted cluster memberships from the fitted model using the accessor function clusters(), i.e., clusters(newModel). If you want to predict the cluster memberships of new observations, you can then use clusters(newModel, newdata = data_frame_with_new_observations) HTH, Bettina Thanks -- View this message in context:http://r.789695.n4.nabble.com/Flexmix-new-data-classification-tp4363996p4363996.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 guidehttp://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- --- Bettina Grün Institut für Angewandte Statistik / IFAS Johannes Kepler Universität Linz Altenbergerstraße 69 4040 Linz, Austria Tel: +43 732 2468-6829 Fax: +43 732 2468-6800 E-Mail:bettina.gr...@jku.at www.ifas.jku.at __ 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 to a file
Hi Thank you a lot for answering so fast! but..what do you mean by example? I 've mentioned above the loop I used and I also show how the file looks like I do not see any loop. I do not archive all posts from R help, only those with interesting answers :-) and if you do not keep the context in future mails for those not using nabble it is lost and it would be necessary to dig in r help archive. 'cause its huge. the way i read the file is x=read.table(filename.txt,header=FALSE,sep=\t,fill=TRUE) y=x[1:45,] maybe you can use even smaller fraction for a data example y-x[1:10,] dput(y) and copy the output from dput to your mail is the easiest way. Regards Petr (i use only some rows in order to test if it works ) -- View this message in context: http://r.789695.n4.nabble.com/Writing-to-a- file-tp3070617p4364034.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@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] Table rearranging
Thank you for your help, Bill. From the original table (not the plyr output), I would like to remove all the lines that do not have a corresponding open/closed measurement. For example, if there is a Closed yellow measurement on 0917, but not an Open yellow 0917 measurement, then the Closed yellow should be deleted. How can I make this change? Jeffrey From: wdun...@tibco.com To: johjeff...@hotmail.com; r-help@r-project.org Subject: RE: [R] Table rearranging Date: Tue, 7 Feb 2012 00:43:25 + Install and load the plyr package and try something like: ddply(d, .(date, color), summarize, + meanOpen=mean(measurement[door==open]), nOpen=sum(door==open), + meanClosed=mean(measurement[door==closed]), nClosed=sum(door==closed)) date color meanOpen nOpen meanClosed nClosed 1 420 red 0.9741929 1 NaN 0 2 513 red 0.9352938 1 0.9620535 1 3 917 yellow NaN 0 0.9941022 1 4 1230 blue 0.9639099 1 0.9893108 1 5 1230 green 0.9765203 1 NaN 0 Bill Dunlap Spotfire, TIBCO Software wdunlap tibco.com -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Jeffrey Joh Sent: Monday, February 06, 2012 4:28 PM To: r-help@r-project.org Subject: [R] Table rearranging I have a table that looks like this: measurementdatedoorcolor 0.93529385513openred 0.97419293420openred 0.962053514513closedred 0.9639099371230openblue 0.976520341230opengreen 0.9893107951230closedblue 0.9941022917closedyellow I would like to create a table that has: Open measurement, Closed measurement, date, color. For every date/color combination, there should be two columns to represent the door open/closed measurement. If there are multiple datapoints with a given door/date/color combination, then they should be averaged. I would also like to make two columns to represent the number of datapoints that were averaged in determining the open/closed measurements. Jeffrey __ 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] xtable beta testing wanted
Dear R-users, I've enhanced the xtable package, adding numerous features that have been requested by users. The changes are listed below. The objective throughout has been to avoid any breaking changes. However, as xtable is widely used and is a dependency of many packages I'd like to have others try it out before I post it to CRAN. Both bug reports and design feedback is welcomed. In particular, I'd like to have package authors that have written xtable methods try it out. If you are interested and willing, you can get the source from R-Forge via SVN: svn checkout svn://scm.r-forge.r-project.org/svnroot/xtable/ I can also send source and Windows binary versions to anyone that is interested. Feedback can be sent to me or posted in the issue tracking system on R-Forge: http://r-forge.r-project.org/tracker/?group_id=1228 Note that the R-Forge build system is currently being updated, so the option of installing directly from R-Forge won't get the latest version of the package. Thanks, Charlie Roosen xtable maintainer -- NEWS 1.7-0 (NOT YET RELEASED) * Added some vectorization code to improve performance. * Let caption be length 2, in which case the second value is the short caption used when creating a list of tables. * Added toLatex method. * Included print.xtable in the exported methods in the NAMESPACE file. * Added print.results argument to print that can be used to suppress the printing. * Added format.args argument to print that can be used to pass additional arguments such as big.marks to formatC(). * Added rotate.colnames and rotate.rownames arguments to print.xtable. * Added booktabs argument to use the \toprule, \midrule, and \bottomrule tags from the Latex booktabs package rather than using \hline for all horizontal lines. * Added scalebox argument to include a \scalebox clause around the tabular environment with the specified value used as the scaling factor. * Changed the print.xtable() arguments to use getOption() to check the options for a default value. This was suggested since print.xtable() has a lot of arguments that the user will typically leave unchanged between tables. * Added an is.null() check on the table.placement argument. * Added examples using the new arguments to the vignette. Charles Roosen, PhD Technical Director T: +41 (0)61 206 92 91 M: +41 (0)79 248 70 71 F: +41 (0) 61 206 92 99 www.mango-solutions.com Mango Solutions AG Aeschenvorstadt 36 4051 Basel Switzerland LEGAL NOTICE This message is intended for the use o...{{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] using mclapply (multi core apply) to do matrix multiplication
7-02-2012, 00:29 (-0800); Alaios escriu: Dear all, I am trying to multiply three different matrices and each matrice is of size 16384,16384 the normal %*% multiplciation operator has not finished one day now. As I am running a system with many cores (and it seems that R is using only one of those) I would like to write fast a brief function that converts the typical for loops of a matrix multiplication to a set of lapply sets (mclapply uses the lapply syntax but it applies the work to many cores). If my thinking is correct , in the sense that this will speed up things a lot, I want you to help me covert the first matrix in rows the second in columns and convert those in a format that lapply would like to work with. If I understand correctly, R uses a specialized library called BLAS to do matrix multiplications. I doubt re-implementing the matrix multiplication code at R-level would be any faster. What you can try is replace BLAS with a multicore version of BLAS although it's not easy if you have to compile it yourself. Also, you may try to re-think the problem you're trying to solve. Maybe there's a different approach that is less computation-intensive. -- Cheers, Ernest __ 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] using mclapply (multi core apply) to do matrix multiplication
I would like to thank you Ernest for your answer. I guess that this is gonna be faster as right now R only sees one core. In my work there is a system with 64 cores and you can see only one working. If I understand it right a [m,n][n,k] matrix multiplication can be split into rows (from first matrice) and columns (from the second matrice) and then combine all the local results of each cpu together. Would that be too weird for mclapply to handle? B.R Alex From: Ernest Adrogué nfdi...@gmail.com To: r-help@r-project.org Sent: Tuesday, February 7, 2012 11:02 AM Subject: Re: [R] using mclapply (multi core apply) to do matrix multiplication 7-02-2012, 00:29 (-0800); Alaios escriu: Dear all, I am trying to multiply three different matrices and each matrice is of size 16384,16384 the normal %*% multiplciation operator has not finished one day now. As I am running a system with many cores (and it seems that R is using only one of those) I would like to write fast a brief function that converts the typical for loops of a matrix multiplication to a set of lapply sets (mclapply uses the lapply syntax but it applies the work to many cores). If my thinking is correct , in the sense that this will speed up things a lot, I want you to help me covert the first matrix in rows the second in columns and convert those in a format that lapply would like to work with. If I understand correctly, R uses a specialized library called BLAS to do matrix multiplications. I doubt re-implementing the matrix multiplication code at R-level would be any faster. What you can try is replace BLAS with a multicore version of BLAS although it's not easy if you have to compile it yourself. Also, you may try to re-think the problem you're trying to solve. Maybe there's a different approach that is less computation-intensive. -- Cheers, Ernest __ 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] using mclapply (multi core apply) to do matrix multiplication
On 07/02/12 11:31, Alaios wrote: I would like to thank you Ernest for your answer. I guess that this is gonna be faster as right now R only sees one core. In my work there is a system with 64 cores and you can see only one working. If I understand it right a [m,n][n,k] matrix multiplication can be split into rows (from first matrice) and columns (from the second matrice) and then combine all the local results of each cpu together. You definitaly can go this way, but I would STRONGLY recommend to search for parallel BLAS, check in the R-admin manual the section Linear Algebra which deals with BLAS et al, and e.g. http://www.r-bloggers.com/compiling-64-bit-r-2-10-1-with-mkl-in-linux/ My guess is that a paralelization on the C level in the BLAS et al. library will be MUCH faster then a paralelization on R level. Also, there is a R-sig-hpc mailing list for these kind of questions. Cheers, Rainer Would that be too weird for mclapply to handle? B.R Alex From: Ernest Adroguénfdi...@gmail.com To: r-help@r-project.org Sent: Tuesday, February 7, 2012 11:02 AM Subject: Re: [R] using mclapply (multi core apply) to do matrix multiplication 7-02-2012, 00:29 (-0800); Alaios escriu: Dear all, I am trying to multiply three different matrices and each matrice is of size 16384,16384 the normal %*% multiplciation operator has not finished one day now. As I am running a system with many cores (and it seems that R is using only one of those) I would like to write fast a brief function that converts the typical for loops of a matrix multiplication to a set of lapply sets (mclapply uses the lapply syntax but it applies the work to many cores). If my thinking is correct , in the sense that this will speed up things a lot, I want you to help me covert the first matrix in rows the second in columns and convert those in a format that lapply would like to work with. If I understand correctly, R uses a specialized library called BLAS to do matrix multiplications. I doubt re-implementing the matrix multiplication code at R-level would be any faster. What you can try is replace BLAS with a multicore version of BLAS although it's not easy if you have to compile it yourself. Also, you may try to re-think the problem you're trying to solve. Maybe there's a different approach that is less computation-intensive. __ 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. -- Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation Biology, UCT), Dipl. Phys. (Germany) Centre of Excellence for Invasion Biology Stellenbosch University South Africa Tel : +33 - (0)9 53 10 27 44 Cell: +33 - (0)6 85 62 59 98 Fax : +33 - (0)9 58 10 27 44 Fax (D):+49 - (0)3 21 21 25 22 44 email: rai...@krugs.de Skype: RMkrug __ 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] using mclapply (multi core apply) to do matrix multiplication
Thank you very much for your point... I hope I can find some easy to follow instructions as I do not have root permission for the many cores system and our system administrator want to have easy instructions to follow. Thanks a gain. From: Rainer M Krug r.m.k...@gmail.com Cc: Ernest Adrogué nfdi...@gmail.com; r-help@r-project.org r-help@r-project.org Sent: Tuesday, February 7, 2012 11:44 AM Subject: Re: [R] using mclapply (multi core apply) to do matrix multiplication On 07/02/12 11:31, Alaios wrote: I would like to thank you Ernest for your answer. I guess that this is gonna be faster as right now R only sees one core. In my work there is a system with 64 cores and you can see only one working. If I understand it right a [m,n][n,k] matrix multiplication can be split into rows (from first matrice) and columns (from the second matrice) and then combine all the local results of each cpu together. You definitaly can go this way, but I would STRONGLY recommend to search for parallel BLAS, check in the R-admin manual the section Linear Algebra which deals with BLAS et al, and e.g. http://www.r-bloggers.com/compiling-64-bit-r-2-10-1-with-mkl-in-linux/ My guess is that a paralelization on the C level in the BLAS et al. library will be MUCH faster then a paralelization on R level. Also, there is a R-sig-hpc mailing list for these kind of questions. Cheers, Rainer Would that be too weird for mclapply to handle? B.R Alex From: Ernest Adroguénfdi...@gmail.com To: r-help@r-project.org Sent: Tuesday, February 7, 2012 11:02 AM Subject: Re: [R] using mclapply (multi core apply) to do matrix multiplication 7-02-2012, 00:29 (-0800); Alaios escriu: Dear all, I am trying to multiply three different matrices and each matrice is of size 16384,16384 the normal %*% multiplciation operator has not finished one day now. As I am running a system with many cores (and it seems that R is using only one of those) I would like to write fast a brief function that converts the typical for loops of a matrix multiplication to a set of lapply sets (mclapply uses the lapply syntax but it applies the work to many cores). If my thinking is correct , in the sense that this will speed up things a lot, I want you to help me covert the first matrix in rows the second in columns and convert those in a format that lapply would like to work with. If I understand correctly, R uses a specialized library called BLAS to do matrix multiplications. I doubt re-implementing the matrix multiplication code at R-level would be any faster. What you can try is replace BLAS with a multicore version of BLAS although it's not easy if you have to compile it yourself. Also, you may try to re-think the problem you're trying to solve. Maybe there's a different approach that is less computation-intensive. __ 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. -- Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation Biology, UCT), Dipl. Phys. (Germany) Centre of Excellence for Invasion Biology Stellenbosch University South Africa Tel : +33 - (0)9 53 10 27 44 Cell: +33 - (0)6 85 62 59 98 Fax : +33 - (0)9 58 10 27 44 Fax (D): +49 - (0)3 21 21 25 22 44 email: rai...@krugs.de Skype: RMkrug [[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] Positioning text in top left corner of plot
Dear all, another questions related to zoo plotting. I would like to do as in the subject. Here a reproducible code: library(zoo) par(mfrow=c(2,1) plot(zoo(seq(1:10),as.Date(seq(1:10),origin=1970-01-01)),xlab=,ylab=,main=Value,las=1) mtext(EUR billions,adj=0,cex=0.7) plot(zoo(seq(1:10),as.Date(seq(1:10),origin=1970-01-01)),xlab=,ylab=,main=Value,las=1) par(xpd=T) text(par(usr)[1],par(usr)[3]+10.5,EUR billions,cex=0.7) In the first graph I use the mtext function, which does the trick but it places the text too close too the y-axis. A second possibility is to use the text function, prior specification of the option parameter par(xpd=T). However, this does not look consistent as, depending on the graph, I would need to add different numbers in the y specification of the axis (par(usr)[3]+XXX). Any idea? -- View this message in context: http://r.789695.n4.nabble.com/Positioning-text-in-top-left-corner-of-plot-tp831723p4364355.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] using mclapply (multi core apply) to do matrix multiplication
On 07/02/12 12:02, Alaios wrote: Thank you very much for your point... I hope I can find some easy to follow instructions as I do not have root permission Me neither on pur cluster - but that won't stop you from compiling and installing R in your home directory. By doing this, you have even more control. Cheers and good luck, Rainer for the many cores system and our system administrator want to have easy instructions to follow. Thanks a gain. *From:* Rainer M Krug r.m.k...@gmail.com *To:* Alaios ala...@yahoo.com *Cc:* Ernest Adrogué nfdi...@gmail.com; r-help@r-project.org r-help@r-project.org *Sent:* Tuesday, February 7, 2012 11:44 AM *Subject:* Re: [R] using mclapply (multi core apply) to do matrix multiplication On 07/02/12 11:31, Alaios wrote: I would like to thank you Ernest for your answer. I guess that this is gonna be faster as right now R only sees one core. In my work there is a system with 64 cores and you can see only one working. If I understand it right a [m,n][n,k] matrix multiplication can be split into rows (from first matrice) and columns (from the second matrice) and then combine all the local results of each cpu together. You definitaly can go this way, but I would STRONGLY recommend to search for parallel BLAS, check in the R-admin manual the section Linear Algebra which deals with BLAS et al, and e.g. http://www.r-bloggers.com/compiling-64-bit-r-2-10-1-with-mkl-in-linux/ My guess is that a paralelization on the C level in the BLAS et al. library will be MUCH faster then a paralelization on R level. Also, there is a R-sig-hpc mailing list for these kind of questions. Cheers, Rainer Would that be too weird for mclapply to handle? B.R Alex From: Ernest Adroguénfdi...@gmail.com mailto:nfdi...@gmail.com To: r-help@r-project.org mailto:r-help@r-project.org Sent: Tuesday, February 7, 2012 11:02 AM Subject: Re: [R] using mclapply (multi core apply) to do matrix multiplication 7-02-2012, 00:29 (-0800); Alaios escriu: Dear all, I am trying to multiply three different matrices and each matrice is of size 16384,16384 the normal %*% multiplciation operator has not finished one day now. As I am running a system with many cores (and it seems that R is using only one of those) I would like to write fast a brief function that converts the typical for loops of a matrix multiplication to a set of lapply sets (mclapply uses the lapply syntax but it applies the work to many cores). If my thinking is correct , in the sense that this will speed up things a lot, I want you to help me covert the first matrix in rows the second in columns and convert those in a format that lapply would like to work with. If I understand correctly, R uses a specialized library called BLAS to do matrix multiplications. I doubt re-implementing the matrix multiplication code at R-level would be any faster. What you can try is replace BLAS with a multicore version of BLAS although it's not easy if you have to compile it yourself. Also, you may try to re-think the problem you're trying to solve. Maybe there's a different approach that is less computation-intensive. __ R-help@r-project.org mailto: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. -- Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation Biology, UCT), Dipl. Phys. (Germany) Centre of Excellence for Invasion Biology Stellenbosch University South Africa Tel : +33 - (0)9 53 10 27 44 Cell: +33 - (0)6 85 62 59 98 Fax : +33 - (0)9 58 10 27 44 Fax (D): +49 - (0)3 21 21 25 22 44 email: rai...@krugs.de mailto:rai...@krugs.de Skype: RMkrug -- Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation Biology, UCT), Dipl. Phys. (Germany) Centre of Excellence for Invasion Biology Stellenbosch University South Africa Tel : +33 - (0)9 53 10 27 44 Cell: +33 - (0)6 85 62 59 98 Fax : +33 - (0)9 58 10 27 44 Fax (D):+49 - (0)3 21 21 25 22 44 email: rai...@krugs.de Skype: RMkrug __ 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] using mclapply (multi core apply) to do matrix multiplication
7-02-2012, 02:31 (-0800); Alaios escriu: I would like to thank you Ernest for your answer. I guess that this is gonna be faster as right now R only sees one core. In my work there is a system with 64 cores and you can see only one working. If I understand it right a [m,n][n,k] matrix multiplication can be split into rows (from first matrice) and columns (from the second matrice) and then combine all the local results of each cpu together. Would that be too weird for mclapply to handle? I never used mclapply, but anyway here's a matrix multiplication function that uses lapply. Because the two lapply's are nested I don't think you can parallelize the two... I would only make the second one work with multiple cores mmult - function(a, b) { a - as.matrix(a) b - as.matrix(b) if (ncol(a) != nrow(b)) stop('non-conforming matrices') out - lapply(1:ncol(b), function(j) lapply(1:nrow(a), function(i) sum(a[i,] * b[,j]))) array(unlist(out), c(nrow(a), ncol(b))) } Also, I'm pretty sure that there are better algorithms. If you do this it would be interesting if you measured the execution time of the different alternatives and post the results :) -- Cheers, Ernest __ 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] box.cox
It means the author of that package has decided to start using a different function going forward but is still providing box.cox() for back-compatibility while the transition is made. Specifically, you should look up bcPower and replace box.cox in your scripts. For more information, type help(car-deprecated) at the prompt. Michael On Tue, Feb 7, 2012 at 1:47 AM, Rosario Garcia Gil m.rosario.gar...@slu.se wrote: Hello I am using box.cox() and I get this error message: Warning message: 'box.cox' is deprecated. Use 'bcPower' instead. See help(Deprecated) and help(car-deprecated). I went to help but I did not understand the explanation, I am still wondering what is really happening. Thanks /R __ 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] Lattice - different axis length
Dear all, I have a question about the lattice package, more specifically about the control of the x-axis length in the different panels. I use the following code to make the stacked barchart: barchart(country ~ climatechangefocalpoint + meteorologyservice + adaptationorvulnerability + cdmcarbonmarkets + energy + aviationmaritimetransport + forestry + pressofficer|period, data=graph5, as.table=T, xlim=c(0,150), layout=c(3,1), par.strip.text=list(cex=1.5), scales = list(alternating = 3, cex=1.2, tick.number=5), xlab=, col=c(grey15,grey75,grey30,grey90,grey45,grey0,grey60,grey100), #col=c(grey(100:1/100), grey(50:1/100), grey(0:1/100)), stack=T) Unfortunately, in the first two periods there are much less data, hence the bars are much shorter in those two panels and most of the space is unused, as the axis still run to 150 as for the last panel. Does anyone know how to cut the x-axis for the first two panels at, say, 50, while leaving it for the third panel as it is? Google and looking through old conversations here didn't help me, so I'm not quite sure whether this is possible at all. Thanks and best, Florian -- View this message in context: http://r.789695.n4.nabble.com/Lattice-different-axis-length-tp4364450p4364450.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] Positioning text in top left corner of plot
Regarding example 1: mtext has a line parameter that should help you example 2: using text(par(usr)[1],par(usr)[4],EUR billions,cex=0.7,adj=c(0,-1.5)) does not depend on the axis (values/range) beyond par('usr') With text(par(usr)[1], you get something closer to your second example, will however not write outside your plot for long texts Benno On Feb 7, 2012, at 12:06 PM, Manta wrote: Dear all, another questions related to zoo plotting. I would like to do as in the subject. Here a reproducible code: library(zoo) par(mfrow=c(2,1)) plot(zoo(seq(1:10),as.Date(seq(1:10),origin=1970-01-01)),xlab=,ylab=,main=Value,las=1) mtext(EUR billions,adj=0,cex=0.7) plot(zoo(seq(1:10),as.Date(seq(1:10),origin=1970-01-01)),xlab=,ylab=,main=Value,las=1) par(xpd=T) text(par(usr)[1],par(usr)[3]+10.5,EUR billions,cex=0.7) In the first graph I use the mtext function, which does the trick but it places the text too close too the y-axis. A second possibility is to use the text function, prior specification of the option parameter par(xpd=T). However, this does not look consistent as, depending on the graph, I would need to add different numbers in the y specification of the axis (par(usr)[3]+XXX). Any idea? -- View this message in context: http://r.789695.n4.nabble.com/Positioning-text-in-top-left-corner-of-plot-tp831723p4364355.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. Benno Pütz Statistical Genetics MPI of Psychiatry Kraepelinstr. 2-10 80804 Munich, Germany T: ++49-(0)89-306 22 222 F: ++49-(0)89-306 22 601 [[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] creating more vectors out of one
Thank you very much for your help! I used the split(apkz, cumsum(apkz==1)) command and it worked! I find the command very interesting and after a bit of thinking I found out, why the command cusum is used! (Because TRUE equals 1 and FALSE equals 0) Thanks again! Marion 2012/2/6 Petr PIKAL petr.pi...@precheza.cz Hi Dear R-helpers, I have got a vector which looks like the following: apkz - c(1,2,3,4,5,6,1,1,2,1,2,3,4) now I am trying to find a way to automatically create several vectors, each starting with the number 1, namely: First vector 1 2 3 4 5 6 Second vector 1 Thrid vector 1 2 Fourth vector 1 2 3 4 Does anyone know how to do that? This came to my mind as first, but I believe there are other options. split(apkz, cumsum(apkz==1)) Regards Petr Thank you very much for your help in advance! Marion [[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. [[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] using mclapply (multi core apply) to do matrix multiplication
7-02-2012, 03:32 (-0800); Alaios escriu: I wouldl ike to thank you for your response. The hardest part in the installation is to find a BLAS library to install. If I understand it right once I install BLAS then I only need to change a flag in the ./configure of R installation.. Our system is running opensuse and has intel cores. according to the link here http://cran.r-project.org/doc/manuals/R-admin.html#BLAS I ahve to find a proper BLAS library to insta.. In the explanation for the different alternatives seem that most of those are not implemented any more and other require special configuration :( This article includes an overview of different BLAS libraries along with benchmarks: http://cran.r-project.org/web/packages/gcbd/vignettes/gcbd.pdf It looks like using single-threaded ATLAS is already an improvement over LAPACK in most cases. I use Debian and it's straightforward to replace one with the other: you only have to install the libatlas3gf-base package and remove liblapack3gf and libblas3gf. Unfortunately, Debian does not include a multi-threaded version of ATLAS although they provide instructions on how to recompile the package yourself with multi-threading enabled. I don't know about SUSE, sorry. -- Bye, Ernest __ 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] xtable beta testing wanted
Hi Charles, This looks great, I'll try it out later today. Any chance that more automated multicolumn headers will be added in a future version? I know you can do it using the add.to.row argument, but it would be great to have an automated way using something like cgroup and n.cgroup in Hmisc:latex. Best, Ista On Tuesday, February 07, 2012 09:40:28 AM Charles Roosen wrote: Dear R-users, I've enhanced the xtable package, adding numerous features that have been requested by users. The changes are listed below. The objective throughout has been to avoid any breaking changes. However, as xtable is widely used and is a dependency of many packages I'd like to have others try it out before I post it to CRAN. Both bug reports and design feedback is welcomed. In particular, I'd like to have package authors that have written xtable methods try it out. If you are interested and willing, you can get the source from R-Forge via SVN: svn checkout svn://scm.r-forge.r-project.org/svnroot/xtable/ I can also send source and Windows binary versions to anyone that is interested. Feedback can be sent to me or posted in the issue tracking system on R-Forge: http://r-forge.r-project.org/tracker/?group_id=1228 Note that the R-Forge build system is currently being updated, so the option of installing directly from R-Forge won't get the latest version of the package. Thanks, Charlie Roosen xtable maintainer -- NEWS 1.7-0 (NOT YET RELEASED) * Added some vectorization code to improve performance. * Let caption be length 2, in which case the second value is the short caption used when creating a list of tables. * Added toLatex method. * Included print.xtable in the exported methods in the NAMESPACE file. * Added print.results argument to print that can be used to suppress the printing. * Added format.args argument to print that can be used to pass additional arguments such as big.marks to formatC(). * Added rotate.colnames and rotate.rownames arguments to print.xtable. * Added booktabs argument to use the \toprule, \midrule, and \bottomrule tags from the Latex booktabs package rather than using \hline for all horizontal lines. * Added scalebox argument to include a \scalebox clause around the tabular environment with the specified value used as the scaling factor. * Changed the print.xtable() arguments to use getOption() to check the options for a default value. This was suggested since print.xtable() has a lot of arguments that the user will typically leave unchanged between tables. * Added an is.null() check on the table.placement argument. * Added examples using the new arguments to the vignette. Charles Roosen, PhD Technical Director T: +41 (0)61 206 92 91 M: +41 (0)79 248 70 71 F: +41 (0) 61 206 92 99 www.mango-solutions.com Mango Solutions AG Aeschenvorstadt 36 4051 Basel Switzerland LEGAL NOTICE This message is intended for the use o...{{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. __ 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] Positioning text in top left corner of plot
Thanks, although I still have a couple of questions: 1. What is the line parameter? I could not find it in the manual... 2. How does exactly work the adj parameter when giving two different values? -- View this message in context: http://r.789695.n4.nabble.com/Positioning-text-in-top-left-corner-of-plot-tp831723p4364757.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] Lattice - different axis length
Hi, set scales=list(x=list(relation='free')) may help. Weidong Gu On Tue, Feb 7, 2012 at 6:42 AM, Florian Weiler fweile...@jhubc.it wrote: Dear all, I have a question about the lattice package, more specifically about the control of the x-axis length in the different panels. I use the following code to make the stacked barchart: barchart(country ~ climatechangefocalpoint + meteorologyservice + adaptationorvulnerability + cdmcarbonmarkets + energy + aviationmaritimetransport + forestry + pressofficer|period, data=graph5, as.table=T, xlim=c(0,150), layout=c(3,1), par.strip.text=list(cex=1.5), scales = list(alternating = 3, cex=1.2, tick.number=5), xlab=, col=c(grey15,grey75,grey30,grey90,grey45,grey0,grey60,grey100), #col=c(grey(100:1/100), grey(50:1/100), grey(0:1/100)), stack=T) Unfortunately, in the first two periods there are much less data, hence the bars are much shorter in those two panels and most of the space is unused, as the axis still run to 150 as for the last panel. Does anyone know how to cut the x-axis for the first two panels at, say, 50, while leaving it for the third panel as it is? Google and looking through old conversations here didn't help me, so I'm not quite sure whether this is possible at all. Thanks and best, Florian -- View this message in context: http://r.789695.n4.nabble.com/Lattice-different-axis-length-tp4364450p4364450.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@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] Table rearranging
On Feb 7, 2012, at 4:21 AM, Jeffrey Joh wrote: Thank you for your help, Bill. From the original table (not the plyr output), I would like to remove all the lines that do not have a corresponding open/closed measurement. For example, if there is a Closed yellow measurement on 0917, but not an Open yellow 0917 measurement, then the Closed yellow should be deleted. How can I make this change? In R you need to assign the results of a function to an object name so you code would look like: modified_data - ddply(d, .(date, color), summarize, meanClosed=mean(measurement[door==closed]), nClosed=sum(door==closed)) -- David Jeffrey From: wdun...@tibco.com To: johjeff...@hotmail.com; r-help@r-project.org Subject: RE: [R] Table rearranging Date: Tue, 7 Feb 2012 00:43:25 + Install and load the plyr package and try something like: ddply(d, .(date, color), summarize, + ddply(d, .(date, color), summarize + meanClosed=mean(measurement[door==closed]), nClosed=sum(door==closed)) date color meanOpen nOpen meanClosed nClosed 1 420 red 0.9741929 1 NaN 0 2 513 red 0.9352938 1 0.9620535 1 3 917 yellow NaN 0 0.9941022 1 4 1230 blue 0.9639099 1 0.9893108 1 5 1230 green 0.9765203 1 NaN 0 Bill Dunlap Spotfire, TIBCO Software wdunlap tibco.com -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org ] On Behalf Of Jeffrey Joh Sent: Monday, February 06, 2012 4:28 PM To: r-help@r-project.org Subject: [R] Table rearranging I have a table that looks like this: measurementdatedoorcolor 0.93529385513openred 0.97419293420openred 0.962053514513closedred 0.9639099371230openblue 0.976520341230opengreen 0.9893107951230closedblue 0.9941022917closedyellow I would like to create a table that has: Open measurement, Closed measurement, date, color. For every date/color combination, there should be two columns to represent the door open/closed measurement. If there are multiple datapoints with a given door/date/color combination, then they should be averaged. I would also like to make two columns to represent the number of datapoints that were averaged in determining the open/closed measurements. Jeffrey __ 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 Winsemius, MD West Hartford, CT __ 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] Lattice - different axis length
On Feb 7, 2012, at 6:42 AM, Florian Weiler wrote: Dear all, I have a question about the lattice package, more specifically about the control of the x-axis length in the different panels. I use the following code to make the stacked barchart: barchart(country ~ climatechangefocalpoint + meteorologyservice + adaptationorvulnerability + cdmcarbonmarkets + energy + aviationmaritimetransport + forestry + pressofficer|period, data=graph5, as.table=T, xlim=c(0,150), layout=c(3,1), par.strip.text=list(cex=1.5), scales = list(alternating = 3, cex=1.2, tick.number=5), xlab=, col = c (grey15 ,grey75,grey30,grey90,grey45,grey0,grey60,grey100), #col=c(grey(100:1/100), grey(50:1/100), grey(0:1/100)), stack=T) Unfortunately, in the first two periods there are much less data, hence the bars are much shorter in those two panels and most of the space is unused, as the axis still run to 150 as for the last panel. Does anyone know how to cut the x-axis for the first two panels at, say, 50, while leaving it for the third panel as it is? Google and looking through old conversations here didn't help me, so I'm not quite sure whether this is possible at all. in the scales list you could add:y=list(relation=free) ?xyplot and look at the relation sub-section in the scales section for more information. Thanks and best, Florian -- View this message in context: http://r.789695.n4.nabble.com/Lattice-different-axis-length-tp4364450p4364450.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. David Winsemius, MD West Hartford, CT __ 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] box.cox
Dear Tal and Rosario, The bcPower() function in the car package computes Box-Cox power transformations; boxcox() in the MASS package estimates the response transformation parameter for the Box-Cox regression model (actually the profile log-likelihood for a range of variables of the parameter). The two functions do different things. The boxCox() function in the car package is a slight generalization of boxcox(), allowing for other families of transformations than the Box-Cox powers; the powerTransform() function in the car package is more general still, in that it will handle multivariate transformations, both conditional (like the Box-Cox regression model) and unconditional. bcPower() was named box.cox() in earlier versions of the car package, associated with the first edition of the R [and S-PLUS] Companion to Applied Regression. In newer versions of the package, associated with the second edition of the book, written with Sandy Weisberg, this and other functions have been renamed to remove periods in function names. The older, deprecated, names of the functions are retained as a courtesy to readers of the first edition. Best, John John Fox Senator William McMaster Professor of Social Statistics Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of Tal Galili Sent: February-07-12 2:50 AM To: Rosario Garcia Gil Cc: r-help@R-project.org Subject: Re: [R] box.cox Use: library(MASS) ?boxcox Contact Details:--- Contact me: tal.gal...@gmail.com | 972-52-7275845 Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) | www.r- statistics.com (English) -- On Tue, Feb 7, 2012 at 8:47 AM, Rosario Garcia Gil m.rosario.gar...@slu.sewrote: Hello I am using box.cox() and I get this error message: Warning message: 'box.cox' is deprecated. Use 'bcPower' instead. See help(Deprecated) and help(car-deprecated). I went to help but I did not understand the explanation, I am still wondering what is really happening. Thanks /R __ 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-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] Positioning text in top left corner of plot
Both are described in the documentation (?(m)text works for me): Usage mtext(text, side = 3, line = 0, outer = FALSE, at = NA, adj = NA, padj = NA, cex = NA, col = NA, font = NA, ...) Arguments texta character or expression vector specifying the text to be written. Other objects are coerced by as.graphicsAnnot. sideon which side of the plot (1=bottom, 2=left, 3=top, 4=right). lineon which MARgin line, starting at 0 counting outwards. […] actually, line does accept fractional values , e.g. line=.5 should work for you and Usage text(x, ...) ## Default S3 method: text(x, y = NULL, labels = seq_along(x), adj = NULL, pos = NULL, offset = 0.5, vfont = NULL, cex = 1, col = NULL, font = NULL, ...) Arguments x, ynumeric vectors of coordinates where the text labels should be written. If the length of x and y differs, the shorter one is recycled. labels a character vector or expression specifying the text to be written. An attempt is made to coerce other language objects (names and calls) to expressions, and vectors and other classed objects to character vectors by as.character. If labels is longer than x and y, the coordinates are recycled to the length of labels. adj one or two values in [0, 1] which specify the x (and optionally y) adjustment of the labels. On most devices values outside that interval will also work. […] Hope the formatting comes out OK… Benno On Feb 7, 2012, at 2:49 PM, Manta wrote: Thanks, although I still have a couple of questions: 1. What is the line parameter? I could not find it in the manual... 2. How does exactly work the adj parameter when giving two different values? -- View this message in context: http://r.789695.n4.nabble.com/Positioning-text-in-top-left-corner-of-plot-tp831723p4364757.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. Benno Pütz Statistical Genetics MPI of Psychiatry Kraepelinstr. 2-10 80804 Munich, Germany T: ++49-(0)89-306 22 222 F: ++49-(0)89-306 22 601 __ 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] Logistic Regression
I'm surprised not to see the simple answer: glm models return the MLE estimate. fit - glm(y ~ x1 + x2 + , family='binomial') There is no need for special packages, this is a standard part of R. Terry Therneau begin included message -- On 02/06/2012 03:08 PM, Ana wrote: I am looking for R packages that can make a Logistic Regression model with parameter estimation by Maximum Likelihood Estimation. Many thanks for helping 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. Hi Ana, I worked out some maximum likelihood estimates for logistic regression in R on my website: http://www.netstorm.be/home/lrm Best regards Thomas __ 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] box.cox
I stand corrected - thank you for the clarification John. Contact Details:--- Contact me: tal.gal...@gmail.com | 972-52-7275845 Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) | www.r-statistics.com (English) -- On Tue, Feb 7, 2012 at 4:12 PM, John Fox j...@mcmaster.ca wrote: Dear Tal and Rosario, The bcPower() function in the car package computes Box-Cox power transformations; boxcox() in the MASS package estimates the response transformation parameter for the Box-Cox regression model (actually the profile log-likelihood for a range of variables of the parameter). The two functions do different things. The boxCox() function in the car package is a slight generalization of boxcox(), allowing for other families of transformations than the Box-Cox powers; the powerTransform() function in the car package is more general still, in that it will handle multivariate transformations, both conditional (like the Box-Cox regression model) and unconditional. bcPower() was named box.cox() in earlier versions of the car package, associated with the first edition of the R [and S-PLUS] Companion to Applied Regression. In newer versions of the package, associated with the second edition of the book, written with Sandy Weisberg, this and other functions have been renamed to remove periods in function names. The older, deprecated, names of the functions are retained as a courtesy to readers of the first edition. Best, John John Fox Senator William McMaster Professor of Social Statistics Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of Tal Galili Sent: February-07-12 2:50 AM To: Rosario Garcia Gil Cc: r-help@R-project.org Subject: Re: [R] box.cox Use: library(MASS) ?boxcox Contact Details:--- Contact me: tal.gal...@gmail.com | 972-52-7275845 Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) | www.r- statistics.com (English) -- On Tue, Feb 7, 2012 at 8:47 AM, Rosario Garcia Gil m.rosario.gar...@slu.sewrote: Hello I am using box.cox() and I get this error message: Warning message: 'box.cox' is deprecated. Use 'bcPower' instead. See help(Deprecated) and help(car-deprecated). I went to help but I did not understand the explanation, I am still wondering what is really happening. Thanks /R __ 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. [[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] na.action in stats::factanal()
Does factanal() force the user to use the formula interface if they wish to specify an na.action? Short answer: yes. Long answer: The handling of na.action is a built in part of the formula processing functions, so it's automatic when dealing with a formula. There are also downstream effects on predict() and resid() that are worked out for the formula case, but aren't clear otherwise. So- a. it would require extra programming and thought to work it out for matrix vector input, and the right answer isn't clear (it's harder than you might think). b. the usual assumption when a matrix/vector is given directly is the user knows what he's doing, or wouldn't have called it this way. For many routines, the matrix input is a speedup for simulations. c. factanal is unusual -- most routines split the two inputs. glm=formula interface glm.fit=matrix interface, lm lm.fit, coxph coxph.fit, 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.
[R] Weighted mad
Greetings UseRs, Pls advise if there is a way to write a func that can be supplied to aggregate to compute weighted MeanAbsolute Dev (MAD). I am having trouble passing the correct weights from each group level and cannot see the code behind aggregate. But maybe 'aggregate' is not the best way to do that. m1 - aggregate(pool[,c(SMM)],by=list(time=pool$ym),weighted.mean,w=pool$wght) Error in weighted.mean.default(X[[1L]], ...) : 'x' and 'w' must have the same length Apparently the grouping does not work on the additional argument. I am using weighted mean here just to be explicit and avoid supplying a custom function gor weighted MAD, which is not difficult to write by itself. It's making it work with aggreagte that is the problem. aggregate function (x, ...) UseMethod(aggregate) environment: namespace:stats Does not show anything... Stephen B __ 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] Lattice - different axis length
Thanks a lot for your help! I found that the relation=free command was very helpful, but then I had a little space between the axis and where the bars actually begun. I tried to deal with it in a panel function but was unable to do so. However, I found a way providing a list to xlim: barchart(country ~ new + cdmcarbonmarkets + energy + aviationmaritimetransport + forestry|period, data=graph5, as.table=T, layout=c(3,1), par.strip.text=list(cex=1.5), * xlim=list(c(.7,10),c(.8,12),c(8,130)), scales = list(alternating = 3, cex=1.2,x=list(relation=free)),* xlab=, col=c(grey25,grey75,grey0,grey100,grey50), stack=T) Oddly enough, I had to play around with the initial value to make the gap go away. The above code now produces exactly what I want (i.e. the bar start at 0 and have the right length). I do not really know why this is the case though. I would have expected the bars to be too short??? Still, I am happy with the result and thank you again for you help -- View this message in context: http://r.789695.n4.nabble.com/Lattice-different-axis-length-tp4364450p4364948.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] using mclapply (multi core apply) to do matrix multiplication
What is the nature of the matrices? Are they sparse or derived from sparse matrices? If they are sparse, have you looked at the packages available in R for sparse matrices? library(sos) summary(sp - findFn('sparse', 999)) will identify help pages in contributed packages containing sparse. The primary one is Matrix, but there are others. If they are not sparse but are derived from sparse matrices, you might be able to do some theoretical work. Of course, this only makes sense if you have a specific class of problems that generates the matrices, which seems plausible since you said you had square matrices of dimension 2^14. Hope this helps. Spencer On 2/7/2012 4:36 AM, Ernest Adrogué wrote: 7-02-2012, 03:32 (-0800); Alaios escriu: I wouldl ike to thank you for your response. The hardest part in the installation is to find a BLAS library to install. If I understand it right once I install BLAS then I only need to change a flag in the ./configure of R installation.. Our system is running opensuse and has intel cores. according to the link here http://cran.r-project.org/doc/manuals/R-admin.html#BLAS I ahve to find a proper BLAS library to insta.. In the explanation for the different alternatives seem that most of those are not implemented any more and other require special configuration :( This article includes an overview of different BLAS libraries along with benchmarks: http://cran.r-project.org/web/packages/gcbd/vignettes/gcbd.pdf It looks like using single-threaded ATLAS is already an improvement over LAPACK in most cases. I use Debian and it's straightforward to replace one with the other: you only have to install the libatlas3gf-base package and remove liblapack3gf and libblas3gf. Unfortunately, Debian does not include a multi-threaded version of ATLAS although they provide instructions on how to recompile the package yourself with multi-threading enabled. I don't know about SUSE, sorry. -- Spencer Graves, PE, PhD President and Chief Technology Officer Structure Inspection and Monitoring, Inc. 751 Emerson Ct. San José, CA 95126 ph: 408-655-4567 web: www.structuremonitoring.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] Weighted mad
You could try something like this. library(plyr) ddply(pool, .(ym), function(z){ weighted.mean(x= z$SMM, w = z$wght) }) ir. Thierry Onkelinx Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest team Biometrie Kwaliteitszorg / team Biometrics Quality Assurance Kliniekstraat 25 1070 Anderlecht Belgium + 32 2 525 02 51 + 32 54 43 61 85 thierry.onkel...@inbo.be www.inbo.be To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher The plural of anecdote is not data. ~ Roger Brinner The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. ~ John Tukey -Oorspronkelijk bericht- Van: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] Namens Bond, Stephen Verzonden: dinsdag 7 februari 2012 15:55 Aan: r-help@r-project.org Onderwerp: [R] Weighted mad Greetings UseRs, Pls advise if there is a way to write a func that can be supplied to aggregate to compute weighted MeanAbsolute Dev (MAD). I am having trouble passing the correct weights from each group level and cannot see the code behind aggregate. But maybe 'aggregate' is not the best way to do that. m1 - aggregate(pool[,c(SMM)],by=list(time=pool$ym),weighted.mean,w=pool$wght) Error in weighted.mean.default(X[[1L]], ...) : 'x' and 'w' must have the same length Apparently the grouping does not work on the additional argument. I am using weighted mean here just to be explicit and avoid supplying a custom function gor weighted MAD, which is not difficult to write by itself. It's making it work with aggreagte that is the problem. aggregate function (x, ...) UseMethod(aggregate) environment: namespace:stats Does not show anything... Stephen B __ 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] Positioning text in top left corner of plot
Thanks, mtext with option 'line' does the trick. I could not see that option in my version of R, but it appeared when I installed it again. -- View this message in context: http://r.789695.n4.nabble.com/Positioning-text-in-top-left-corner-of-plot-tp831723p4365088.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] Weighted mad
Dear Stephen The names of methods of aggregate can be viewed with the methods argument. Typing aggregate.data.frame shows you the method used for data frames. methods(aggregate) [1] aggregate.data.frame aggregate.defaultaggregate.formula* [4] aggregate.ts Non-visible functions are asterisked However, you need the function by, which passes a data frame to the function. pool - data.frame(SMM = rnorm(10), ym = rep(1:2, each = 5), wght = 1) pool SMM ym wght 1 0.904640493 11 2 1.021857174 11 3 1.234153020 11 4 -0.697625918 11 5 0.073146470 11 6 1.438553786 21 7 -0.83118 21 8 -0.055872825 21 9 0.858622161 21 10 0.001968816 21 mad - by(pool, pool$ym, + function(pool) { weighted.mean(x = pool$SMM, w = pool$wght) }) mad[1:2] pool$ym 1 2 0.5072342 0.2824177 Regards, Chris Campbell MANGO SOLUTIONS Data Analysis that Delivers +44 1249 767700 -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Bond, Stephen Sent: 07 February 2012 14:55 To: r-help@r-project.org Subject: [R] Weighted mad Greetings UseRs, Pls advise if there is a way to write a func that can be supplied to aggregate to compute weighted MeanAbsolute Dev (MAD). I am having trouble passing the correct weights from each group level and cannot see the code behind aggregate. But maybe 'aggregate' is not the best way to do that. m1 - aggregate(pool[,c(SMM)],by=list(time=pool$ym),weighted.mean,w=pool$wght) Error in weighted.mean.default(X[[1L]], ...) : 'x' and 'w' must have the same length Apparently the grouping does not work on the additional argument. I am using weighted mean here just to be explicit and avoid supplying a custom function gor weighted MAD, which is not difficult to write by itself. It's making it work with aggreagte that is the problem. aggregate function (x, ...) UseMethod(aggregate) environment: namespace:stats Does not show anything... Stephen B __ 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. Kind Regards Chris Campbell [Mango Solutions] T:+44 (0)1249 767700 F: +44 (0)1249 767707 Mango Solutions Websitehttp://www.mango-solutions.com Unit 2 Greenways Business Park Bellinger Close Chippenham Wilts SN15 1BN UK With effect from Friday 24th February 2012 Mango Business Solutions Ltd main office will be moving to: 2 Methuen Park Chippenham Wilts SN14 0GB UK Please note, there may be some disruption to services; including emails, from 1730 hrs on Thursday 23rd February. For any urgent enquiries, please use the following mobile numbers: 07966062462, 07900580808, 07967808091 LEGAL NOTICE This message is intended for the use o...{{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.
[R] Best way to remove all objects but leave the functions in a workspace.
I think the subject says it all. Thanks in advance, KW -- [[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] using mclapply (multi core apply) to do matrix multiplication
I wouldl ike to thank you for your response. The hardest part in the installation is to find a BLAS library to install. If I understand it right once I install BLAS then I only need to change a flag in the ./configure of R installation.. Our system is running opensuse and has intel cores. according to the link here http://cran.r-project.org/doc/manuals/R-admin.html#BLAS I ahve to find a proper BLAS library to insta.. In the explanation for the different alternatives seem that most of those are not implemented any more and other require special configuration :( Are there not any rpm package to do the work nice and transparently? B.R Alex From: Rainer M Krug r.m.k...@gmail.com Cc: Ernest Adrogué nfdi...@gmail.com; r-help@r-project.org r-help@r-project.org Sent: Tuesday, February 7, 2012 12:08 PM Subject: Re: [R] using mclapply (multi core apply) to do matrix multiplication On 07/02/12 12:02, Alaios wrote: Thank you very much for your point... I hope I can find some easy to follow instructions as I do not have root permission Me neither on pur cluster - but that won't stop you from compiling and installing R in your home directory. By doing this, you have even more control. Cheers and good luck, Rainer for the many cores system and our system administrator want to have easy instructions to follow. Thanks a gain. *From:* Rainer M Krug r.m.k...@gmail.com *Cc:* Ernest Adrogué nfdi...@gmail.com; r-help@r-project.org r-help@r-project.org *Sent:* Tuesday, February 7, 2012 11:44 AM *Subject:* Re: [R] using mclapply (multi core apply) to do matrix multiplication On 07/02/12 11:31, Alaios wrote: I would like to thank you Ernest for your answer. I guess that this is gonna be faster as right now R only sees one core. In my work there is a system with 64 cores and you can see only one working. If I understand it right a [m,n][n,k] matrix multiplication can be split into rows (from first matrice) and columns (from the second matrice) and then combine all the local results of each cpu together. You definitaly can go this way, but I would STRONGLY recommend to search for parallel BLAS, check in the R-admin manual the section Linear Algebra which deals with BLAS et al, and e.g. http://www.r-bloggers.com/compiling-64-bit-r-2-10-1-with-mkl-in-linux/ My guess is that a paralelization on the C level in the BLAS et al. library will be MUCH faster then a paralelization on R level. Also, there is a R-sig-hpc mailing list for these kind of questions. Cheers, Rainer Would that be too weird for mclapply to handle? B.R Alex From: Ernest Adroguénfdi...@gmail.com mailto:nfdi...@gmail.com To: r-help@r-project.org mailto:r-help@r-project.org Sent: Tuesday, February 7, 2012 11:02 AM Subject: Re: [R] using mclapply (multi core apply) to do matrix multiplication 7-02-2012, 00:29 (-0800); Alaios escriu: Dear all, I am trying to multiply three different matrices and each matrice is of size 16384,16384 the normal %*% multiplciation operator has not finished one day now. As I am running a system with many cores (and it seems that R is using only one of those) I would like to write fast a brief function that converts the typical for loops of a matrix multiplication to a set of lapply sets (mclapply uses the lapply syntax but it applies the work to many cores). If my thinking is correct , in the sense that this will speed up things a lot, I want you to help me covert the first matrix in rows the second in columns and convert those in a format that lapply would like to work with. If I understand correctly, R uses a specialized library called BLAS to do matrix multiplications. I doubt re-implementing the matrix multiplication code at R-level would be any faster. What you can try is replace BLAS with a multicore version of BLAS although it's not easy if you have to compile it yourself. Also, you may try to re-think the problem you're trying to solve. Maybe there's a different approach that is less computation-intensive. __ R-help@r-project.org mailto: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. -- Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation Biology, UCT), Dipl. Phys. (Germany) Centre of Excellence for Invasion Biology Stellenbosch University South Africa Tel : +33 - (0)9 53 10 27 44 Cell: +33 - (0)6 85 62 59 98 Fax : +33 - (0)9 58 10 27 44 Fax (D): +49 - (0)3 21 21 25 22 44 email: rai...@krugs.de mailto:rai...@krugs.de Skype: RMkrug -- Rainer M.
[R] fixed effects with clustered standard errors
Dear R-helpers, I have a very simple question and I really hope that someone could help me I would like to estimate a simple fixed effect regression model with clustered standard errors by individuals. For those using Stata, the counterpart would be xtreg with the fe option, or areg with the absorb option and in both case the clustering is achieved with vce(cluster id) My question is : how could I do that with R ? An important point is that I have too many individuals, therefore I cannot include dummies and should use the demeaning usual procedure. I tried with the plm package with the within option, but R quikcly tells me that the memory limits are attained (I have over 10go ram!) while the dataset is only 700mo (about 50 000 individuals, highly unbalanced) I dont understand... plm do indeed demean the data so the computation should be fast and light enough... ?! Are there any other solutions ? Many thanks in advance ! ;) John __ 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] GLM Quasibinomial - 48 models
I've originally made 48 GLM binomial models and compare the AIC values. But dispersion was very large: Example: Residual deviance: 8811.6 on 118 degrees of freedom I was suggested to do a quasibinomial afterwards but found that it did not help the dispersion factor of models and received a warning: Residual deviance: 3005.7 on 67 degrees of freedom AIC: NA Number of Fisher Scoring iterations: 13 Warning message: In summary.glm(qModel48.glm) : observations with zero weight not used for calculating dispersion 1) What does this warning message mean? 2) Is this quasibinomial step necessary or can I just compare the AIC of 48 glm models binomial even though dispersion is very large, well over 2. Jean -- View this message in context: http://r.789695.n4.nabble.com/GLM-Quasibinomial-48-models-tp4364137p4364137.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] Writing to a file
Thanks a lot for the interest :) My loop is the following counter = 0 for (i in 1:nrow(y)) { for (j in 1:ncol(y)) { if (y[i,j]==Func_0005634) { counter = counter + 1 } if(y[i,j]==Func_0005737){ counter = counter + 1 } if(y[i,j]==Func_0005515){ counter = counter + 1 } } if(counter == 3 ){ cat(y[i,1], file = foo.csv, \n) } counter = 0 } and after read.table(foo.csv) I get V1 1 45 which is the last result why does it overwrite? how can I have all the results? Eager to a reply from you! -- View this message in context: http://r.789695.n4.nabble.com/Writing-to-a-file-tp3070617p4364149.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] glht (multicomparisons) with a binomial response variable
Thanks very much, Frank, I will try, first I have to adjust your solution to my problem cause i have just for months : april, may, june and july. I will tell you!! Thanks again - Mario Garrido Escudero PhD student Dpto. de Biología Animal, Ecología, Parasitología, Edafología y Qca. Agrícola Universidad de Salamanca -- View this message in context: http://r.789695.n4.nabble.com/glht-multicomparisons-with-a-binomial-response-variable-tp4360898p4364501.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] Actual vs. predicted survival times
Dear R-help, I am using R 2.14.1 on Windows 7. I would like to produce a plot like the attached - although simplified to actual vs. Predicted survival time with distinguishing marks for censored and observed points. I have a dataset and have fitted a Cox model to it. In an attempt to visualise how accurate the model is it would be ideal if I could plot the actual survival times against the predicted survival times. I have been looking on the internet to see if there are ways to do this in R. The only post I found (https://stat.ethz.ch/pipermail/r-help/2009-February/189888.html) that seemed directly relevant suggested that I shouldn't be generating survival times at all. Given that, I was concerned about proceeding but I would like to have access to a plot to make a decision on its usefulness. I appreciate that there are predict.coxph and predict.cph options available to me. My first attempt was as follows: # fit Cox model # coxfita = coxph(Surv(tsecond,seccens)~stroke(smess)+rels(smess)+asleep(smess)+eeg1(smess)+eeg2(smess)+ct1(smess)+ct2(smess)+treat(smess),data=smess) # Find censored and observed groups # messcens - subset(smess,seccens==1) messobs - subset(smess,seccens==0) # Obtain predicted survival times # explp - exp(predict(coxfita,type=lp)) explp2 - mean(ssmess$tsecond,na.rm=TRUE)*explp smess2 - data.frame(ssmess,explp2) # Find censored and observed groups # smesscens - subset(smess2,seccens==1) smessobs - subset(smess2,seccens==0) # Produce plot # plot(smesscens$explp2,messcens$tsecond,pch=4,col=blue,ylab=Actual Survival Time,xlab=Predicted Survival Time,main=Survival Times,xlim=c(0,3500),ylim=c(0,3500)) points(smessobs$explp2,messobs$tsecond,pch=4,col=red) This leads to the attached plot. It doesn't seem correct though as the predicted times all start over 500 days. Any suggestions would be very welcome. Many thanks, Laura Actual vs. Survival LJB.pdf Description: Actual vs. Survival LJB.pdf __ 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] nice report generator?
Hi Michael, Hi all, an alternative for reading/writing tables from/to Excel is the package XLConnect. It is platform-independent and therefore runs under Windows, UNIX/Linux Mac. It supports reading/writing worksheets as well as named ranges. In addition, for reporting purposes, there is support for cell styles. You can find an example here: http://miraisolutions.wordpress.com/2011/08/31/xlconnect-a-platform-independent-interface-to-excel/ XLConnect – A platform-independent interface to Excel . With respect to cell styles: you either completely code up the styling according to your needs using /setCellStyle/ (which can be cumbersome if you need to apply a lot of styling) or you can use so-called style actions controlled via /setStyleAction/. See the XLConnect reference manual for more information. Best regards, Martin -- View this message in context: http://r.789695.n4.nabble.com/nice-report-generator-tp4169939p4365144.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] binomial vs quasibinomial
After looking at 48 glm binomial models I decided to try the quasibinomial with the top model 25 (lowest AIC). To try to account for overdispersion (residual deviance 2679.7/68 d.f.) After doing so the dispersion factor is the same for the quasibinomial and less sectors of the beach were significant by p-value. While the p-values in the binomial were more significant for each section of the beach. -- telling me more about the beach. Is this ok? Can I just look at the binomial glm model 25 and look at its p-values for beach sections and forget about the quasibinomial model 25? J Call: glm(formula = cbind(Shells, TotalEggs - Shells) ~ Sector:Veg:Aeventexhumed, family = quasibinomial, data = data.to.analyze) -- View this message in context: http://r.789695.n4.nabble.com/binomial-vs-quasibinomial-tp4364371p4364371.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] Writing to a file
As I said to you a while back, use append = TRUE. Michael On Tue, Feb 7, 2012 at 4:18 AM, Felicity felicity...@hotmail.com wrote: Thanks a lot for the interest :) My loop is the following counter = 0 for (i in 1:nrow(y)) { for (j in 1:ncol(y)) { if (y[i,j]==Func_0005634) { counter = counter + 1 } if(y[i,j]==Func_0005737){ counter = counter + 1 } if(y[i,j]==Func_0005515){ counter = counter + 1 } } if(counter == 3 ){ cat(y[i,1], file = foo.csv, \n) } counter = 0 } and after read.table(foo.csv) I get V1 1 45 which is the last result why does it overwrite? how can I have all the results? Eager to a reply from you! -- View this message in context: http://r.789695.n4.nabble.com/Writing-to-a-file-tp3070617p4364149.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@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] Creating time series (ts) object
Thank you for your help. I'll simply remove leap days from my serie. Regards, Phil -- View this message in context: http://r.789695.n4.nabble.com/Creating-time-series-ts-object-tp4362762p4364664.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] Problem with plotting a square 1 x 3 plot and placement of outer margin text
Dear R-helpers, Please see the attached plot. The problem is that I have too much space between the x-axis label (which is mtext in an outer margin) and the plots. My par settings for this plot are: par(mfrow=c(1,3),oma=c(2,2,2,2),mar=c(5.1,4.1,4.1,2.1),pty=s) #here is the code that produces the three plots, which I have deleted for simplicity mtext(Log Wetland Area,side=1,outer=TRUE) It works fine (less space between plots and outer margin text)) when I set pty=m but then I get very long and skinny rectangular plots. I would like to keep the square plots. Any help would be much appreciated! Many thanks, Mark __ 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] R2WinBUGS error message
Hi! I am new to BUGS and running BUGS from R. I am trying to run a regression model from R, however I have this error message: Error in file(con, wb) : cannot open the connection In addition: Warning messages: 1: In file.create(to[okay]) : cannot create file 'c:/Program Files/WinBUGS14//System/Rsrc/Registry_Rsave.odc', reason 'Permission denied' 2: In file(con, wb) : cannot open file 'c:/Program Files/WinBUGS14//System/Rsrc/Registry.odc': Permission denied Can anyone help me out ? I am running R2WinBUGS from a windows 7 machine. -- View this message in context: http://r.789695.n4.nabble.com/R2WinBUGS-error-message-tp4364838p4364838.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] Vectorizing a loop
Hello Folks, I'm trying to vectorize a loop that processes rows of a dataframe. It involves lots of conditionals, such as If column 10 == 3, and if column 3 is True, and both column 5 and 6 are False, then set column 4 to True. So, for example, any ideas about vectorizing the following? df = data.frame( list(a=c(1,2,3,4), b=c(a,b,c,d), c=c(T,F,T,F), d=NA, e=c(F,F,T,T)) ) for (i in 1:nrow(df)) { if (df[i,3] %in% c(FALSE,NA) (df[i,1] 2 | df[i,5]) ) { df[i,4] = 1 } if (df[i,5] %in% c(TRUE, NA) df[i,2] == b) { df[i,4] = 2 df[i,5] = T } } Thanks, Allie __ 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] Reading a Fixed-Width File with SAS import instructions directly into an R data frame
Hi everyone, I'm wondering if anyone has written (or knows of) an R function that takes the SAS import code to read in an ASCII / fixed-width data file.. and then parses out the SAS code to figure out how to structure a (foreign package) read.fwf command so that fixed-width data file can be read directly into an R data frame? I'm envisioning a function that returns an R data frame from ASCII data set, using two parameters: the location of the ASCII data and the location of the SAS import code. I have done some work on the topic (https://gist.github.com/1253937), but it's got a lot of bugs. Before I spend a bunch of time cleaning this up, I wanted to check that I'm not re-inventing the wheel. Any leads on related (and probably much further along) work would be appreciated. Thanks!! [[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] Actual vs. predicted survival times
On Feb 7, 2012, at 10:32 AM, Bonnett, Laura wrote: Dear R-help, I am using R 2.14.1 on Windows 7. I would like to produce a plot like the attached - although simplified to actual vs. Predicted survival time with distinguishing marks for censored and observed points. I have a dataset and have fitted a Cox model to it. In an attempt to visualise how accurate the model is it would be ideal if I could plot the actual survival times against the predicted survival times. I have been looking on the internet to see if there are ways to do this in R. The only post I found (https://stat.ethz.ch/pipermail/r-help/2009-February/189888.html ) that seemed directly relevant suggested that I shouldn't be generating survival times at all. Given that, I was concerned about proceeding but I would like to have access to a plot to make a decision on its usefulness. I appreciate that there are predict.coxph and predict.cph options available to me. My first attempt was as follows: # fit Cox model # coxfita = coxph(Surv(tsecond,seccens)~stroke(smess)+rels(smess) +asleep(smess)+eeg1(smess)+eeg2(smess)+ct1(smess)+ct2(smess) +treat(smess),data=smess) # Find censored and observed groups # messcens - subset(smess,seccens==1) messobs - subset(smess,seccens==0) # Obtain predicted survival times # explp - exp(predict(coxfita,type=lp)) That gives you relative risks, not survival times. explp2 - mean(ssmess$tsecond,na.rm=TRUE)*explp Why are you multiplying times by relative risks? That makes no sense. -- David. smess2 - data.frame(ssmess,explp2) # Find censored and observed groups # smesscens - subset(smess2,seccens==1) smessobs - subset(smess2,seccens==0) # Produce plot # plot(smesscens$explp2,messcens$tsecond,pch=4,col=blue,ylab=Actual Survival Time,xlab=Predicted Survival Time,main=Survival Times,xlim=c(0,3500),ylim=c(0,3500)) points(smessobs$explp2,messobs$tsecond,pch=4,col=red) This leads to the attached plot. It doesn't seem correct though as the predicted times all start over 500 days. Any suggestions would be very welcome. Many thanks, Laura Actual vs. Survival LJB.pdf__ 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 Winsemius, MD West Hartford, CT __ 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] lmer with spatial and temporal random factors, not nested
Marte Lilleeng mlilleeng at gmail.com writes: Hi, I am new to this list. The r-sig-mixed-mod...@r-project.org mailing list would be more appropriate for this question -- please direct any further questions there ... I have a question regarding including both spatial and temporal random factors in lmer. These two are not nested, and an example of model I try to fit is model1-lmer(Richness~Y+Canopy+Veg_cm+Treatment+(1|Site/Block/Plot)+ (1|Year), family=poisson, REML=FALSE), where richness = integer Y Treatment = factor Canopy Veg_cm = numerical, continous Site/Block/Plot= factor Year = integer Fine, but REML=FALSE is unnecessary/irrelevant for generalized linear mixed model (family!=gaussian) fits. I get the following warning message: Warning messages: 1: In mer_finalize(ans) : Cholmod warning 'not positive definite' at file:../Cholesky/t_cholmod_rowfac.c, line 432 2: In mer_finalize(ans) : singular convergence (7) Is this due to the nature of my fixed/random factors or the way I put up the random factors? Hard to tell exactly. It's probably due to overfitting and/or lack of balance (glmer handles lack of balance, but extreme lack of balance can lead to technical difficulties like this one). In lme I could include a component for autocorrelation, ex:cor=corAR1(form=~Year|Site/Block/ID). Does the equivalent exist for lmer? No, sorry. Crossed random effects are possible in lme (see p. 165?) of Pinheiro and Bates 2000, and glmmPQL in the MASS package can handle a Poisson response, so that might be the best way to go. However, I would also strongly encourage you to do some graphical exploration of your data and make sure there aren't outliers, almost-empty blocks, etc. Ben Bolker __ 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] na.action in stats::factanal()
Thanks for the clear explanation Terry! It gets ugly for many factanal applications, where you are dealing with 300 variables… One question: what would be wrong with auto generating the formula from a matrix call? That way the matrix call gets the benefit of returning scores. Also: you say a person using the non-formula version might know what they are doing. My guess is that is not the case; Most non-experts do just this. And an expert could still not get scores back, no? best, tim PS: If anyone who cares about documentation is reading, it would be lovely to include a valid example for getting scores in a realistic dataset with NAs… where the na.action has to be set. On 7 Feb 2012, at 2:40 PM, Terry Therneau wrote: Does factanal() force the user to use the formula interface if they wish to specify an na.action? Short answer: yes. Long answer: The handling of na.action is a built in part of the formula processing functions, so it's automatic when dealing with a formula. There are also downstream effects on predict() and resid() that are worked out for the formula case, but aren't clear otherwise. So- a. it would require extra programming and thought to work it out for matrix vector input, and the right answer isn't clear (it's harder than you might think). b. the usual assumption when a matrix/vector is given directly is the user knows what he's doing, or wouldn't have called it this way. For many routines, the matrix input is a speedup for simulations. c. factanal is unusual -- most routines split the two inputs. glm=formula interface glm.fit=matrix interface, lm lm.fit, coxph coxph.fit, 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.
[R] predict.naiveBayes() bug in e1071 package
Hi, I'm currently using the R package e1071 to train naive bayes classifiers and came across a bug: When the posterior probabilities of all classes are small, the result from the predict.naiveBayes function become NaNs. This is an issue with the treatment of the log-transformed probabilities inside the predict.naiveBayes function. Here is an example to demonstrate the problem (you might need to increase 'nvar' depending on your machine): 8 N - 100 nvar - 60 varnames - paste(v, 1:nvar, sep=) dat - sapply(1:nvar, function(dummy) {c(rnorm(N/2, 0, 1), rnorm(N/2, 10, 1))}) colnames(dat) - varnames out - rep(c(a,b), each=N/2) names(dat) - varnames nb - naiveBayes(x=dat, y=out) new.dat - t(rnorm(nvar, 5, 0.1)) colnames(new.dat) - varnames predict(nb, new.dat, type=raw) 8 the results of the last line is usually NaNs. As for the solution: To protect agains very small numbers, the e1071:::predict.naiveBayes function takes the probabilities into log-space and adds instead of multiplying probabilities. However, when calculating the posterior probabilities of each class (when type = raw), the log of the probabilities are exponentiated, which defeats the purpose of the logspace transformation. I suggest the following change to the code: Towards the end of the predict.naiveBayes function, you currently do: L - exp(L) L / sum(L) # this is what is returned you can instead use sapply(L, function(lp) {1 / sum(exp(L - lp))}) the above comes from the following equality: x / (x + y + z) = 1 / (1 + exp(log(y) - log(x)) + exp(log(z) - log(x))) Best wishes, /Ali Tofigh __ 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] Vectorizing a loop
On Tue, 7 Feb 2012, Alexander Shenkin wrote: Hello Folks, I'm trying to vectorize a loop that processes rows of a dataframe. It involves lots of conditionals, such as If column 10 == 3, and if column 3 is True, and both column 5 and 6 are False, then set column 4 to True. So, for example, any ideas about vectorizing the following? df = data.frame( list(a=c(1,2,3,4), b=c(a,b,c,d), c=c(T,F,T,F), d=NA, e=c(F,F,T,T)) ) for (i in 1:nrow(df)) { if (df[i,3] %in% c(FALSE,NA) (df[i,1] 2 | df[i,5]) ) { df[i,4] = 1 } if (df[i,5] %in% c(TRUE, NA) df[i,2] == b) { df[i,4] = 2 df[i,5] = T } } Your code attempts to do some things with NA that won't behave the way you expect them to. Specifically, you cannot use %in% to test for NA, and you cannot give the if function an NA. It only appears to work because you don't actually give it a complete set of test values consistent with your tests in the loop. My guess at your intent is: df - data.frame( list( a=c(1,2,3,4,5) , b=c(a,b,c,d,e) , c=c(TRUE,FALSE,TRUE,FALSE,NA) , d=NA , e=c(FALSE,FALSE,TRUE,TRUE,NA) ) ) tmpdf - df for (i in 1:nrow(df)) { if ( ( is.na(df[i,3]) || !df[i,3] ) ( df[i,1] 2 || ( is.na( df[i,5] ) || df[i,5] ) ) ) { df[i,4] - 1 } if ( ( is.na( df[i,5] ) || df[i,5] ) df[i,2] == b ) { df[i,4] - 2 df[i,5] - TRUE } } df2 - df df - tmpdf # intermediate logical vectors for clarity tmp - ( is.na(df[[3]]) | !df[[3]] ) ( df[[1]] 2 | df[[5]] ) tmp2 - ( is.na(df[[5]]) | df[[5]] ) df[[2]] == b df[ tmp, d ] - 1 df[ tmp2, d ] - 2 df[ tmp2, e ] - TRUE --- Jeff NewmillerThe . . Go Live... DCN:jdnew...@dcn.davis.ca.usBasics: ##.#. ##.#. Live Go... Live: OO#.. Dead: OO#.. Playing Research Engineer (Solar/BatteriesO.O#. #.O#. with /Software/Embedded Controllers) .OO#. .OO#. rocks...1k __ 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] Vectorizing a loop
On Tue, Feb 07, 2012 at 11:39:42AM -0500, Alexander Shenkin wrote: Hello Folks, I'm trying to vectorize a loop that processes rows of a dataframe. It involves lots of conditionals, such as If column 10 == 3, and if column 3 is True, and both column 5 and 6 are False, then set column 4 to True. So, for example, any ideas about vectorizing the following? df = data.frame( list(a=c(1,2,3,4), b=c(a,b,c,d), c=c(T,F,T,F), d=NA, e=c(F,F,T,T)) ) for (i in 1:nrow(df)) { if (df[i,3] %in% c(FALSE,NA) (df[i,1] 2 | df[i,5]) ) { df[i,4] = 1 } if (df[i,5] %in% c(TRUE, NA) df[i,2] == b) { df[i,4] = 2 df[i,5] = T } } Hi. Try the following. cond1 - (df[,3] %in% c(FALSE,NA)) (df[,1] 2 | df[,5]) df[,4] - ifelse(cond1, 1, df[,4]) cond2 - (df[,5] %in% c(TRUE, NA)) (df[,2] == b) df[,4] - ifelse(cond2, 2, df[,4]) df[,5] - ifelse(cond2, TRUE, df[,5]) Hope this helps. Petr Savicky. __ 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] Vectorizing a loop
On Feb 7, 2012, at 12:56 PM, Jeff Newmiller wrote: On Tue, 7 Feb 2012, Alexander Shenkin wrote: Hello Folks, I'm trying to vectorize a loop that processes rows of a dataframe. It involves lots of conditionals, such as If column 10 == 3, and if column 3 is True, and both column 5 and 6 are False, then set column 4 to True. So, for example, any ideas about vectorizing the following? df = data.frame( list(a=c(1,2,3,4), b=c(a,b,c,d), c=c(T,F,T,F), d=NA, e=c(F,F,T,T)) ) for (i in 1:nrow(df)) { if (df[i,3] %in% c(FALSE,NA) (df[i,1] 2 | df[i,5]) ) { df[i,4] = 1 } if (df[i,5] %in% c(TRUE, NA) df[i,2] == b) { df[i,4] = 2 df[i,5] = T } } Your code attempts to do some things with NA that won't behave the way you expect them to. Specifically, you cannot use %in% to test for NA, Huh? NA %in% NA [1] TRUE NA %in% c(5, NA) [1] TRUE NA %in% c(5, 6) [1] FALSE -- David. and you cannot give the if function an NA. It only appears to work because you don't actually give it a complete set of test values consistent with your tests in the loop. My guess at your intent is: df - data.frame( list( a=c(1,2,3,4,5) , b=c(a,b,c,d,e) , c=c(TRUE,FALSE,TRUE,FALSE,NA) , d=NA , e=c(FALSE,FALSE,TRUE,TRUE,NA) ) ) tmpdf - df for (i in 1:nrow(df)) { if ( ( is.na(df[i,3]) || !df[i,3] ) ( df[i,1] 2 || ( is.na( df[i,5] ) || df[i,5] ) ) ) { df[i,4] - 1 } if ( ( is.na( df[i,5] ) || df[i,5] ) df[i,2] == b ) { df[i,4] - 2 df[i,5] - TRUE } } df2 - df df - tmpdf # intermediate logical vectors for clarity tmp - ( is.na(df[[3]]) | !df[[3]] ) ( df[[1]] 2 | df[[5]] ) tmp2 - ( is.na(df[[5]]) | df[[5]] ) df[[2]] == b df[ tmp, d ] - 1 df[ tmp2, d ] - 2 df[ tmp2, e ] - TRUE --- Jeff NewmillerThe . . Go Live... DCN:jdnew...@dcn.davis.ca.usBasics: ##.#. ##.#. Live Go... Live: OO#.. Dead: OO#.. Playing Research Engineer (Solar/BatteriesO.O#. #.O#. with /Software/Embedded Controllers) .OO#. .OO#. rocks...1k __ 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 Winsemius, MD West Hartford, CT __ 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] predict.naiveBayes() bug in e1071 package
On Feb 7, 2012, at 12:43 PM, Ali Tofigh wrote: Hi, I'm currently using the R package e1071 to train naive bayes classifiers and came across a bug: When the posterior probabilities of all classes are small, the result from the predict.naiveBayes function become NaNs. This should be sent to the maintainer of the package. The name of the maintainer can always be found in the DESCRIPTION file. Several of the authors are regular readers of rhelp, but I do not know whether David Meyer is. I'm sure a well-documented bug report, as this appears to be, will be welcomed. -- David. This is an issue with the treatment of the log-transformed probabilities inside the predict.naiveBayes function. Here is an example to demonstrate the problem (you might need to increase 'nvar' depending on your machine): 8 N - 100 nvar - 60 varnames - paste(v, 1:nvar, sep=) dat - sapply(1:nvar, function(dummy) {c(rnorm(N/2, 0, 1), rnorm(N/ 2, 10, 1))}) colnames(dat) - varnames out - rep(c(a,b), each=N/2) names(dat) - varnames nb - naiveBayes(x=dat, y=out) new.dat - t(rnorm(nvar, 5, 0.1)) colnames(new.dat) - varnames predict(nb, new.dat, type=raw) 8 the results of the last line is usually NaNs. As for the solution: To protect agains very small numbers, the e1071:::predict.naiveBayes function takes the probabilities into log-space and adds instead of multiplying probabilities. However, when calculating the posterior probabilities of each class (when type = raw), the log of the probabilities are exponentiated, which defeats the purpose of the logspace transformation. I suggest the following change to the code: Towards the end of the predict.naiveBayes function, you currently do: L - exp(L) L / sum(L) # this is what is returned you can instead use sapply(L, function(lp) {1 / sum(exp(L - lp))}) the above comes from the following equality: x / (x + y + z) = 1 / (1 + exp(log(y) - log(x)) + exp(log(z) - log(x))) Best wishes, /Ali Tofigh __ 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 Winsemius, MD West Hartford, CT __ 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] Best way to remove all objects but leave the functions in a workspace.
Best? Probably not. No money back if this deletes everything by mistake: remove(list=ls()[sapply(ls(),function(n){!is.function(get(n))})]) As a function, maybe: clearNonF=function(e=.GlobalEnv){remove(list=ls(e)[sapply(ls(e),function(n){!is.function(get(n))})],envir=e)} clearNonF() will clear out the current workspace. Maybe. On Tue, Feb 7, 2012 at 4:13 PM, Keith Weintraub kw1...@gmail.com wrote: I think the subject says it all. Barry __ 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] UTF8
Hi, I just found the solution to the problem in the following posts: https://stat.ethz.ch/pipermail/r-sig-mac/2008-October/005395.html AND https://stat.ethz.ch/pipermail/r-sig-mac/2010-February/007106.html Now everything works. cheers Norma -- View this message in context: http://r.789695.n4.nabble.com/UTF8-tp873280p4365624.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] .Random.seed not found
Hi All, I have a user who is receiving this error after running the following: library(rjags) library(R2jags) x = c(1.6907, 1.7242, 1.7552, 1.7842, 1.8113, 1.8369, 1.8610, 1.8839) n = c(59, 60, 62, 56, 63, 59, 62, 60) r = c(6, 13, 18, 28, 52, 53, 61, 60) N = 8 data=list(x,n,r,N) inits=function(){ list(alpha.star=0,beta=0) } parameters=c(alpha,beta,rhat) result=jags(data,inits,model.file=test.bug,parameters,n.chains=3,n.iter=1,n.burnin=1000,n.thin=1) module glm loaded Error: object '.Random.seed not found I don't think it has anything to do with rjags or R2jags because when I run the same code on the same system, I receive no error messages. The R version is 2.14.0 on a Linux 64 machine. Any help would be appreciated. Sylvia Wilkerson [[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] replace elements of a data frame
Hello, I am getting mad at finding a simple way to replace elements of a df. Here is a short df : names-c(BO,C,CL) price-c(10,25,20) df-data.frame(names,price) I want to replace BO by BOB, C by CR, CL by CLO, and the list is more long. I can do that for each element: df[df==BO]-BOB But my df is bigger indeed with other elements. I was thinking using replace(), but can't get any clean result ( NA or all elements replaced with only one), neither with sapply(). TY for any help, and sorry for the n00b question. Arnaud Gaboury A2CT2 Ltd. __ 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] 3D confidence ellipsoid with ellipse projections onto 2D plane
I have a 3xN matrix of parameters obtained from N regressions where the 3 parameters are jointly statistically significant. I would like to reproduce a 3D confidence ellipsoid projecting 2D ellipses onto the XY plane as in Figure 5.2 in this texthttp://books.google.com/books?id=PMeJGeXA09ECpg=PA172dq=confidence+ellipsoidhl=ensa=Xei=oWQxT4OgBYLSiALW85WtAwved=0CEIQ6AEwAg#v=onepageq=confidence%20ellipsoidf=false . Is this possible using some combination of ellipse3d() and ellipse()? Any insights would be greatly appreciated. Thanks. [[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] replace elements of a data frame
Hi Arnaud, Take a look at require(car) ?recode HTH, Jorge.- On Tue, Feb 7, 2012 at 1:05 PM, Arnaud Gaboury wrote: Hello, I am getting mad at finding a simple way to replace elements of a df. Here is a short df : names-c(BO,C,CL) price-c(10,25,20) df-data.frame(names,price) I want to replace BO by BOB, C by CR, CL by CLO, and the list is more long. I can do that for each element: df[df==BO]-BOB But my df is bigger indeed with other elements. I was thinking using replace(), but can't get any clean result ( NA or all elements replaced with only one), neither with sapply(). TY for any help, and sorry for the n00b question. Arnaud Gaboury A2CT2 Ltd. __ 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] .Random.seed not found
You need to use set.seed() or one of the random number generators to have .Random.seed defined. Also, if you load a previously saved workspace and that workspace had a .Random.seed defined, you will have it defined. HTH Peter On Tue, Feb 7, 2012 at 9:58 AM, Wilkerson, Sylvia (NIH/CIT) [E] syl...@mail.nih.gov wrote: Hi All, I have a user who is receiving this error after running the following: library(rjags) library(R2jags) x = c(1.6907, 1.7242, 1.7552, 1.7842, 1.8113, 1.8369, 1.8610, 1.8839) n = c(59, 60, 62, 56, 63, 59, 62, 60) r = c(6, 13, 18, 28, 52, 53, 61, 60) N = 8 data=list(x,n,r,N) inits=function(){ list(alpha.star=0,beta=0) } parameters=c(alpha,beta,rhat) result=jags(data,inits,model.file=test.bug,parameters,n.chains=3,n.iter=1,n.burnin=1000,n.thin=1) module glm loaded Error: object '.Random.seed not found I don't think it has anything to do with rjags or R2jags because when I run the same code on the same system, I receive no error messages. The R version is 2.14.0 on a Linux 64 machine. Any help would be appreciated. Sylvia Wilkerson [[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] replace elements of a data frame
I did indeed have a look at recode(), and was able to replace, but an error warning : recode(names,BO,BOO,df) Warning message: In recode.default(names, BO, BOO, df) : Name(s) of vars duplicates with an object outside the dataFrame. df names price 1 BOO10 2 C25 3CL20 As you can see, BO has been replaced by BOO, but with a warning! Arnaud Gaboury A2CT2 Ltd. Trade: +41 22 849 88 63 Fax: +41 22 849 88 66 arnaud.gabo...@a2ct2.com This email and any files transmitted with it are confidential and intended solely for the use of the individual or entity to whom they are addressed. Access to this email by anyone else is unauthorized. If you are not the intended recipient, any disclosure, copying, distribution or any action taken or omitted to be taken in reliance on it, is prohibited and may be unlawful. If you have received this email in error please notify the sender. From: Jorge I Velez [mailto:jorgeivanve...@gmail.com] Sent: mardi 7 février 2012 19:55 To: Arnaud Gaboury Cc: r-help@r-project.org Subject: Re: [R] replace elements of a data frame Hi Arnaud, Take a look at require(car) ?recode HTH, Jorge.- On Tue, Feb 7, 2012 at 1:05 PM, Arnaud Gaboury wrote: Hello, I am getting mad at finding a simple way to replace elements of a df. Here is a short df : names-c(BO,C,CL) price-c(10,25,20) df-data.frame(names,price) I want to replace BO by BOB, C by CR, CL by CLO, and the list is more long. I can do that for each element: df[df==BO]-BOB But my df is bigger indeed with other elements. I was thinking using replace(), but can't get any clean result ( NA or all elements replaced with only one), neither with sapply(). TY for any help, and sorry for the n00b question. Arnaud Gaboury A2CT2 Ltd. __ 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] replace elements of a data frame
I used in fact recode() from epilcac package, not the one you mentioned! Arnaud Gaboury A2CT2 Ltd. -Original Message- From: Arnaud Gaboury Sent: mardi 7 février 2012 20:25 To: Jorge I Velez Cc: r-help@r-project.org; Arnaud Gaboury Subject: RE: [R] replace elements of a data frame I did indeed have a look at recode(), and was able to replace, but an error warning : recode(names,BO,BOO,df) Warning message: In recode.default(names, BO, BOO, df) : Name(s) of vars duplicates with an object outside the dataFrame. df names price 1 BOO10 2 C25 3CL20 As you can see, BO has been replaced by BOO, but with a warning! Arnaud Gaboury A2CT2 Ltd. From: Jorge I Velez [mailto:jorgeivanve...@gmail.com] Sent: mardi 7 février 2012 19:55 To: Arnaud Gaboury Cc: r-help@r-project.org Subject: Re: [R] replace elements of a data frame Hi Arnaud, Take a look at require(car) ?recode HTH, Jorge.- On Tue, Feb 7, 2012 at 1:05 PM, Arnaud Gaboury wrote: Hello, I am getting mad at finding a simple way to replace elements of a df. Here is a short df : names-c(BO,C,CL) price-c(10,25,20) df-data.frame(names,price) I want to replace BO by BOB, C by CR, CL by CLO, and the list is more long. I can do that for each element: df[df==BO]-BOB But my df is bigger indeed with other elements. I was thinking using replace(), but can't get any clean result ( NA or all elements replaced with only one), neither with sapply(). TY for any help, and sorry for the n00b question. Arnaud Gaboury A2CT2 Ltd. __ 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] Need to Write a Code that can find the molecular weight of various compounds
Thanks Paul, I was trying to get the function for getMass to input into another column into a data set using R code: FakeCompounds$Molecule - with(FakeCompounds, getMolecule(Molecular.Formula)) showData(FakeCompounds, placement='-20+200', font=getRcmdr('logFont'), maxwidth=80, maxheight=30) FakeCompounds$Molecular.Mass - with(FakeCompounds, getMass(Molecular.Formula)) which gave me the response error: [24] ERROR: $ operator is invalid for atomic vectors Which I didn't quite understand in the context of Rdisop, does it only allow for one variable at a time? So I tried to create a loop, ending up with: c(getFormula(getMolecule(CH3)),getFormula(getMolecule(CH2)))-B i-1 for (i in 1:length(B)){G[i]-getMass(B[i])} but I still get the same error: [25] ERROR: $ operator is invalid for atomic vectors The final goal is to have a list that can calculate molecular mass from a list of chemical formulas. Any help? I'm not sure if this is a problem only related to Rdisop or not. Thanks a lot! Matt -- View this message in context: http://r.789695.n4.nabble.com/Need-to-Write-a-Code-that-can-find-the-molecular-weight-of-various-compounds-tp4342874p4365836.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] lme, lmer, convergence
Hello, all, I am running some simulations to estimate power for a complicated epidemiological study, and am using lme and lmer to get these estimates. I have to run a few thousand iterations, and once in a great while, an iteration will create fake data such that the model won't converge. I see from Google searches that this is not an uncommon situation. My question: is there a way to extract the convergence value from an lme or lmer model? It prints on the screen, but how can I get hold of it to evaluate it with some function? What I'd like to do is build a failsafe into my program so that if a particular model in an iteration doesn't converge, I call a redo on that iteration. This way the program will keep running and not stop in a fit of pique in the middle of my long simulation. If I can't do this, my fallback will be to try setting lmeControl options such that even bad models return parameter estimates etc -- once or twice in 10,000 iterations should not ruin things too badly -- but I'd like to try it the cleaner way first. Erin Jonaitis, Ph.D. Assistant Scientist, Wisconsin Alzheimer's Institute 7818 Big Sky Drive Madison, WI 53719 __ 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] replace elements of a data frame
On 07-02-2012, at 20:24, Arnaud Gaboury wrote: I did indeed have a look at recode(), and was able to replace, but an error warning : recode(names,BO,BOO,df) Warning message: In recode.default(names, BO, BOO, df) : Name(s) of vars duplicates with an object outside the dataFrame. df names price 1 BOO10 2 C25 3CL20 As you can see, BO has been replaced by BOO, but with a warning! library(car) names-c(BO,C,CL) price-c(10,25,20) df-data.frame(names,price) recode(df$names,'BO'='BOO'; 'CL'='CLO'; 'C'='CR') results in [1] BOO CR CLO Levels: BOO CLO CR Note the single quotes. Berend __ 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, lmer, convergence
Erin McMullen Jonaitis jonaitis at wisc.edu writes: Hello, all, I am running some simulations to estimate power for a complicated epidemiological study, and am using lme and lmer to get these estimates. I have to run a few thousand iterations, and once in a great while, an iteration will create fake data such that the model won't converge. I see from Google searches that this is not an uncommon situation. My question: is there a way to extract the convergence value from an lme or lmer model? It prints on the screen, but how can I get hold of it to evaluate it with some function? What I'd like to do is build a failsafe into my program so that if a particular model in an iteration doesn't converge, I call a redo on that iteration. This way the program will keep running and not stop in a fit of pique in the middle of my long simulation. If I can't do this, my fallback will be to try setting lmeControl options such that even bad models return parameter estimates etc -- once or twice in 10,000 iterations should not ruin things too badly -- but I'd like to try it the cleaner way first. There's a somewhat hack-ish solution, which is to use options(warn=2) to 'upgrade' warnings to errors, and then use try() or tryCatch() to catch them. More fancily, I used code that looked something like this to save warnings as I went along (sorry about the - ) in a recent simulation study. You could also check w$message to do different things in the case of different warnings. withCallingHandlers(tryCatch(fun(n=nvec[j],tau=tauvec[i],...), error = function(e) { warn[k,i,j] - paste(ERROR:,e$message) NA_ans}), warning = function(w) { warn[k,i,j] - w$message invokeRestart(muffleWarning) }) In the slightly longer run, we are working on getting the development (lme4Eigen) version of lme4 to save convergence warnings in an accessible slot. I don't know if I would hold my breath for this to be back-ported to nlme, though ... Some of these discussions might be better suited for r-sig-mixed-models at r-project.org Ben Bolker __ 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] survfit is too slow! Looking for an alternative
Hi All I found survfit function was very slow for a large dataset and I am looking for an alternative way to quickly get the predicted survival probabilities. My historical data set is a pool of loans with monthly observed default status for 24 months. I would like to fit the proportional hazard model with time varying covariate such as unemployment rates and time constant variables at loan application in a counting process format, and then use the model to predict the probability of default in each month during next 2 years for a pool of new loans. I have read some posts from other R users. It sounds like using (average survival probability)^exp((X-means(X)*Beta) can quickly get the predicted survival probabilities. My predictors for the model include both continuous variables and categorical variables and my dataset is in counting process format with both time varying and time constant predictors. So how should I take the mean? I guess its the mean of training data? And the denominator for the mean is the number of observations (i.e, the number of rows of training data in the counting process format)? What if the predictor is a categorical variable? Any comments and suggestions are greatly appreciated. Thanks! Ying [[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] Positioning text in top left corner of plot
You might want to look at the grconvertX and grconvertY functions. You can use them to convert a coordinate relative to the plotting device (screen, paper) or plotting region to user coordinates to use with text or mtext. On Tue, Feb 7, 2012 at 4:06 AM, Manta mantin...@libero.it wrote: Dear all, another questions related to zoo plotting. I would like to do as in the subject. Here a reproducible code: library(zoo) par(mfrow=c(2,1) plot(zoo(seq(1:10),as.Date(seq(1:10),origin=1970-01-01)),xlab=,ylab=,main=Value,las=1) mtext(EUR billions,adj=0,cex=0.7) plot(zoo(seq(1:10),as.Date(seq(1:10),origin=1970-01-01)),xlab=,ylab=,main=Value,las=1) par(xpd=T) text(par(usr)[1],par(usr)[3]+10.5,EUR billions,cex=0.7) In the first graph I use the mtext function, which does the trick but it places the text too close too the y-axis. A second possibility is to use the text function, prior specification of the option parameter par(xpd=T). However, this does not look consistent as, depending on the graph, I would need to add different numbers in the y specification of the axis (par(usr)[3]+XXX). Any idea? -- View this message in context: http://r.789695.n4.nabble.com/Positioning-text-in-top-left-corner-of-plot-tp831723p4364355.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. -- Gregory (Greg) L. Snow Ph.D. 538...@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] truncated regression
Hi all, I'm having problems with the command treg (truncated regression) in FEAR 1.15. In the manual, there is a simple explication of the argument that have to be put in the commandbut in the examples there are more passages, like the product between matrix and the command rnorm.trunc: I don't understand why I have to these passages, and above all I don't understand how to apply it with my dataset. Thank you very much to all who is going to help me. Regards Gilda [[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] gmodels error: no method for coercing this S4 class to a vector
Dear All, I'm having a problem using functions in the gmodels library on an object of class mer from the lmer package. Code for a reproducible example is below. # Load lme4 library and sample data library(lme4) library(faraway) library(gmodels) data(penicillin) # Fit a linear mixed effects model fit4 - lmer(yield ~ treat + (1|blend), penicillin) # Extract confidence intervals for fixed effects using gmodels library ci(fit4) The call to 'ci' results in the following error: Error in as.vector(data) : no method for coercing this S4 class to a vector This error does not occur when running the same code on R version 2.13.1. Has anyone else run into this problem? Best regards, Tim Clough --- Session info below --- sessionInfo() R version 2.14.1 (2011-12-22) Platform: x86_64-pc-mingw32/x64 (64-bit) locale: [1] LC_COLLATE=English_United States.1252 [2] LC_CTYPE=English_United States.1252 [3] LC_MONETARY=English_United States.1252 [4] LC_NUMERIC=C [5] LC_TIME=English_United States.1252 attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] gmodels_2.15.1 lme4_0.999375-42 Matrix_1.0-3 lattice_0.20-0 [5] faraway_1.0.5 loaded via a namespace (and not attached): [1] gdata_2.8.2 grid_2.14.1 gtools_2.6.2 MASS_7.3-16 nlme_3.1-103 [6] stats4_2.14.1 tools_2.14.1 __ 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] binomial vs quasibinomial
Not really an R question, now is it ? more like pure stats. I'm guessing you didn't get an answer because this list can't tell you how to analyze your data (or in your case, approve an incorrect analysis). Regarding the part of your question that is R related, I think you may be confused on what the dispersion parameter is. In summary.glm it is reported in the line above the null deviance: (Dispersion parameter for quasibinomial family taken to be ...). You said the dispersion factor is the same for the quasibinomial I doubt that for any real data it gets estimated at 1 (which was assumed for the binomial family model). Also, since you say your p-values are higher on the second model, that leads me to believe it is indeed taken to be greater than 1. Check the two models again. Regards On Tue, Feb 7, 2012 at 4:11 AM, Jhope jeanwaij...@gmail.com wrote: After looking at 48 glm binomial models I decided to try the quasibinomial with the top model 25 (lowest AIC). To try to account for overdispersion (residual deviance 2679.7/68 d.f.) After doing so the dispersion factor is the same for the quasibinomial and less sectors of the beach were significant by p-value. While the p-values in the binomial were more significant for each section of the beach. -- telling me more about the beach. Is this ok? Can I just look at the binomial glm model 25 and look at its p-values for beach sections and forget about the quasibinomial model 25? J Call: glm(formula = cbind(Shells, TotalEggs - Shells) ~ Sector:Veg:Aeventexhumed, family = quasibinomial, data = data.to.analyze) -- View this message in context: http://r.789695.n4.nabble.com/binomial-vs-quasibinomial-tp4364371p4364371.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@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] fixed effects linear model in R
Dear R-helpers, First of all, sorry for those who have (eventually) already received that request. The mail has been bumped several times, so I am not sure the list has received it... and I need help (if you have time)! ;-) I have a very simple question and I really hope that someone could help me I would like to estimate a simple fixed effect regression model with clustered standard errors by individuals. For those using Stata, the counterpart would be xtreg with the fe option, or areg with the absorb option and in both case the clustering is achieved with vce(cluster id) My question is : how could I do that with R ? An important point is that I have too many individuals, therefore I cannot include dummies and should use the demeaning usual procedure. I tried with the plm package with the within option, but R quikcly tells me that the memory limits are attained (I have over 10go ram!) while the dataset is only 700mo (about 50 000 individuals, highly unbalanced) I dont understand... plm do indeed demean the data so the computation should be fast and light enough... and instead it saturates my memory and do not converge... Do you have an idea ? Moreover, it is possible to obtain cluster robust standard errors with plm ? Are there any other solutions for fixed effects linear models (with the demean trick in order not to create as many dummies as individuals) ? Many thanks in advance ! ;) John __ 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 need
I have mad a for loop to try and output values which i have named spectrum. However, I cannot seem to get the answers to come out as a vector which is what i need. They come out as separate values which I am then unable to join together. Thank you for(f in seq(0,0.5,0.1)) { sigmasqaured - 1 i = complex(real = 0, imaginary = 1) spectrum - (sigmasqaured)/(abs(1-2.7607*exp(2*pi*i*f)+3.8106*exp(4*pi*i*f)-2.6535*exp(6*pi*i*f)+0.9258*exp(8*pi*i*f))^2) print(spectrum) } [[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] fixed effects linear model in R
Based on Paul Allison's booklet Fixed Effect Regression Models (2009), the FE model can be estimated by person-mean centering all of your variables (but not the outcome), and then including a random intercept for each person. The centering gives you the FE model estimates, and the random intercept adjusts the standard errors for clustering by individuals. Note that your data must be in person-period (or long) format to do this. In case you are unfamiliar with person-mean centering, that simply means taking the mean of each person's values for a given variable for all of the periods in your data, and then calculating a deviation from that mean at each time period. For example, a person's average income over four years might be $50,000, but in each year their actual income would be slightly higher or lower than this (these would be the person-mean deviations). In symbolic form, your code might look something like this: library(lme4) variable_pmcentered = variable - person_mean mod = lmer(outcome ~ variable_pmcentered + person_mean + other predictors + (1|personID)) The advantage of this method (which Allison calls a hybrid method) over traditional FE models is that you get the benefits of a FE model (subtracting out time-invariant omitted variables) along with the benefits of random effect models (e.g., estimating coefficients for time-invariant variables, estimating interactions with time, letting intercepts and slopes varying randomly, etc.) See Allison's booklet for more details on this method. Allison, Paul D. 2009. Fixed Effects Regression Models. Los Angeles, C.A.: Sage. Andrew Miles On Feb 7, 2012, at 5:00 PM, caribou...@gmx.fr wrote: Dear R-helpers, First of all, sorry for those who have (eventually) already received that request. The mail has been bumped several times, so I am not sure the list has received it... and I need help (if you have time)! ;-) I have a very simple question and I really hope that someone could help me I would like to estimate a simple fixed effect regression model with clustered standard errors by individuals. For those using Stata, the counterpart would be xtreg with the fe option, or areg with the absorb option and in both case the clustering is achieved with vce(cluster id) My question is : how could I do that with R ? An important point is that I have too many individuals, therefore I cannot include dummies and should use the demeaning usual procedure. I tried with the plm package with the within option, but R quikcly tells me that the memory limits are attained (I have over 10go ram!) while the dataset is only 700mo (about 50 000 individuals, highly unbalanced) I dont understand... plm do indeed demean the data so the computation should be fast and light enough... and instead it saturates my memory and do not converge... Do you have an idea ? Moreover, it is possible to obtain cluster robust standard errors with plm ? Are there any other solutions for fixed effects linear models (with the demean trick in order not to create as many dummies as individuals) ? Many thanks in advance ! ;) John __ 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] Help need
Instead of a for loop, why not use the vectorization inherent in R? sigmasqaured - 1 i - complex(real = 0, imaginary =1) f - seq(0,0.5,0.1) spectrum - (sigmasqaured)/(abs(1-2.7607*exp(2*pi*i*f)+3.8106*exp(4*pi*i*f)-2.6535*exp(6*pi*i*f)+0.9258*exp(8*pi*i*f))^2) spectrum [1] 9.632720e+00 1.411130e+03 2.947753e+00 6.479994e-02 1.295175e-02 8.042731e-03 On Tue, Feb 7, 2012 at 1:08 PM, Jaymin Shah jayminsh...@live.com wrote: I have mad a for loop to try and output values which i have named spectrum. However, I cannot seem to get the answers to come out as a vector which is what i need. They come out as separate values which I am then unable to join together. Thank you for(f in seq(0,0.5,0.1)) { sigmasqaured - 1 i = complex(real = 0, imaginary = 1) spectrum - (sigmasqaured)/(abs(1-2.7607*exp(2*pi*i*f)+3.8106*exp(4*pi*i*f)-2.6535*exp(6*pi*i*f)+0.9258*exp(8*pi*i*f))^2) print(spectrum) } [[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. [[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] Vectorizing a loop
On Tue, 7 Feb 2012, David Winsemius wrote: On Feb 7, 2012, at 12:56 PM, Jeff Newmiller wrote: On Tue, 7 Feb 2012, Alexander Shenkin wrote: Hello Folks, I'm trying to vectorize a loop that processes rows of a dataframe. It involves lots of conditionals, such as If column 10 == 3, and if column 3 is True, and both column 5 and 6 are False, then set column 4 to True. So, for example, any ideas about vectorizing the following? df = data.frame( list(a=c(1,2,3,4), b=c(a,b,c,d), c=c(T,F,T,F), d=NA, e=c(F,F,T,T)) ) for (i in 1:nrow(df)) { if (df[i,3] %in% c(FALSE,NA) (df[i,1] 2 | df[i,5]) ) { df[i,4] = 1 } if (df[i,5] %in% c(TRUE, NA) df[i,2] == b) { df[i,4] = 2 df[i,5] = T } } Your code attempts to do some things with NA that won't behave the way you expect them to. Specifically, you cannot use %in% to test for NA, Huh? NA %in% NA [1] TRUE NA %in% c(5, NA) [1] TRUE NA %in% c(5, 6) [1] FALSE Sorry, SQL rules bleeding through... %in% is clearly more forgiving in R than IN is in SQL. However, the second if did check whether df[i,5] was NA, yet the first if did not. Since comparisons with NA are neither false nor true that test failed. NA | 1 [1] TRUE NA 1 [1] NA NA 1 [1] NA # intermediate logical vectors for clarity tmp - ( is.na(df[[3]]) | !df[[3]] ) ( df[[1]] 2 | df[[5]] ) tmp2 - ( is.na(df[[5]]) | df[[5]] ) df[[2]] == b df[ tmp, d ] - 1 df[ tmp2, d ] - 2 df[ tmp2, e ] - TRUE --- Jeff NewmillerThe . . Go Live... DCN:jdnew...@dcn.davis.ca.usBasics: ##.#. ##.#. Live Go... Live: OO#.. Dead: OO#.. Playing Research Engineer (Solar/BatteriesO.O#. #.O#. with /Software/Embedded Controllers) .OO#. .OO#. rocks...1k __ 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(0.9.0): could not find function ==
I'm curious if you have a guess whether the issue I was having is a result of the problems with the 0.9.0 version or if they're due to fundamental changes in ggplot2? It looks like a bug - I'll add it to the to do list. Hadley -- Assistant Professor / Dobelman Family Junior Chair Department of Statistics / Rice University 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] Help need
You need to store the values of each iteration in a vector, and then display the vector after you the loop terminates. I made a few updates to your code, and it seems to do what you want now. And thanks for including the code. That made is easy to know how to help. spectrum = c() for(f in seq(0,0.5,0.1)) { sigmasqaured - 1 i = complex(real = 0, imaginary = 1) spectrum - c(spectrum, (sigmasqaured)/(abs(1-2.7607*exp(2*pi*i*f)+3.8106*exp(4*pi*i*f)-2.6535*exp(6*pi*i*f)+0.9258*exp(8*pi*i*f))^2)) } spectrum Andrew Miles On Feb 7, 2012, at 4:08 PM, Jaymin Shah wrote: I have mad a for loop to try and output values which i have named spectrum. However, I cannot seem to get the answers to come out as a vector which is what i need. They come out as separate values which I am then unable to join together. Thank you for(f in seq(0,0.5,0.1)) { sigmasqaured - 1 i = complex(real = 0, imaginary = 1) spectrum - (sigmasqaured)/(abs(1-2.7607*exp(2*pi*i*f)+3.8106*exp(4*pi*i*f)-2.6535*exp(6*pi*i*f)+0.9258*exp(8*pi*i*f))^2) print(spectrum) } [[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. [[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] rJava load failure, 64-bit R, 64-bit Java, R 2.14.1 (works fine on R 2.12.2 64-bit, same computer)
Hello R users, I have encountered difficulty in attempting to load the package rJava for 64-bit R 2.14.1. It appears that the problem is specific (for my system) to this version of R. I had no trouble loading rJava on an earlier version (2.12.2) of R. N:\java -version java version 1.6.0_25 Java(TM) SE Runtime Environment (build 1.6.0_25-b06) Java HotSpot(TM) 64-Bit Server VM (build 20.0-b11, mixed mode) I used the following to install rJava: install.packages('rJava', .libPaths()[1], 'http://www.rforge.net/') library(rJava) Error : .onLoad failed in loadNamespace() for 'rJava', details: call: inDL(x, as.logical(local), as.logical(now), ...) error: unable to load shared object 'C:/Program Files/R/R-2.14.1/library/rJava/libs/x64/rJava.dll': LoadLibrary failure: The specified path is invalid. Error: package/namespace load failed for 'rJava' The specified file is in fact present, and my R console says (64-bit) at the top I am aware of the following discussion on this topic: http://stackoverflow.com/questions/2399027/cannot-load-rjava-because-cannot-load-a-shared-library I have updated my path to include the following: C:\Program Files\R\R-2.14.1\bin\x64 And C:\Program Files\Java\jre6\bin\ When I download the package and install from a local zip I get the following error: library(rJava) Error in get(Info[i, 1], envir = env) : cannot allocate memory block of size 2.9 Gb Error: package/namespace load failed for 'rJava' There should not be any memory issues: memory.limit() [1] 16381 Thank you! Jacob L Strunk [[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] How to figure out the number of data values in a list?
Hello, I am working with the following data: nb10 - read.table(http://www.adjoint-functors.net/su/web/314/R/NB10;) I need to figure out how many data values are in this list. -- View this message in context: http://r.789695.n4.nabble.com/How-to-figure-out-the-number-of-data-values-in-a-list-tp4366679p4366679.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] Setting up infile for R CMD BATCH
Suppose I create an R program called myTest.R with only one line like the following: type - as.integer(readline(input type (1: type1; 2: type2)? )) Then I'd like to run myTest.R in batch mode by constructing an input file called answers.R with the following: source(myTest.R) 1 When I ran the following at the terminal: R CMD BATCH answer.R output.Rout it failed to pick up the answer '1' from the 2nd line in answers.R as shown inside output.Rout: source(myTest.R) input type (0: quit; 1: type1; 2: type2)? 1 [1] 1 What am I missing here? Thanks in advance, Gang __ 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 figure out the number of data values in a list?
-Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of Ajata Paul Sent: Tuesday, February 07, 2012 2:49 PM To: r-help@r-project.org Subject: [R] How to figure out the number of data values in a list? Hello, I am working with the following data: nb10 - read.table(http://www.adjoint-functors.net/su/web/314/R/NB10;) I need to figure out how many data values are in this list. Maybe, nrow(nb10) ? You might want to read the introductory manual that comes with every R installation. Dan Daniel J. Nordlund Washington State Department of Social and Health Services Planning, Performance, and Accountability Research and Data Analysis Division Olympia, WA 98504-5204 __ 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] Setting up infile for R CMD BATCH
You're not missing anything. In your output.Rout: the 1 right after the source('test') is the 1 inputed from answers.R. the [1] 1 is the result of test. Remove the second line from answers.R and see what happens (hint: script ends after the readline prompt). Just out of curiosity, why will you use a script that requires user input (readlines) in batch mode ? Cheers On Tue, Feb 7, 2012 at 4:05 PM, Gang Chen gangch...@gmail.com wrote: Suppose I create an R program called myTest.R with only one line like the following: type - as.integer(readline(input type (1: type1; 2: type2)? )) Then I'd like to run myTest.R in batch mode by constructing an input file called answers.R with the following: source(myTest.R) 1 When I ran the following at the terminal: R CMD BATCH answer.R output.Rout it failed to pick up the answer '1' from the 2nd line in answers.R as shown inside output.Rout: source(myTest.R) input type (0: quit; 1: type1; 2: type2)? 1 [1] 1 What am I missing here? Thanks in advance, Gang __ 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] units for mapproject() function result
Does anyone know what the units are for projected coordinates obtained using mapproj's mapproject function with an Albers projection? Thanks for any and all help! Buck Stockhausen *** * Dr. William T. Stockhausen * *** * Resource Ecology and Fisheries Management * * Alaska Fisheries Science Center * * National Marine Fisheries Service * * National Oceanic and Atmospheric Administration * * 7600 Sand Point Way N.E.* * Seattle, Washington 98115-6349 * *** * email: william.stockhau...@noaa.gov * * voice: 206-526-4241 fax: 206-526-6723 * * web : http://www.afsc.noaa.gov * *** All models are wrong, some are useful.--G.E.P. Box Beware of geeks bearing equations.--W. Buffett *** Disclaimer: The opinions expressed above are personal and do not necessarily reflect official NOAA policy. [[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] R equivalent of Python str()?
Hi, I was wondering if there's a function in R that is meant to return a string representation of an object. Basically, it's like print() but it doesn't print anything, it only returns a string. I know there's a str() function but it's not quite the same. I mean a function that returns the same string that print() would display. -- Bye, Ernest __ 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 equivalent of Python str()?
?dump ?dput 2012/2/7 Ernest Adrogué nfdi...@gmail.com: Hi, I was wondering if there's a function in R that is meant to return a string representation of an object. Basically, it's like print() but it doesn't print anything, it only returns a string. I know there's a str() function but it's not quite the same. I mean a function that returns the same string that print() would display. -- Bye, Ernest __ 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. -- Jim Holtman Data Munger Guru What is the problem that you are trying to solve? Tell me what you want to do, not how you want to do it. __ 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] Zoomable time series plots
Not sure if the question is appropos, but I have multiple csv's which are read into an xts object, corresponding to telemetry data (accelerometer, magnetometer/compass, and gyroscope). For examination, it would be über useful if plot.zoo (or something similar) allowed me to zoom in and out of a subset of the time axes. -- Sent from my mobile device Envoyait de mon portable __ 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] Zoomable time series plots
On Tue, Feb 7, 2012 at 8:02 PM, Hasan Diwan hasan.di...@gmail.com wrote: Not sure if the question is appropos, but I have multiple csv's which are read into an xts object, corresponding to telemetry data (accelerometer, magnetometer/compass, and gyroscope). For examination, it would be über useful if plot.zoo (or something similar) allowed me to zoom in and out of a subset of the time axes. In the examples section of ?xyplot.zoo (in the zoo package) see the first example in the ## Not run: portion using playwith. -- Statistics Software Consulting GKX Group, GKX Associates Inc. tel: 1-877-GKX-GROUP email: ggrothendieck at 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 equivalent of Python str()?
Possibly as.character() is what the OP was seeking Michael On Feb 7, 2012, at 7:15 PM, jim holtman jholt...@gmail.com wrote: ?dump ?dput 2012/2/7 Ernest Adrogué nfdi...@gmail.com: Hi, I was wondering if there's a function in R that is meant to return a string representation of an object. Basically, it's like print() but it doesn't print anything, it only returns a string. I know there's a str() function but it's not quite the same. I mean a function that returns the same string that print() would display. -- Bye, Ernest __ 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. -- Jim Holtman Data Munger Guru What is the problem that you are trying to solve? Tell me what you want to do, not how you want to do it. __ 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] How to figure out the number of data values in a list?
Also, this is the 3rd post in the last ~36 hours asking elementary questions about that data set and my homework senses are tingling. Perhaps a moratorium on them? M On Feb 7, 2012, at 6:09 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.gov wrote: -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of Ajata Paul Sent: Tuesday, February 07, 2012 2:49 PM To: r-help@r-project.org Subject: [R] How to figure out the number of data values in a list? Hello, I am working with the following data: nb10 - read.table(http://www.adjoint-functors.net/su/web/314/R/NB10;) I need to figure out how many data values are in this list. Maybe, nrow(nb10) ? You might want to read the introductory manual that comes with every R installation. Dan Daniel J. Nordlund Washington State Department of Social and Health Services Planning, Performance, and Accountability Research and Data Analysis Division Olympia, WA 98504-5204 __ 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] units for mapproject() function result
On Wed, 08 Feb 2012, William Stockhausen wrote: Does anyone know what the units are for projected coordinates obtained using mapproj's mapproject function with an Albers projection? Thanks for any and all help! Buck Stockhausen I don't know for sure, but it looks like radians to me, with some unspecified origin(depending on the parameters specified). Certainly the maps package data is specified in radians internally. Hope this helps, Ray Brownrigg __ 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] Table rearranging
Hi David, I am not sure how ddply/summarize solves my issue. I have the following table: ID measurement date door color 1 0.93529385 513 open red 2 0.97419293 420 open red 3 0.962053514 513 closed red 4 0.963909937 1230 open blue 5 0.97652034 1230 open green 6 0.989310795 1230 closed blue 7 0.9941022 917 closed yellow 8 0.8945757 1230 open blue I only want to keep the lines that have corresponding open/closed measurements. For example, I want to keep lines 4,6,8 because for the 1230 blue condition, there exists both open and closed measurements. However, the 513 red condition has an open measurement, but no closed measurement. Therefore, line 1 should be deleted. Jeffrey CC: r-help@r-project.org; wdun...@tibco.com From: dwinsem...@comcast.net To: johjeff...@hotmail.com Subject: Re: [R] Table rearranging Date: Tue, 7 Feb 2012 09:08:00 -0500 On Feb 7, 2012, at 4:21 AM, Jeffrey Joh wrote: Thank you for your help, Bill. From the original table (not the plyr output), I would like to remove all the lines that do not have a corresponding open/closed measurement. For example, if there is a Closed yellow measurement on 0917, but not an Open yellow 0917 measurement, then the Closed yellow should be deleted. How can I make this change? In R you need to assign the results of a function to an object name so you code would look like: modified_data - ddply(d, .(date, color), summarize, meanClosed=mean(measurement[door==closed]), nClosed=sum(door==closed)) -- David Jeffrey From: wdun...@tibco.com To: johjeff...@hotmail.com; r-help@r-project.org Subject: RE: [R] Table rearranging Date: Tue, 7 Feb 2012 00:43:25 + Install and load the plyr package and try something like: ddply(d, .(date, color), summarize, + ddply(d, .(date, color), summarize + meanClosed=mean(measurement[door==closed]), nClosed=sum(door==closed)) date color meanOpen nOpen meanClosed nClosed 1 420 red 0.9741929 1 NaN 0 2 513 red 0.9352938 1 0.9620535 1 3 917 yellow NaN 0 0.9941022 1 4 1230 blue 0.9639099 1 0.9893108 1 5 1230 green 0.9765203 1 NaN 0 Bill Dunlap Spotfire, TIBCO Software wdunlap tibco.com -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org ] On Behalf Of Jeffrey Joh Sent: Monday, February 06, 2012 4:28 PM To: r-help@r-project.org Subject: [R] Table rearranging I have a table that looks like this: measurement date door color 0.93529385 513 open red 0.97419293 420 open red 0.962053514 513 closed red 0.963909937 1230 open blue 0.97652034 1230 open green 0.989310795 1230 closed blue 0.9941022 917 closed yellow I would like to create a table that has: Open measurement, Closed measurement, date, color. For every date/color combination, there should be two columns to represent the door open/closed measurement. If there are multiple datapoints with a given door/date/color combination, then they should be averaged. I would also like to make two columns to represent the number of datapoints that were averaged in determining the open/closed measurements. Jeffrey __ 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 Winsemius, MD West Hartford, CT __ 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] Setting up infile for R CMD BATCH
Thanks for the help. You're not missing anything. In your output.Rout: the 1 right after the source('test') is the 1 inputed from answers.R. the [1] 1 is the result of test. Remove the second line from answers.R and see what happens (hint: script ends after the readline prompt). That number '1' at the 2nd line was meant to be the answer for the readline() part, but apparently it does not worked as I intended. If the script ends right after the readline prompt as you suggested, my question is, how can I feed in this answer '1' in the infile 'answers.R'? Just out of curiosity, why will you use a script that requires user input (readlines) in batch mode ? I know this sounds contradictory, but my intention is that, in case the user has all the answers available from a previous run in interactive mode, s/he may try out the batch mode the next time. On Tue, Feb 7, 2012 at 4:05 PM, Gang Chen gangch...@gmail.com wrote: Suppose I create an R program called myTest.R with only one line like the following: type - as.integer(readline(input type (1: type1; 2: type2)? )) Then I'd like to run myTest.R in batch mode by constructing an input file called answers.R with the following: source(myTest.R) 1 When I ran the following at the terminal: R CMD BATCH answer.R output.Rout it failed to pick up the answer '1' from the 2nd line in answers.R as shown inside output.Rout: source(myTest.R) input type (0: quit; 1: type1; 2: type2)? 1 [1] 1 What am I missing here? Thanks in advance, Gang __ 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.