Re: [R] zero line
Hi -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of bgnumis Sent: Sunday, March 17, 2013 4:12 PM To: r-help@r-project.org Subject: [R] zero line Hi all, I´m plotting cf (with two axis) and addind a shaded color up and down on the 0 line x axis (tfr1 is the time). The thing is that when I plot this graph adds a line up on the first plot. I hope you can understand what I mean. Not much plot(cf[,1],type=l,col=black, xlab=Time, ylab=Correlation, axes=F) Error in plot(cf[, 1], type = l, col = black, xlab = Time, ylab = Correlation, : object 'cf' not found No object cf here. It would be good if you make your code reproducible by putting a result of dput(cf) and make available other objects used for plotting. Or you can try to present some fake data illustrating your points. How should I erase this sencond line, it is suposed they have the same data? Why zero line doesn´t complete the shaded area? plot(cf[,1],type=l,col=black, xlab=Time, ylab=Correlation, axes=F) grid() cord.x - c(0,1:length(cf[,1]), length(cf[,1])) cord.y - c(0,cf[,1],0) polygon(cord.x,cord.y,col='skyblue') par(new=T) From below code you have some object named T. Hence your above code probably do not do what you expect as your T is put in par not reserved word TRUE. Regards Petr variable-matrix(0,length(cf[,1]),1) plot(tfr1[(T+Tcor):(nrow(tfr1)),1],variable,type=l,ylab=,col=black ,xlab=) mtext(Correlation,side=4,line=2,col=4) title(Dinamic Correlation,font=4) [[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 plot u-v wind by R?
hi R users: I have a dataset including u wind in x-axis and v wind in y-axis. How can I plot the u,v wind data in vector or barb figure? which command ? thank you . -- TANG Jie [[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] how to plot u-v wind by R?
On 03/18/2013 08:50 PM, Jie Tang wrote: hi R users: I have a dataset including u wind in x-axis and v wind in y-axis. How can I plot the u,v wind data in vector or barb figure? which command ? thank you . Hi Jie, Maybe the vectorField function (plotrix)? Jim __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting CIE chromaticity diagram?
As it turns out, I am working on this right now, and it's a bit of a challenge to get the color gradient correct. I don't think this is really the right kind of question for this list, so let's see if anyone has attempted it or has ideas, and then take the conversation private. I can send you my not quite perfect attempt a little later today. Bryan Prof. Bryan Hanson Dept of Chemistry Biochemistry DePauw University Greencastle IN 46135 USA academic.depauw.edu/~hanson/deadpezsociety.html github.com/bryanhanson academic.depauw.edu/~hanson/UMP/Index.html On Mar 18, 2013, at 6:12 AM, ishi soichi soichi...@gmail.com wrote: Has anyone plotted or is it possible to plot CIE *xy* chromaticity diagram http://en.wikipedia.org/wiki/File:CIE1931xy_blank.svg I need this plot in color. ishida [[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] plotting CIE chromaticity diagram?
ishi soichi soichi777 at gmail.com writes: Has anyone plotted or is it possible to plot CIE *xy* chromaticity diagram http://en.wikipedia.org/wiki/File:CIE1931xy_blank.svg I need this plot in color. ishida I think that plotting the spectral locus and the line of purples is trivial, but that's not actually what you are asking. Though filling in the interior with a color gradient would be a challenging and interesting problem, there are those who would consider this a misguided attempt at representation (even though it is frequently done) as the CIE coordinates provide no information on color appearance, per se. They specify what lights look alike, not what they look like. -- Kenneth Knoblauch Inserm U846 Stem-cell and Brain Research Institute Department of Integrative Neurosciences 18 avenue du Doyen Lépine 69500 Bron France tel: +33 (0)4 72 91 34 77 fax: +33 (0)4 72 91 34 61 portable: +33 (0)6 84 10 64 10 http://www.sbri.fr/members/kenneth-knoblauch.html __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting CIE chromaticity diagram?
ishi soichi soichi777 at gmail.com writes: Has anyone plotted or is it possible to plot CIE *xy* chromaticity diagram http://en.wikipedia.org/wiki/File:CIE1931xy_blank.svg I need this plot in color. ishida sermon And following up on my previous mail (diatribe), after having contemplated the image at the link that you provide, there is another significant distortion to take into account in the representation that you propose. Suppose that you are considering the appearance of isolated, centrally viewed lights in isolation, in neutral adaptation by an observer whose vision corresponds to the average observer of CIE1931. Then, the set of primaries that you will be using to generate those colors (on a screen or on paper) will have a gamut that excludes a significant portion of the diagram, so that if you take this limited gamut and stretch it out to fill the diagram, then the coordinates of the colors indicated will deviate significantly from the positions in which they are represented (not to even speak of their appearance). So, in general, while this leads to a pretty diagram, there is a lot of potential for misunderstanding in such a (mis)representation. /sermon -- Kenneth Knoblauch Inserm U846 Stem-cell and Brain Research Institute Department of Integrative Neurosciences 18 avenue du Doyen Lépine 69500 Bron France tel: +33 (0)4 72 91 34 77 fax: +33 (0)4 72 91 34 61 portable: +33 (0)6 84 10 64 10 http://www.sbri.fr/members/kenneth-knoblauch.html __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting CIE chromaticity diagram?
I am, unfortunately, well-aware of the limitations that Ken points out (and I do appreciate him making these points). One can readily demonstrate the gamut limitations by printing the diagram Ishida links to on different devices. My hope is to get something close and include a disclaimer. Bryan On Mar 18, 2013, at 7:08 AM, Ken Knoblauch ken.knobla...@inserm.fr wrote: ishi soichi soichi777 at gmail.com writes: Has anyone plotted or is it possible to plot CIE *xy* chromaticity diagram http://en.wikipedia.org/wiki/File:CIE1931xy_blank.svg I need this plot in color. ishida sermon And following up on my previous mail (diatribe), after having contemplated the image at the link that you provide, there is another significant distortion to take into account in the representation that you propose. Suppose that you are considering the appearance of isolated, centrally viewed lights in isolation, in neutral adaptation by an observer whose vision corresponds to the average observer of CIE1931. Then, the set of primaries that you will be using to generate those colors (on a screen or on paper) will have a gamut that excludes a significant portion of the diagram, so that if you take this limited gamut and stretch it out to fill the diagram, then the coordinates of the colors indicated will deviate significantly from the positions in which they are represented (not to even speak of their appearance). So, in general, while this leads to a pretty diagram, there is a lot of potential for misunderstanding in such a (mis)representation. /sermon -- Kenneth Knoblauch Inserm U846 Stem-cell and Brain Research Institute Department of Integrative Neurosciences 18 avenue du Doyen Lépine 69500 Bron France tel: +33 (0)4 72 91 34 77 fax: +33 (0)4 72 91 34 61 portable: +33 (0)6 84 10 64 10 http://www.sbri.fr/members/kenneth-knoblauch.html __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting CIE chromaticity diagram?
Unfortunately, a bad idea that looks good never dies. That said, I still think that it is an interesting challenge to create this graph by programming. I can imagine that if you specified the coordinates at a certain number of points in the interior and along the spectrum locus, you could use an interpolation algorithm, to reproduce what is on the wikipedia page. Of course, I'm kicking myself for encouraging this be done given the opportunity it provides for misguidance and misrepresentation. Ken Quoting Bryan Hanson han...@depauw.edu: I am, unfortunately, well-aware of the limitations that Ken points out (and I do appreciate him making these points). One can readily demonstrate the gamut limitations by printing the diagram Ishida links to on different devices. My hope is to get something close and include a disclaimer. Bryan On Mar 18, 2013, at 7:08 AM, Ken Knoblauch ken.knobla...@inserm.fr wrote: ishi soichi soichi777 at gmail.com writes: Has anyone plotted or is it possible to plot CIE *xy* chromaticity diagram http://en.wikipedia.org/wiki/File:CIE1931xy_blank.svg I need this plot in color. ishida sermon And following up on my previous mail (diatribe), after having contemplated the image at the link that you provide, there is another significant distortion to take into account in the representation that you propose. Suppose that you are considering the appearance of isolated, centrally viewed lights in isolation, in neutral adaptation by an observer whose vision corresponds to the average observer of CIE1931. Then, the set of primaries that you will be using to generate those colors (on a screen or on paper) will have a gamut that excludes a significant portion of the diagram, so that if you take this limited gamut and stretch it out to fill the diagram, then the coordinates of the colors indicated will deviate significantly from the positions in which they are represented (not to even speak of their appearance). So, in general, while this leads to a pretty diagram, there is a lot of potential for misunderstanding in such a (mis)representation. /sermon -- Kenneth Knoblauch Inserm U846 Stem-cell and Brain Research Institute Department of Integrative Neurosciences 18 avenue du Doyen Lépine 69500 Bron France tel: +33 (0)4 72 91 34 77 fax: +33 (0)4 72 91 34 61 portable: +33 (0)6 84 10 64 10 http://www.sbri.fr/members/kenneth-knoblauch.html __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Kenneth Knoblauch Inserm U846 Stem-cell and Brain Research Institute Department of Integrative Neurosciences 18 avenue du Doyen Lépine 69500 Bron France tel: +33 (0)4 72 91 34 77 fax: +33 (0)4 72 91 34 61 portable: +33 (0)6 84 10 64 10 http://www.sbri.fr/members/kenneth-knoblauch.html __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Q-Q plot for chi-quadrat w/o knowing the contingency table
At 14:56 17/03/2013, Wim Kreinen wrote: Hello, I have chi2-values from 2x2 contingency tables, each degree of freedom =1. But I dont know the values of the contingency table. All I have in addtion: The allelic odds ratio. Now I would like to do a Q-Q-plot. Is that possible? If I understand you correctly ?qqplot should help you here. If yes, please let me know. If no, please let me know as well. Thanks Wim [[alternative HTML version deleted]] Michael Dewey i...@aghmed.fsnet.co.uk http://www.aghmed.fsnet.co.uk/home.html __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] How many samples ACTUALLY used in regression?
Dear All, is there a simple way that covers all regression models to extract the number of samples from a data frame/matrix actually used in a regression model? For instance I might have a data of 100 rows and 4 colums (1 response + 3 explanatory variables). If 3 samples have one or more NAs in the explanatory variable columns these samples will be dropped in any model: my.model = lm(y ~ x + w + z, my.data) my.model = glm(y ~ x + w + z, my.data, family = binomial) my.model = polr(y ~ x + w + z, my.data) … I don't seem to be able to find one single method that works in the exact same way -- irrespective of the model type -- to interrogate my.model to see how many samples of my.data were actually used. Is there such function or do I need to hack something together? Best wishes Federico __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Counting confidence intervals
Hi, I have a 2 x 1 matrix of confidence intervals. The first column is the lower and the next column is the upper. I want to cont how many times a number say 12 lies in the interval. Can anyone assist? -- Thanks, Jim. [[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] !0 + !0 == !0 - !0
Bert Gunter gunter.berton at gene.com writes: But this has nothing to do with 7.31 and everything to do with operator precedence and automatic casting from integers to logical and vice-versa. I also think it fair to say that all (??) languages have these sorts of malapropisms due to operator precedence. Maybe FAQ 7.31 was referred to not for its direct relevance but as a measure of the old-hand-ness of the people who will get the joke. I have to admit it took me a minute ... 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] How many samples ACTUALLY used in regression?
Federico Calboli f.calboli at imperial.ac.uk writes: is there a simple way that covers all regression models to extract the number of samples from a data frame/matrix actually used in a regression model? For instance I might have a data of 100 rows and 4 colums (1 response + 3 explanatory variables). If 3 samples have one or more NAs in the explanatory variable columns these samples will be dropped in any model: my.model = lm(y ~ x + w + z, my.data) my.model = glm(y ~ x + w + z, my.data, family = binomial) my.model = polr(y ~ x + w + z, my.data) I don't seem to be able to find one single method that works in the exact same way -- irrespective of the model type -- to interrogate my.model to see how many samples of my.data were actually used. Is there such function or do I need to hack something together? I haven't tested it (don't want to bother to put together the test data), but does nrow(model.frame(my.model)) work ? __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Counting confidence intervals
Hi Jim, Try either of the following (untested): sum( x[1, ] 12 x[2, ] 12) sum(apply(x, 2, function(x) x[1] 12 x[2] 12)) where x is your 2x1000 matrix. HTH, Jorge.- On Tue, Mar 19, 2013 at 12:03 AM, Jim Silverton wrote: Hi, I have a 2 x 1 matrix of confidence intervals. The first column is the lower and the next column is the upper. I want to cont how many times a number say 12 lies in the interval. Can anyone assist? -- Thanks, Jim. [[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] Counting confidence intervals
sum(M[1]12 12=M[2]) untested, no data --- 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 --- Sent from my phone. Please excuse my brevity. Jim Silverton jim.silver...@gmail.com wrote: Hi, I have a 2 x 1 matrix of confidence intervals. The first column is the lower and the next column is the upper. I want to cont how many times a number say 12 lies in the interval. Can anyone assist? __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] header intact
-- Thanks Regards, Kamal Kishore Research Scholar Department of Biostatistics National Institute of Mental Health and Neuro Sciences (NIMHANS) Banglore-560029 Mob-+91-9591349768 Email:kamalkishorest...@gmail.com [[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] try/tryCatch
Hi All, I have tried every fix on my try or tryCatch that I have found on the internet, but so far have not been able to get my R code to continue with the for loop after the lmer model results in an error. Here is two attemps of my code, the input is a 3D array file, but really any function would do metatrialstry-function(mydata){ a-matrix(data=NA, nrow=dim(mydata)[3], ncol=5) #colnames(a)=c(sens, spec, corr, sens_se, spec_se, counter)#colnames(a)=c(sens, spec, corr, sens_se, spec_se) k=1 for(ii in 1:dim(mydata)[3]){ tmp-mydata[,,ii] tmp1-as.data.frame(tmp) names(tmp1)=c(persons, d1, tp, fn, fp, fn, detect, d0, outcome) lm1-try(lmer(outcome~0+d1+d0+(0+d1+d0 | persons), family=binomial, data=tmp1, nAGQ=3), silent=T) if(class(lm1)[1]!='try-error'){ a[ii,1]=lm1@fixef[1] a[ii,2]=lm1@fixef[2] a[ii,3]=vcov(lm1)[1,2]/prod(sqrt(diag(vcov(lm1 a[ii,4:5]=sqrt(diag(vcov(lm1))) } } #k=k+1 #a[ii,6]=k return(a) } # # try / try catch ### # metatrialstry-function(mydata){ a-matrix(data=NA, nrow=dim(mydata)[3], ncol=5) #colnames(a)=c(sens, spec, corr, sens_se, spec_se, counter) colnames(a)=c(sens, spec, corr, sens_se, spec_se) #a[,6]=rep(0, length(a[,6])) for(ii in 1:dim(mydata)[3]){ tmp-mydata[,,ii] tmp1-as.data.frame(tmp) names(tmp1)=c(persons, d1, tp, fn, fp, fn, detect, d0, outcome) lm1-tryCatch(lmer(outcome~0+d1+d0+(0+d1+d0 | persons), family=binomial, data=tmp1, nAGQ=3), error=function(e) e) a[ii,1]=lm1@fixef[1] a[ii,2]=lm1@fixef[2] a[ii,3]=vcov(lm1)[1,2]/prod(sqrt(diag(vcov(lm1 a[ii,4:5]=sqrt(diag(vcov(lm1))) } return(a) } Any guidance would be greatly appreciated... thanks! Lisa -- Lisa Wang email: wang.li...@gmail.com cell: +49 -0176-87786557 Tübingen, Germany, 72070 [[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] biclustering error in drawHeatmap
Hi Everyone, I have a problem that keeps occurring using the biclust function (info below). Does anyone have an idea what is going wrong? Best Wishes and thank you, Bart library(biclust) data-read.table(file.choose(), header = T, sep = \t) data2-as.matrix(data[,-1]) data3-discretize(data2) dim(data3) [1] 266 30 clust-biclust(data3,method = BCCC(), number = 10) clust An object of class Biclust call: biclust(x = data3, method = BCCC(), number = 10) Number of Clusters found: 10 First 5 Cluster sizes: BC 1 BC 2 BC 3 BC 4 BC 5 Number of Rows: 82 58 39 19 15 Number of Columns: 20 15 19 16 17 drawHeatmap(x = data3, bicResult = clust, bicluster = 10) Error in drawHeatmap(x = data3, bicResult = clust, bicluster = 10) : Error: the bicluster does not exist in the result set [[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] Problem with generated parameter estimates
Dear All, I would be very grateful for your help concerning the following question: Below mentioned programme is available on net to generate longitudinal data. Usually we get almost same parameter estimates as used to generate the data. The problem here is I am not able to get it for data used here, despite increasing sample size and number of simulations. Is it normal to expect this type of behavior from mixed effect models? The problem is this code running fine for data(#data) provided by the author. I am not good with programming, can anyone comments on the code or possible explanation for this behavior? Any help would be greatly appreciated. sim.lmer - function(n,p,error,B0,varb0,B1,varb1,cor) { #start function ### ### INPUTS ### #n is number of subjects #p is number of time points #error is Residuals of measurement #B0 is fixed intercept effect (average group intercept) #B1 is fixed slope effect (average group slope) #varb0 is the variance of individual intercepts #varb1 is the variance of individual slopes #cor is the correlation between ind. intercepts and ind. slopes # ASSUMPTIONS### # Random effects have a bivariate normal distribution with zero mean # Random errors have a normal distribution with zero mean # The responses are generated based on a linear mixed effect regression model # Y = B0 + B1*Weeks + b0 + b1 + e # # # # Start function cov - cor*sqrt(varb0*varb1) #correlation between random effects d - matrix(c(varb0,cov,cov,varb1),nrow=2,ncol=2) #var-cov matrix of random effects ind - mvrnorm(n,c(0,0),d) #generate bivariate normal random effects b0 - ind[,1] # individual intercepts' deviation from fixed intercept b1 - ind[,2] # individual slopes' deviation from fixed slope d2 - (error^2)*diag(p) # var-cov matrix of error terms at time points err - mvrnorm(n,rep(0,p),d2) # generate multivariate normal error terms with zero mean ind.slo - B1+b1 ind.int - B0+b0 rand.eff - cbind(ind.slo,ind.int) ## data - matrix(nrow=n,ncol=p) for(i in 1:n) { for(k in 1:p) { data[i,k]= B0 + B1*(k-1) + (b0[i] + b1[i]*(k-1)) }} data2 - data + err ## data2 - as.data.frame(data2) names - c() for(i in 1:p) names[i]=paste(Score,i,sep=) colnames(data2) - names ### Add column names to data mynames - colnames(data2) data2$ID - 1:n d - reshape(data2,varying=mynames,idvar=ID, v.names=Score,timevar=Weeks,times=1:p,direction=long) list(data=d,rand.eff=rand.eff) } ### End of Function require(MASS) require(plyr) require(lme4) data - sim.lmer(n=1000,p=6,error=10,B0=38.93,varb0=30.92,B1=-2.31,varb1=0.56,cor=-0.7) #data - sim.lmer(n=1000,p=10,error=5,B0=40,varb0=1050,B1=.2,varb1=.4,cor=.7) d - data$data LMER.1 - lmer(Score ~ Weeks + (1 + Weeks | ID), data = d, REML=F) summary(LMER.1) -- Thanks Regards, Kamal Kishore PhD Student [[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] date time manipulation- R 2.15.1 windows 7
Hi Arun, Thank u so much for your help. It helped me solve the problem THANK YOU -Yashvardhan Kajaria On Thu, Mar 14, 2013 at 11:28 PM, arun smartpink...@yahoo.com wrote: HI, 1. date1-c(5 jan 2013, 1 jan 2013) date1-as.Date(date1,format=%d %b %Y) date1[1]-date1[2] #Time difference of 4 days 2. If you only have the week number of year without any other information, it would be difficult to predict which day that would be. You could get the week number from the date: library(lubridate) date2-2013-01-26 wday(ymd(date2),label=TRUE) # 1 parsed with %Y-%m-%d #[1] Sat week(ymd(date2)) # 1 parsed with %Y-%m-%d #[1] 4 A.K. - Original Message - From: yash kajaria yash.kaja...@gmail.com To: r-help@r-project.org Cc: Sent: Thursday, March 14, 2013 9:45 AM Subject: [R] date time manipulation- R 2.15.1 windows 7 Hi, I wanted to learn how to solve a date and time manipulation where i can do the following two 1. difference of two dates eg (differnce between 5th jan 2013 and 1st jan 2013) 2.Suppose i have week number of the year, i want to know if i can find out the day it refers to eg( say week 2 of 2013 would be 6th jan 2013 and the day is sunday) i need my result to tell me that its the 6th of jan 2013 as well as the day (sunday) Can u please help me out? Thanks, Yashvardhan Kajaria [[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] save scores from sem
I'm not aware of any routine that those the job, although I think that it could be relatively easily done by multiplication the manifest variable vector with the estimates for the specific effect. To make an example: v1; v2; v3; v4 are manifest variables that loads on one y latent variablein a data frame called A the code for the model should be like: model -specifymodel( y - v1, lam1, NA y - v2, lam2, NA y - v3, lam3, NA y - v4, lam4, NA After fitting the model with sem model.sem - sem(model, data=A) you should be able to compute the y variable like: attach(data) data$y-v1*lam1+v2*lam2+v3*lam3+v4*lam4 #change the loading name with the actual loading (number) or extract them from the objectiveML object (they are located in model.sem[[15]]) Note that those loadings are unstandardized and that the resulting variable will not be standardized. Hope it helps Regards, Marko -- Marko Tončić Assistant Researcher University of Rijeka Faculty of Humanities and Social Sciences Department of Psychology Sveučilišna Avenija 4, 51000 Rijeka, Croatia __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Assign Observations to Tree Nodes
Hi, is there a way to assign observations from a dataset to nodes of a tree?  For example, if I have a data frame like the following: ID Var1 Var2 Var3 ⦠VarX 1      .      .      .      n      And then run a classification tree with this form: tree - ctree(Var1 ~ Var2 + Var3 + ... + VarX, data=set) I would like to produce a vector for each node containing the observations (IDs) that are contained within that node.  Or really, any output that would tell me what observations are contained within nodes.  My goal is to run logits/other analyses on certain nodes and be able to easily pinpoint subsets of the original data frame who meet the criteria of certain nodes. Thanks in advance for any help/suggestions. Jeff [[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] Counting confidence intervals
Hi, Try this: set.seed(25) mat1- matrix(cbind(sample(1:15,20,replace=TRUE),sample(16:30,20,replace=TRUE)),ncol=2) nrow(mat1[sapply(seq_len(nrow(mat1)),function(i) any(seq(mat1[i,1],mat1[i,2])==12)),]) #[1] 17 set.seed(25) mat2- matrix(cbind(sample(1:15,1e5,replace=TRUE),sample(16:30,1e5,replace=TRUE)),ncol=2) system.time(res-nrow(mat2[sapply(seq_len(nrow(mat2)),function(i) any(seq(mat2[i,1],mat2[i,2])==12)),])) # user system elapsed # 1.552 0.000 1.549 res #[1] 80070 head(mat2[sapply(seq_len(nrow(mat2)),function(i) any(seq(mat2[i,1],mat2[i,2])==12)),]) # [,1] [,2] #[1,] 7 29 #[2,] 11 30 #[3,] 3 30 #[4,] 2 26 #[5,] 10 22 #[6,] 6 22 A.K. From: Jim Silverton jim.silver...@gmail.com To: r-help@r-project.org Sent: Monday, March 18, 2013 9:03 AM Subject: Re: [R] Counting confidence intervals Hi, I have a 2 x 1 matrix of confidence intervals. The first column is the lower and the next column is the upper. I want to cont how many times a number say 12 lies in the interval. Can anyone assist? -- Thanks, Jim. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Problem with write.table
Hi everyone, I'm trying to create unique filenames and then write a data frame into these files. This is the code I am using. str1-fname str2-.ext.txt; FILENAME-paste0(str1,str2); write.table(df, file=FILENAME,col.names=FALSE,row.names=FALSE,quote=FALSE,sep=\t); Ideally, a file called fname.ext.txt should be created, containing the data frame df. However, the file is not getting created at all. Am I going wrong somewhere? [[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] Counting confidence intervals
Thanks. Jeff On Mon, Mar 18, 2013 at 9:30 AM, Jeff Newmiller jdnew...@dcn.davis.ca.uswrote: sum(M[1]12 12=M[2]) untested, no data --- 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 --- Sent from my phone. Please excuse my brevity. Jim Silverton jim.silver...@gmail.com wrote: Hi, I have a 2 x 1 matrix of confidence intervals. The first column is the lower and the next column is the upper. I want to cont how many times a number say 12 lies in the interval. Can anyone assist? -- Thanks, Jim. [[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] Counting confidence intervals
Thats cumbersome, Arun. sum(mat1[,1] 12 mat1[,2] 12) [1] 17 will do the job and even faster: system.time(replicate(1, sum(mat1[,1] 12 mat1[,2] 12))) # user system elapsed # 0.067 0.001 0.078 HTH, Jorge.- On Tue, Mar 19, 2013 at 1:06 AM, arun wrote: Hi, Try this: set.seed(25) mat1- matrix(cbind(sample(1:15,20,replace=TRUE),sample(16:30,20,replace=TRUE)),ncol=2) nrow(mat1[sapply(seq_len(nrow(mat1)),function(i) any(seq(mat1[i,1],mat1[i,2])==12)),]) #[1] 17 set.seed(25) mat2- matrix(cbind(sample(1:15,1e5,replace=TRUE),sample(16:30,1e5,replace=TRUE)),ncol=2) system.time(res-nrow(mat2[sapply(seq_len(nrow(mat2)),function(i) any(seq(mat2[i,1],mat2[i,2])==12)),])) # user system elapsed # 1.552 0.000 1.549 res #[1] 80070 head(mat2[sapply(seq_len(nrow(mat2)),function(i) any(seq(mat2[i,1],mat2[i,2])==12)),]) # [,1] [,2] #[1,]7 29 #[2,] 11 30 #[3,]3 30 #[4,]2 26 #[5,] 10 22 #[6,]6 22 A.K. From: Jim Silverton To: r-help@r-project.org Sent: Monday, March 18, 2013 9:03 AM Subject: Re: [R] Counting confidence intervals Hi, I have a 2 x 1 matrix of confidence intervals. The first column is the lower and the next column is the upper. I want to cont how many times a number say 12 lies in the interval. Can anyone assist? -- Thanks, Jim. [[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. [[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] Counting confidence intervals
Hi, Jorge's method will be faster. #system.time(res1-sum(apply(mat2,1,function(x) x[1]12 x[2]12))) #instead of 2, it should be 1 # user system elapsed # 0.440 0.000 0.445 system.time(res1-sum(apply(mat2,1,function(x) x[1]=12 x[2]12))) # # user system elapsed # 0.500 0.000 0.502 res1 #[1] 80070 A.K. From: Jim Silverton jim.silver...@gmail.com To: arun smartpink...@yahoo.com Sent: Monday, March 18, 2013 10:08 AM Subject: Re: [R] Counting confidence intervals thanks arun!! On Mon, Mar 18, 2013 at 10:06 AM, arun smartpink...@yahoo.com wrote: Hi, Try this: set.seed(25) mat1- matrix(cbind(sample(1:15,20,replace=TRUE),sample(16:30,20,replace=TRUE)),ncol=2) nrow(mat1[sapply(seq_len(nrow(mat1)),function(i) any(seq(mat1[i,1],mat1[i,2])==12)),]) #[1] 17 set.seed(25) mat2- matrix(cbind(sample(1:15,1e5,replace=TRUE),sample(16:30,1e5,replace=TRUE)),ncol=2) system.time(res-nrow(mat2[sapply(seq_len(nrow(mat2)),function(i) any(seq(mat2[i,1],mat2[i,2])==12)),])) # user system elapsed # 1.552 0.000 1.549 res #[1] 80070 head(mat2[sapply(seq_len(nrow(mat2)),function(i) any(seq(mat2[i,1],mat2[i,2])==12)),]) # [,1] [,2] #[1,] 7 29 #[2,] 11 30 #[3,] 3 30 #[4,] 2 26 #[5,] 10 22 #[6,] 6 22 A.K. From: Jim Silverton jim.silver...@gmail.com To: r-help@r-project.org Sent: Monday, March 18, 2013 9:03 AM Subject: Re: [R] Counting confidence intervals Hi, I have a 2 x 1 matrix of confidence intervals. The first column is the lower and the next column is the upper. I want to cont how many times a number say 12 lies in the interval. Can anyone assist? -- Thanks, Jim. [[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. -- Thanks, Jim. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Counting confidence intervals
If you don't use apply() it would be even faster: system.time(sum(mat2[,1] 12 mat2[,2] 12)) user system elapsed 0.004 0.000 0.003 Regards, Jorge.- On Tue, Mar 19, 2013 at 1:21 AM, arun wrote: Hi, Jorge's method will be faster. #system.time(res1-sum(apply(mat2,1,function(x) x[1]12 x[2]12))) #instead of 2, it should be 1 # user system elapsed # 0.440 0.000 0.445 system.time(res1-sum(apply(mat2,1,function(x) x[1]=12 x[2]12))) # # user system elapsed # 0.500 0.000 0.502 res1 #[1] 80070 A.K. From: Jim Silverton jim.silver...@gmail.com To: arun smartpink...@yahoo.com Sent: Monday, March 18, 2013 10:08 AM Subject: Re: [R] Counting confidence intervals thanks arun!! On Mon, Mar 18, 2013 at 10:06 AM, arun smartpink...@yahoo.com wrote: Hi, Try this: set.seed(25) mat1- matrix(cbind(sample(1:15,20,replace=TRUE),sample(16:30,20,replace=TRUE)),ncol=2) nrow(mat1[sapply(seq_len(nrow(mat1)),function(i) any(seq(mat1[i,1],mat1[i,2])==12)),]) #[1] 17 set.seed(25) mat2- matrix(cbind(sample(1:15,1e5,replace=TRUE),sample(16:30,1e5,replace=TRUE)),ncol=2) system.time(res-nrow(mat2[sapply(seq_len(nrow(mat2)),function(i) any(seq(mat2[i,1],mat2[i,2])==12)),])) # user system elapsed # 1.552 0.000 1.549 res #[1] 80070 head(mat2[sapply(seq_len(nrow(mat2)),function(i) any(seq(mat2[i,1],mat2[i,2])==12)),]) # [,1] [,2] #[1,]7 29 #[2,] 11 30 #[3,]3 30 #[4,]2 26 #[5,] 10 22 #[6,]6 22 A.K. From: Jim Silverton jim.silver...@gmail.com To: r-help@r-project.org Sent: Monday, March 18, 2013 9:03 AM Subject: Re: [R] Counting confidence intervals Hi, I have a 2 x 1 matrix of confidence intervals. The first column is the lower and the next column is the upper. I want to cont how many times a number say 12 lies in the interval. Can anyone assist? -- Thanks, Jim. [[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. -- Thanks, Jim. [[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] How many samples ACTUALLY used in regression?
On Mar 18, 2013, at 7:36 AM, Federico Calboli f.calb...@imperial.ac.uk wrote: Dear All, is there a simple way that covers all regression models to extract the number of samples from a data frame/matrix actually used in a regression model? For instance I might have a data of 100 rows and 4 colums (1 response + 3 explanatory variables). If 3 samples have one or more NAs in the explanatory variable columns these samples will be dropped in any model: my.model = lm(y ~ x + w + z, my.data) my.model = glm(y ~ x + w + z, my.data, family = binomial) my.model = polr(y ~ x + w + z, my.data) … I don't seem to be able to find one single method that works in the exact same way -- irrespective of the model type -- to interrogate my.model to see how many samples of my.data were actually used. Is there such function or do I need to hack something together? Best wishes Federico I don't know that this would be universal to all possible R model implementations, but should work for those that at least abide by certain standards[1] relative to the internal use of ?model.frame. In the case where model functions use 'model = TRUE' as the default in their call (eg. lm(), glm() and MASS::polr()), the returned model object will have a component called 'model', such that: nrow(my.model$model) returns the number of rows in the internally created data frame. Note that 'model = TRUE' is not the default for many functions, for example Terry's coxph() in survival or Frank's lrm() in rms. Note also that the value of 'na.action' in the modeling function call may have a potential effect on whether the number of rows in the retained 'model' data frame is really the correct value. You can also use model.frame(), independently matching arguments passed to the model function, to replicate what takes place internally in many modeling functions. The result of model.frame() will be a data frame, again, subject to similar limitations as above. Regards, Marc Schwartz [1]: http://developer.r-project.org/model-fitting-functions.txt __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] save scores from sem
Dear Marko, I can't quite tell whether this is an original question (as suggested by no Re: in the subject field and no text from an earlier message) or an answer to a question already posed (as suggested by the phrasing of the message); if the latter, I apologize for missing the original question. In any event, see the fscore() function in the sem package (?fscore) for computation of factor scores. I hope this helps, John John Fox Sen. William McMaster Prof. of Social Statistics Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox/ On Mon, 18 Mar 2013 13:21:48 +0100 marKo mton...@ffri.hr wrote: I'm not aware of any routine that those the job, although I think that it could be relatively easily done by multiplication the manifest variable vector with the estimates for the specific effect. To make an example: v1; v2; v3; v4 are manifest variables that loads on one y latent variablein a data frame called A the code for the model should be like: model -specifymodel( y - v1, lam1, NA y - v2, lam2, NA y - v3, lam3, NA y - v4, lam4, NA After fitting the model with sem model.sem - sem(model, data=A) you should be able to compute the y variable like: attach(data) data$y-v1*lam1+v2*lam2+v3*lam3+v4*lam4 #change the loading name with the actual loading (number) or extract them from the objectiveML object (they are located in model.sem[[15]]) Note that those loadings are unstandardized and that the resulting variable will not be standardized. Hope it helps Regards, Marko -- Marko Ton?i? Assistant Researcher University of Rijeka Faculty of Humanities and Social Sciences Department of Psychology Sveu?ilina Avenija 4, 51000 Rijeka, Croatia __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Problem with write.table
Hi -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of Sahana Srinivasan Sent: Monday, March 18, 2013 3:04 PM To: r-help@r-project.org Subject: [R] Problem with write.table Hi everyone, I'm trying to create unique filenames and then write a data frame into these files. This is the code I am using. str1-fname str2-.ext.txt; FILENAME-paste0(str1,str2); write.table(df, file=FILENAME,col.names=FALSE,row.names=FALSE,quote=FALSE,sep=\t); Ideally, a file called fname.ext.txt should be created, containing the data frame df. However, the file is not getting created at all. Am I going wrong somewhere? Maybe OS and write permissions? I tried your code with some data.frame and I got a file fname.ext.txt without problems. Check your objecd df and FILENAME. Regards Petr [[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] How many samples ACTUALLY used in regression?
Perhaps a crude but reliable way is to check the number of residuals, e.g., length(my.model$resid). Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: ca...@usgs.gov brian_c...@usgs.gov tel: 970 226-9326 On Mon, Mar 18, 2013 at 8:39 AM, Marc Schwartz marc_schwa...@me.com wrote: On Mar 18, 2013, at 7:36 AM, Federico Calboli f.calb...@imperial.ac.uk wrote: Dear All, is there a simple way that covers all regression models to extract the number of samples from a data frame/matrix actually used in a regression model? For instance I might have a data of 100 rows and 4 colums (1 response + 3 explanatory variables). If 3 samples have one or more NAs in the explanatory variable columns these samples will be dropped in any model: my.model = lm(y ~ x + w + z, my.data) my.model = glm(y ~ x + w + z, my.data, family = binomial) my.model = polr(y ~ x + w + z, my.data) I don't seem to be able to find one single method that works in the exact same way -- irrespective of the model type -- to interrogate my.model to see how many samples of my.data were actually used. Is there such function or do I need to hack something together? Best wishes Federico I don't know that this would be universal to all possible R model implementations, but should work for those that at least abide by certain standards[1] relative to the internal use of ?model.frame. In the case where model functions use 'model = TRUE' as the default in their call (eg. lm(), glm() and MASS::polr()), the returned model object will have a component called 'model', such that: nrow(my.model$model) returns the number of rows in the internally created data frame. Note that 'model = TRUE' is not the default for many functions, for example Terry's coxph() in survival or Frank's lrm() in rms. Note also that the value of 'na.action' in the modeling function call may have a potential effect on whether the number of rows in the retained 'model' data frame is really the correct value. You can also use model.frame(), independently matching arguments passed to the model function, to replicate what takes place internally in many modeling functions. The result of model.frame() will be a data frame, again, subject to similar limitations as above. Regards, Marc Schwartz [1]: http://developer.r-project.org/model-fitting-functions.txt __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Confirmatory factor analysis using the sem package. TLI CFI and RMSEA absent from model summary.
Hi R-help, I am using the sem package to run confirmatory factor analysis (cfa) on some questionnaire data collected from 307 participants. I have been running R-2.15.3 in Windows in conjunction with R studio. The model I am using was developed from exploratory factor analysis of a separate dataset (n=439); it includes 18 items that load onto 3 factors. I have used the sem package documentation and this video (http://vimeo.com/38941937) to run the cfa and obtain a chi-square statistic for the model. However, when I use the summary() function, the model does not provide indices of fit; I need these as part of my analysis output. In particular, I am looking for the Tucker Lewis Index (TLI), Comparative Fit Index (CFI), the Root Mean Square of Approximation (RMSEA). I have checked the documentation and cannot seem to find any reason for this; none of the arguments listed with the sem command indicate that I have to specify these as part of the output. In addition, the analysis example f! rom the video includes these indices as part of the output, but my analysis does not provide them. I have included my code with comments below: library(sem) validation.data - structure(list(V1 = c(5L, 4L, 2L, 4L, 5L, 6L, 6L, 4L, 5L, 3L, 6L, 5L, 4L, 5L, 3L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 4L, 5L, 4L, 4L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 2L, 6L, 5L, 6L, 4L, 5L, 6L, 5L, 5L, 4L, 5L, 5L, 3L, 5L, 5L, 5L, 5L, 4L, 2L, 5L, 5L, 5L, 4L, 6L, 4L, 6L, 5L, 5L, 5L, 4L, 5L, 5L, 4L, 5L, 6L, 4L, 5L, 4L, 5L, 5L, 5L, 3L, 5L, 5L, 3L, 5L, 4L, 5L, 2L, 6L, 4L, 4L, 4L, 5L, 5L, 4L, 6L, 2L, 4L, 5L, 4L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 2L, 4L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 6L, 5L, 5L, 4L, 4L, 5L, 4L, 4L, 5L, 4L, 4L, 4L, 5L, 5L, 6L, 5L, 5L, 5L, 4L, 4L, 3L, 4L, 5L, 5L, 5L, 2L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 6L, 4L, 5L, 5L, 6L, 5L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 4L, 4L, 4L, 5L, 2L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 4L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 4L, 4L, 2L, 6L, 5L, 5L, 5L, 6L, 4L, 3L, 5L, 5L, 5L, 5L, 4L, 4L, 6L, 6L, 5L, 5L, 5L, 3L, 6L, 5L, 5L, 5L, 4L, 5L, 5L, 4L, 5L, 4L, 6L, 5L, 4L, 4L, 5L, 4L, 3L, 4L, 5L, 4L, 5L, 6L, 2L, 4L, 4L, 5L, 4L, 4L, 4L, 5L, 6L, 5L, 4L, 4L, 4L, 6L, 6L, 4L, 5L, 5L, 5L, 2L, 4L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 3L, 5L, 6L, 5L, 4L, 4L, 3L, 3L, 4L, 5L, 5L, 1L, 4L, 5L, 3L, 5L, 1L, 6L, 5L, 4L, 4L, 5L, 5L, 4L, 5L, 5L, 6L, 5L, 5L, 5L), V2 = c(5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 5L, 6L, 6L, 4L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 4L, 6L, 5L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 5L, 4L, 6L, 5L, 6L, 5L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 5L, 4L, 6L, 4L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 4L, 2L, 5L, 6L, 4L, 5L, 5L, 6L, 6L, 5L, 6L, 4L, 6L, 5L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 4L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 3L, 5L, 6L, 6L, 6L, 5L, 6L, 5L, 5L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 4L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 5L, 6L, 5L, 5L, 6L, 5L, 5L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 6L, 6L), V3 = c(5L, 5L, 3L, 6L, 5L, 2L, 4L, 4L, 4L, 4L, 6L, 5L, 5L, 4L, 4L, 4L, 4L, 5L, 4L, 4L, 1L, 3L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 1L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 4L, 3L, 2L, 5L, 6L, 6L, 5L, 5L, 5L, 6L, 5L, 6L, 5L, 4L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 3L, 5L, 4L, 5L, 5L, 5L, 3L, 5L, 5L, 5L, 3L, 5L, 4L, 5L, 4L, 5L, 4L, 3L, 5L, 3L, 5L, 3L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 2L, 5L, 4L, 5L, 4L, 6L, 4L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 6L, 5L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 4L, 4L, 5L, 5L, 4L, 6L, 4L, 6L, 5L, 4L, 4L, 5L, 5L, 5L, 4L, 4L, 1L, 5L, 4L, 4L, 5L, 5L, 4L, 6L, 3L, 4L, 4L, 5L, 4L, 4L, 5L, 3L, 5L, 6L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 3L, 5L, 4L, 6L, 5L, 4L, 6L, 3L, 4L, 2L, 4L, 4L, 5L, 4L, 4L, 5L, 3L, 4L, 5L, 5L, 4L, 3L, 3L, 5L, 5L, 5L, 5L, 4L, 3L, 4L, 2L, 5L, 5L, 6L, 4L, 5L, 4L, 4L, 5L, 4L, 6L, 6L, 4L, 4L, 4L, 5L, 5L, 6L, 5L, 4L, 6L, 5L, 5L, 5L, 4L, 6L, 6L, 3L, 2L, 3L, 6L, 4L, 5L, 3L, 6L, 3L, 4L, 4L, 5L, 4L, 6L, 4L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 4L, 5L, 5L, 6L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 4L, 3L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 6L, 5L, 4L, 4L, 3L, 4L, 5L, 5L, 4L, 1L, 6L, 4L, 4L, 4L, 2L, 6L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 4L, 6L, 5L, 5L, 6L), V4 = c(5L, 3L, 4L, 6L, 5L, 4L, 6L, 4L,
Re: [R] How many samples ACTUALLY used in regression?
On 18/03/2013 14:51, Cade, Brian wrote: Perhaps a crude but reliable way is to check the number of residuals, e.g., length(my.model$resid). Not very reliable (what about zero weights, for example?), and the component is usually 'residuals'. No one has so far mentioned nobs(), which seems to me to be the closest. Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: ca...@usgs.gov brian_c...@usgs.gov tel: 970 226-9326 On Mon, Mar 18, 2013 at 8:39 AM, Marc Schwartz marc_schwa...@me.com wrote: On Mar 18, 2013, at 7:36 AM, Federico Calboli f.calb...@imperial.ac.uk wrote: Dear All, is there a simple way that covers all regression models to extract the number of samples from a data frame/matrix actually used in a regression model? For instance I might have a data of 100 rows and 4 colums (1 response + 3 explanatory variables). If 3 samples have one or more NAs in the explanatory variable columns these samples will be dropped in any model: my.model = lm(y ~ x + w + z, my.data) my.model = glm(y ~ x + w + z, my.data, family = binomial) my.model = polr(y ~ x + w + z, my.data) … I don't seem to be able to find one single method that works in the exact same way -- irrespective of the model type -- to interrogate my.model to see how many samples of my.data were actually used. Is there such function or do I need to hack something together? Best wishes Federico I don't know that this would be universal to all possible R model implementations, but should work for those that at least abide by certain standards[1] relative to the internal use of ?model.frame. In the case where model functions use 'model = TRUE' as the default in their call (eg. lm(), glm() and MASS::polr()), the returned model object will have a component called 'model', such that: nrow(my.model$model) returns the number of rows in the internally created data frame. Note that 'model = TRUE' is not the default for many functions, for example Terry's coxph() in survival or Frank's lrm() in rms. Note also that the value of 'na.action' in the modeling function call may have a potential effect on whether the number of rows in the retained 'model' data frame is really the correct value. You can also use model.frame(), independently matching arguments passed to the model function, to replicate what takes place internally in many modeling functions. The result of model.frame() will be a data frame, again, subject to similar limitations as above. Regards, Marc Schwartz [1]: http://developer.r-project.org/model-fitting-functions.txt __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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. -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How many samples ACTUALLY used in regression?
On 18 Mar 2013, at 15:07, Prof Brian Ripley rip...@stats.ox.ac.uk wrote: On 18/03/2013 14:51, Cade, Brian wrote: Perhaps a crude but reliable way is to check the number of residuals, e.g., length(my.model$resid). Not very reliable (what about zero weights, for example?), and the component is usually 'residuals'. No one has so far mentioned nobs(), which seems to me to be the closest. Given a my.data where 3 out of 100 rows will be discarded due to NAs test = lm(formula = y ~ x + w, my.data, model = T) nobs(test) [1] 97 # as expected But if I substitute 1 NA in one of the row of the housing data: house.plr = polr(formula = Sat ~ Infl + Type + Cont, data = housing, weights = Freq) nobs(house.plr) [1] 1661 because of weights (which would not be take into account in a glm() fit). Because I only care about the raw number of observations, is there a (hopefully) trivial way of getting nobs(poor.fit) to behave like a nobs(vlm.fit)? BW Federico Brian Brian S. Cade, PhD U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: ca...@usgs.gov brian_c...@usgs.gov tel: 970 226-9326 On Mon, Mar 18, 2013 at 8:39 AM, Marc Schwartz marc_schwa...@me.com wrote: On Mar 18, 2013, at 7:36 AM, Federico Calboli f.calb...@imperial.ac.uk wrote: Dear All, is there a simple way that covers all regression models to extract the number of samples from a data frame/matrix actually used in a regression model? For instance I might have a data of 100 rows and 4 colums (1 response + 3 explanatory variables). If 3 samples have one or more NAs in the explanatory variable columns these samples will be dropped in any model: my.model = lm(y ~ x + w + z, my.data) my.model = glm(y ~ x + w + z, my.data, family = binomial) my.model = polr(y ~ x + w + z, my.data) … I don't seem to be able to find one single method that works in the exact same way -- irrespective of the model type -- to interrogate my.model to see how many samples of my.data were actually used. Is there such function or do I need to hack something together? Best wishes Federico I don't know that this would be universal to all possible R model implementations, but should work for those that at least abide by certain standards[1] relative to the internal use of ?model.frame. In the case where model functions use 'model = TRUE' as the default in their call (eg. lm(), glm() and MASS::polr()), the returned model object will have a component called 'model', such that: nrow(my.model$model) returns the number of rows in the internally created data frame. Note that 'model = TRUE' is not the default for many functions, for example Terry's coxph() in survival or Frank's lrm() in rms. Note also that the value of 'na.action' in the modeling function call may have a potential effect on whether the number of rows in the retained 'model' data frame is really the correct value. You can also use model.frame(), independently matching arguments passed to the model function, to replicate what takes place internally in many modeling functions. The result of model.frame() will be a data frame, again, subject to similar limitations as above. Regards, Marc Schwartz [1]: http://developer.r-project.org/model-fitting-functions.txt __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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. -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] try/tryCatch
It would help if you told us what type of error you are getting and to also provide sample data so that we could run it to see what happens. I use 'try' a lot to catch errors and have not had any problems with it. On Mon, Mar 18, 2013 at 6:11 AM, lisa wang.li...@gmail.com wrote: Hi All, I have tried every fix on my try or tryCatch that I have found on the internet, but so far have not been able to get my R code to continue with the for loop after the lmer model results in an error. Here is two attemps of my code, the input is a 3D array file, but really any function would do metatrialstry-function(mydata){ a-matrix(data=NA, nrow=dim(mydata)[3], ncol=5) #colnames(a)=c(sens, spec, corr, sens_se, spec_se, counter)#colnames(a)=c(sens, spec, corr, sens_se, spec_se) k=1 for(ii in 1:dim(mydata)[3]){ tmp-mydata[,,ii] tmp1-as.data.frame(tmp) names(tmp1)=c(persons, d1, tp, fn, fp, fn, detect, d0, outcome) lm1-try(lmer(outcome~0+d1+d0+(0+d1+d0 | persons), family=binomial, data=tmp1, nAGQ=3), silent=T) if(class(lm1)[1]!='try-error'){ a[ii,1]=lm1@fixef[1] a[ii,2]=lm1@fixef[2] a[ii,3]=vcov(lm1)[1,2]/prod(sqrt(diag(vcov(lm1 a[ii,4:5]=sqrt(diag(vcov(lm1))) } } #k=k+1 #a[ii,6]=k return(a) } # # try / try catch ### # metatrialstry-function(mydata){ a-matrix(data=NA, nrow=dim(mydata)[3], ncol=5) #colnames(a)=c(sens, spec, corr, sens_se, spec_se, counter) colnames(a)=c(sens, spec, corr, sens_se, spec_se) #a[,6]=rep(0, length(a[,6])) for(ii in 1:dim(mydata)[3]){ tmp-mydata[,,ii] tmp1-as.data.frame(tmp) names(tmp1)=c(persons, d1, tp, fn, fp, fn, detect, d0, outcome) lm1-tryCatch(lmer(outcome~0+d1+d0+(0+d1+d0 | persons), family=binomial, data=tmp1, nAGQ=3), error=function(e) e) a[ii,1]=lm1@fixef[1] a[ii,2]=lm1@fixef[2] a[ii,3]=vcov(lm1)[1,2]/prod(sqrt(diag(vcov(lm1 a[ii,4:5]=sqrt(diag(vcov(lm1))) } return(a) } Any guidance would be greatly appreciated... thanks! Lisa -- Lisa Wang email: wang.li...@gmail.com cell: +49 -0176-87786557 Tübingen, Germany, 72070 [[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. -- 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. [[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] Confirmatory factor analysis using the sem package. TLI CFI and RMSEA absent from model summary.
Dear Kevin, See ?summary.objectiveML, and in particular the description of the fit.indices argument. By default, the summary() method doesn't print many fit indices, but many are available optionally. I hope this helps, John John Fox Sen. William McMaster Prof. of Social Statistics Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox/ On Mon, 18 Mar 2013 15:00:06 + Kevin Cheung k.che...@derby.ac.uk wrote: Hi R-help, I am using the sem package to run confirmatory factor analysis (cfa) on some questionnaire data collected from 307 participants. I have been running R-2.15.3 in Windows in conjunction with R studio. The model I am using was developed from exploratory factor analysis of a separate dataset (n=439); it includes 18 items that load onto 3 factors. I have used the sem package documentation and this video (http://vimeo.com/38941937) to run the cfa and obtain a chi-square statistic for the model. However, when I use the summary() function, the model does not provide indices of fit; I need these as part of my analysis output. In particular, I am looking for the Tucker Lewis Index (TLI), Comparative Fit Index (CFI), the Root Mean Square of Approximation (RMSEA). I have checked the documentation and cannot seem to find any reason for this; none of the arguments listed with the sem command indicate that I have to specify these as part of the output. In addition, the analysis example! f! rom the video includes these indices as part of the output, but my analysis does not provide them. I have included my code with comments below: library(sem) validation.data - structure(list(V1 = c(5L, 4L, 2L, 4L, 5L, 6L, 6L, 4L, 5L, 3L, 6L, 5L, 4L, 5L, 3L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 4L, 5L, 4L, 4L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 2L, 6L, 5L, 6L, 4L, 5L, 6L, 5L, 5L, 4L, 5L, 5L, 3L, 5L, 5L, 5L, 5L, 4L, 2L, 5L, 5L, 5L, 4L, 6L, 4L, 6L, 5L, 5L, 5L, 4L, 5L, 5L, 4L, 5L, 6L, 4L, 5L, 4L, 5L, 5L, 5L, 3L, 5L, 5L, 3L, 5L, 4L, 5L, 2L, 6L, 4L, 4L, 4L, 5L, 5L, 4L, 6L, 2L, 4L, 5L, 4L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 2L, 4L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 6L, 5L, 5L, 4L, 4L, 5L, 4L, 4L, 5L, 4L, 4L, 4L, 5L, 5L, 6L, 5L, 5L, 5L, 4L, 4L, 3L, 4L, 5L, 5L, 5L, 2L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 6L, 4L, 5L, 5L, 6L, 5L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 4L, 4L, 4L, 5L, 2L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 4L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 4L, 4L, 2L, 6L, 5L, 5L, 5L, 6L, 4L, 3L, 5L, 5L, 5L, 5L, 4L, 4L, 6L, 6L, 5L, 5L, 5L, 3L, 6L, 5L, 5L, 5L, 4L, 5L, 5L, 4L, 5L, 4L, 6L, 5L, 4L, 4L, 5L, 4L, 3L, 4L, 5L, 4L, 5L, 6L, 2L, 4L, 4L, 5L, 4L, 4L, 4L, 5L, 6L, 5L, 4L, 4L, 4L, 6L, 6L, 4L, 5L, 5L, 5L, 2L, 4L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 3L, 5L, 6L, 5L, 4L, 4L, 3L, 3L, 4L, 5L, 5L, 1L, 4L, 5L, 3L, 5L, 1L, 6L, 5L, 4L, 4L, 5L, 5L, 4L, 5L, 5L, 6L, 5L, 5L, 5L), V2 = c(5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 5L, 6L, 6L, 4L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 4L, 6L, 5L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 5L, 4L, 6L, 5L, 6L, 5L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 5L, 4L, 6L, 4L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 4L, 2L, 5L, 6L, 4L, 5L, 5L, 6L, 6L, 5L, 6L, 4L, 6L, 5L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 4L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 3L, 5L, 6L, 6L, 6L, 5L, 6L, 5L, 5L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 4L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 5L, 6L, 5L, 5L, 6L, 5L, 5L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 6L, 6L), V3 = c(5L, 5L, 3L, 6L, 5L, 2L, 4L, 4L, 4L, 4L, 6L, 5L, 5L, 4L, 4L, 4L, 4L, 5L, 4L, 4L, 1L, 3L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 1L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 4L, 3L, 2L, 5L, 6L, 6L, 5L, 5L, 5L, 6L, 5L, 6L, 5L, 4L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 3L, 5L, 4L, 5L, 5L, 5L, 3L, 5L, 5L, 5L, 3L, 5L, 4L, 5L, 4L, 5L, 4L, 3L, 5L, 3L, 5L, 3L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 2L, 5L, 4L, 5L, 4L, 6L, 4L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 6L, 5L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 4L, 4L, 5L, 5L, 4L, 6L, 4L, 6L, 5L, 4L, 4L, 5L, 5L, 5L, 4L, 4L, 1L, 5L, 4L, 4L, 5L, 5L, 4L, 6L, 3L, 4L, 4L, 5L, 4L, 4L, 5L, 3L,
Re: [R] Counting confidence intervals
I want to cont how many times a number say 12 lies in the interval. Can anyone assist? Has anyone else ever wished there was a moderately general 'inside' or 'within' function in R for this problem? For example, something that behaves more or less like within - function(x, interval=NULL, closed=c(TRUE, TRUE), lower=min(interval), upper=max(interval)) { #interval must be a length 2 vector #closed is taken in the order (lower, upper) #lower and upper may be vectors and will be recycled (by etc) if not of length length(x) low.comp - if(closed[1]) = else high.comp - if(closed[2]) = else do.call(low.comp, list(lower, x)) do.call(high.comp, list(upper, x)) } #Examples within(1:5, c(2,4)) within(1:5, c(2,4), closed=c(FALSE, TRUE)) within(1:5, lower=5:1, upper=10:14) S Ellison LGC *** This email and any attachments are confidential. Any use...{{dropped:8}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] !0 + !0 == !0 - !0
Ben Bolker bbolker at gmail.com writes: Maybe FAQ 7.31 was referred to not for its direct relevance but as a measure of the old-hand-ness of the people who will get the joke. !1i|!0 Chuck __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Confirmatory factor analysis using the sem package. TLI CFI and RMSEA absent from model summary.
Dear John, Thank you for taking the time to help me with this, I have been able to obtain fit indices using the information that you provided. Note to users searching archived R-help posts about this issue: The instructional video I looked at (http://vimeo.com/38941937) gave fit indices using the default summary() command without any additional arguments. This may have been due to a different version of R (I noticed that the instructor was using a mac based OS). With regards, Kevin Kevin Yet Fong Cheung, Bsc., MRes., MBPsS. Postgraduate Researcher Centre for Psychological Research University of Derby Kedleston Road Derby DE22 1GB k.che...@derby.ac.ukmailto:k.che...@derby.ac.uk 01332592081 http://derby.academia.edu/KevinCheung -Original Message- From: John Fox [mailto:j...@mcmaster.ca] Sent: 18 March 2013 15:55 To: Kevin Cheung Cc: r-help@r-project.org Subject: Re: [R] Confirmatory factor analysis using the sem package. TLI CFI and RMSEA absent from model summary. Dear Kevin, See ?summary.objectiveML, and in particular the description of the fit.indices argument. By default, the summary() method doesn't print many fit indices, but many are available optionally. I hope this helps, John John Fox Sen. William McMaster Prof. of Social Statistics Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox/ On Mon, 18 Mar 2013 15:00:06 + Kevin Cheung k.che...@derby.ac.uk wrote: Hi R-help, I am using the sem package to run confirmatory factor analysis (cfa) on some questionnaire data collected from 307 participants. I have been running R-2.15.3 in Windows in conjunction with R studio. The model I am using was developed from exploratory factor analysis of a separate dataset (n=439); it includes 18 items that load onto 3 factors. I have used the sem package documentation and this video (http://vimeo.com/38941937) to run the cfa and obtain a chi-square statistic for the model. However, when I use the summary() function, the model does not provide indices of fit; I need these as part of my analysis output. In particular, I am looking for the Tucker Lewis Index (TLI), Comparative Fit Index (CFI), the Root Mean Square of Approximation (RMSEA). I have checked the documentation and cannot seem to find any reason for this; none of the arguments listed with the sem command indicate that I have to specify these as part of the output. In addition, the analysis example! f! rom the video includes these indices as part of the output, but my analysis does not provide them. I have included my code with comments below: library(sem) validation.data - structure(list(V1 = c(5L, 4L, 2L, 4L, 5L, 6L, 6L, 4L, 5L, 3L, 6L, 5L, 4L, 5L, 3L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 4L, 5L, 4L, 4L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 2L, 6L, 5L, 6L, 4L, 5L, 6L, 5L, 5L, 4L, 5L, 5L, 3L, 5L, 5L, 5L, 5L, 4L, 2L, 5L, 5L, 5L, 4L, 6L, 4L, 6L, 5L, 5L, 5L, 4L, 5L, 5L, 4L, 5L, 6L, 4L, 5L, 4L, 5L, 5L, 5L, 3L, 5L, 5L, 3L, 5L, 4L, 5L, 2L, 6L, 4L, 4L, 4L, 5L, 5L, 4L, 6L, 2L, 4L, 5L, 4L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 2L, 4L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 6L, 5L, 5L, 4L, 4L, 5L, 4L, 4L, 5L, 4L, 4L, 4L, 5L, 5L, 6L, 5L, 5L, 5L, 4L, 4L, 3L, 4L, 5L, 5L, 5L, 2L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 6L, 4L, 5L, 5L, 6L, 5L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 4L, 4L, 4L, 5L, 2L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 4L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 4L, 4L, 2L, 6L, 5L, 5L, 5L, 6L, 4L, 3L, 5L, 5L, 5L, 5L, 4L, 4L, 6L, 6L, 5L, 5L, 5L, 3L, 6L, 5L, 5L, 5L, 4L, 5L, 5L, 4L, 5L, 4L, 6L, 5L, 4L, 4L, 5L, 4L, 3L, 4L, 5L, 4L, 5L, 6L, 2L, 4L, 4L, 5L, 4L, 4L, 4L, 5L, 6L, 5L, 4L, 4L, 4L, 6L, 6L, 4L, 5L, 5L, 5L, 2L, 4L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 3L, 5L, 6L, 5L, 4L, 4L, 3L, 3L, 4L, 5L, 5L, 1L, 4L, 5L, 3L, 5L, 1L, 6L, 5L, 4L, 4L, 5L, 5L, 4L, 5L, 5L, 6L, 5L, 5L, 5L), V2 = c(5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 5L, 6L, 6L, 4L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 4L, 6L, 5L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 5L, 4L, 6L, 5L, 6L, 5L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 5L, 4L, 6L, 4L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 4L, 2L, 5L, 6L, 4L, 5L, 5L, 6L, 6L, 5L, 6L, 4L, 6L, 5L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 4L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L,
Re: [R] !0 + !0 == !0 - !0
* Bert Gunter thagre.ore...@trar.pbz [2013-03-17 20:30:56 -0700]: I also think it fair to say that all (??) languages have these sorts of malapropisms due to operator precedence. Except for those languages which do _not_ have operator precedence. Like, e.g., Lisp. -- Sam Steingold (http://sds.podval.org/) on Ubuntu 12.04 (precise) X 11.0.11103000 http://www.childpsy.net/ http://openvotingconsortium.org http://jihadwatch.org http://palestinefacts.org http://mideasttruth.com http://camera.org DRM access management == prison freedom management. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] !0 + !0 == !0 - !0
Sam: Yes. Good point. (which is why my ?? was necessary). -- Bert On Mon, Mar 18, 2013 at 10:05 AM, Sam Steingold s...@gnu.org wrote: * Bert Gunter thagre.ore...@trar.pbz [2013-03-17 20:30:56 -0700]: I also think it fair to say that all (??) languages have these sorts of malapropisms due to operator precedence. Except for those languages which do _not_ have operator precedence. Like, e.g., Lisp. -- Sam Steingold (http://sds.podval.org/) on Ubuntu 12.04 (precise) X 11.0.11103000 http://www.childpsy.net/ http://openvotingconsortium.org http://jihadwatch.org http://palestinefacts.org http://mideasttruth.com http://camera.org DRM access management == prison freedom management. -- Bert Gunter Genentech Nonclinical Biostatistics Internal Contact Info: Phone: 467-7374 Website: http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm [[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 a hyperlink to a csv file
Besides what others have suggested, there are at some packages for writing Excel files directly, and either of those might have a way to write something that Excel will recognize as a hyperlink. xlsx, XLConnect are both Java-based; I would start with these. -Don -- Don MacQueen Lawrence Livermore National Laboratory 7000 East Ave., L-627 Livermore, CA 94550 925-423-1062 On 3/16/13 6:28 AM, Brian Smith bsmith030...@gmail.com wrote: Hi Marc, Thanks for the reply. The question is whether it is possible to write some text (with a hyperlink) to a csv file, such that when you open the file in excel, it shows the text as hyperlinked. I guess it boils down to whether there are any 'tags' that you can put in the csv/txt file so that excel recognizes it as hyperlinked text. Does that make sense? thanks! On Sat, Mar 16, 2013 at 4:16 AM, Marc Girondot marc_...@yahoo.fr wrote: Le 15/03/13 12:53, Brian Smith a écrit : Hi, I was wondering if it is possible to create a hyperlink in a csv file using R code and some package. For example, in the following code: links - cbind(rep('Click for Google',3),google search address goes here) ## R Mailing list blocks if I put the actual web address here write.table(links,'test.csv', sep=',',row.names=F,col.names=**F) the web address should be linked to 'Click for Google'. The browseURL() function open your internet browser with the url indicated as parameter: browseURL(http://www.r-**project.org http://www.r-project.org) But I am not sure how you want call it. You should be more precise about the context you want to use it. For example: links - data.frame(c(' for Google', ' for Bing'), c(http://www.google.com;, http://www.bing.com;), stringsAsFactors = FALSE) cat(Choose an option:\n, paste(1:2, links[,1],\n)) f-scan(nmax=1, quiet=TRUE) browseURL(links[f,2]) Sincerely Marc -- __** Marc Girondot, Pr Laboratoire Ecologie, Systématique et Evolution Equipe de Conservation des Populations et des Communautés CNRS, AgroParisTech et Université Paris-Sud 11 , UMR 8079 Bâtiment 362 91405 Orsay Cedex, France Tel: 33 1 (0)1.69.15.72.30 Fax: 33 1 (0)1.69.15.73.53 e-mail: marc.giron...@u-psud.fr Web: http://www.ese.u-psud.fr/epc/**conservation/Marc.htmlhttp://www.ese.u-ps ud.fr/epc/conservation/Marc.html Skype: girondot [[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] melt with complications
## Can someone suggest a simpler expression than either of these, with the goal ## of taking a long matrix into a wide one with exactly one of the factors converted to ## columns and all the rest retained as factors. I want something that generalizes beyond ## the three factors illustrated here. ## Rich meltTest - data.frame(A=rep(c(B,C), each=12), D=rep(c(E,F,G), each=4, times=2), H=rep(c(I,J,K,L), times=6), M=1:24) meltTest result.melt - do.call(rbind, { tmp - cast(D ~ H | A, value=M, data=meltTest) lapply(names(tmp), function(x) cbind(A=x, tmp[[x]])) ## explicit use of name A }) result.melt result.reshape - reshape(meltTest, direction=wide, timevar=H, idvar=c(A,D)) names(result.reshape)[3:6] - unique(as.character(meltTest$H)) ## explicit use of name H result.reshape [[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] try/tryCatch
here is the error: aa-metatrialstry(beta_5_50) Error in asMethod(object) : matrix is not symmetric [1,2] metatrials, the function that i am attempting to convert with try/tryCatch, gives me back a matrix with as many rows are there are simulations (z) in the aray with dim(x,y,z). with the data i attached, x is 500(number of patients), y is 9 (these are covariates), and z is 500. metatrials-function(mydata){ a-matrix(data=NA, nrow=dim(mydata)[3], ncol=5) colnames(a)=c(sens, spec, corr, sens_se, spec_se) for(ii in 1:dim(mydata)[3]){ tmp-mydata[,,ii] tmp1-as.data.frame(tmp) names(tmp1)=c(persons, d1, tp, fn, fp, fn, detect, d0, outcome) lm1-lmer(outcome~0+d1+d0+(0+d1+d0 | persons), family=binomial, data=tmp1, nAGQ=3) a[ii,1]=lm1@fixef[1] a[ii,2]=lm1@fixef[2] a[ii,3]=vcov(lm1)[1,2]/prod(sqrt(diag(vcov(lm1 a[ii,4:5]=sqrt(diag(vcov(lm1))) } return(a) } ## what i want is for the function to go on to the next data set in the array and simply return an NA for that line in the metatrials results. so basically, just keep going. thanks so much for your help! -Lisa On Mon, Mar 18, 2013 at 4:24 PM, jim holtman jholt...@gmail.com wrote: It would help if you told us what type of error you are getting and to also provide sample data so that we could run it to see what happens. I use 'try' a lot to catch errors and have not had any problems with it. On Mon, Mar 18, 2013 at 6:11 AM, lisa wang.li...@gmail.com wrote: Hi All, I have tried every fix on my try or tryCatch that I have found on the internet, but so far have not been able to get my R code to continue with the for loop after the lmer model results in an error. Here is two attemps of my code, the input is a 3D array file, but really any function would do metatrialstry-function(mydata){ a-matrix(data=NA, nrow=dim(mydata)[3], ncol=5) #colnames(a)=c(sens, spec, corr, sens_se, spec_se, counter)#colnames(a)=c(sens, spec, corr, sens_se, spec_se) k=1 for(ii in 1:dim(mydata)[3]){ tmp-mydata[,,ii] tmp1-as.data.frame(tmp) names(tmp1)=c(persons, d1, tp, fn, fp, fn, detect, d0, outcome) lm1-try(lmer(outcome~0+d1+d0+(0+d1+d0 | persons), family=binomial, data=tmp1, nAGQ=3), silent=T) if(class(lm1)[1]!='try-error'){ a[ii,1]=lm1@fixef[1] a[ii,2]=lm1@fixef[2] a[ii,3]=vcov(lm1)[1,2]/prod(sqrt(diag(vcov(lm1 a[ii,4:5]=sqrt(diag(vcov(lm1))) } } #k=k+1 #a[ii,6]=k return(a) } # # try / try catch ### # metatrialstry-function(mydata){ a-matrix(data=NA, nrow=dim(mydata)[3], ncol=5) #colnames(a)=c(sens, spec, corr, sens_se, spec_se, counter) colnames(a)=c(sens, spec, corr, sens_se, spec_se) #a[,6]=rep(0, length(a[,6])) for(ii in 1:dim(mydata)[3]){ tmp-mydata[,,ii] tmp1-as.data.frame(tmp) names(tmp1)=c(persons, d1, tp, fn, fp, fn, detect, d0, outcome) lm1-tryCatch(lmer(outcome~0+d1+d0+(0+d1+d0 | persons), family=binomial, data=tmp1, nAGQ=3), error=function(e) e) a[ii,1]=lm1@fixef[1] a[ii,2]=lm1@fixef[2] a[ii,3]=vcov(lm1)[1,2]/prod(sqrt(diag(vcov(lm1 a[ii,4:5]=sqrt(diag(vcov(lm1))) } return(a) } Any guidance would be greatly appreciated... thanks! Lisa -- Lisa Wang email: wang.li...@gmail.com cell: +49 -0176-87786557 Tübingen, Germany, 72070 [[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. -- 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. -- Lisa Wang email: wang.li...@gmail.com cell: +49 -0176-87786557 Tübingen, Germany, 72070 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Fit a mixture of lognormal and normal distributions
Hello I am trying to find an automated way of fitting a mixture of normal and log-normal distributions to data which is clearly bimodal. Here's a simulated example: x.1-rnorm(6000, 2.4, 0.6)x.2-rlnorm(1, 1.3,0.1)X-c(x.1, x.2) hist(X,100,freq=FALSE, ylim=c(0,1.5))lines(density(x.1), lty=2, lwd=2)lines(density(x.2), lty=2, lwd=2)lines(density(X), lty=4) Currently i am using mixtools and just run: library(mixtools)mixmdl = normalmixEM(X, k=2, epsilon = 1e-08, maxit = 1000, maxrestarts=20, verb = TRUE, fast=FALSE, ECM = FALSE, arbmean = TRUE, arbvar = TRUE) plot(mixmdl,which=2)lines(density(X), lty=2, lwd=2) This is obviously not the best way of doing this. The estimates it gives me are:mu 3.6595737(x.2 log()=1.29) 2.3113135(x.1) These are not too far off but I was wondering if someone knows of a better R package/way of doing this or has any hints that would help me coding it from scratch (EM+writing down the formulae for ML? but... would this really be better than using a more-advanced-already-available-package-made-by-pros?). The objective is to estimate threshold values at specific FDRs. (some help with that would also be most helpful.) Thanks to all in advance!To' [[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] data.frame with NA
I have this little data.frame http://dl.dropbox.com/u/102669/nanotna.rdata Two column contains NA, so the best thing to do is use na.locf function (with fromLast = T) But locf function doesn't work because NA in my data.frame are not recognized as real NA. Is there a way to substitute fake NA with real NA? In this case na.locf function should work Thank you __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Confirmatory factor analysis using the sem package. TLI CFI and RMSEA absent from model summary.
Hi Kevin, With sem package use : summary(X,fit.indices=c(RMSEA,...)) to get RMSEA or anothers fit indices. See the section ML.methods in sem.pdf Hervé Hervé Le 18/03/2013 16:00, Kevin Cheung a écrit : Hi R-help, I am using the sem package to run confirmatory factor analysis (cfa) on some questionnaire data collected from 307 participants. I have been running R-2.15.3 in Windows in conjunction with R studio. The model I am using was developed from exploratory factor analysis of a separate dataset (n=439); it includes 18 items that load onto 3 factors. I have used the sem package documentation and this video (http://vimeo.com/38941937) to run the cfa and obtain a chi-square statistic for the model. However, when I use the summary() function, the model does not provide indices of fit; I need these as part of my analysis output. In particular, I am looking for the Tucker Lewis Index (TLI), Comparative Fit Index (CFI), the Root Mean Square of Approximation (RMSEA). I have checked the documentation and cannot seem to find any reason for this; none of the arguments listed with the sem command indicate that I have to specify these as part of the output. In addition, the analysis example f! rom the video includes these indices as part of the output, but my analysis does not provide them. I have included my code with comments below: library(sem) validation.data - structure(list(V1 = c(5L, 4L, 2L, 4L, 5L, 6L, 6L, 4L, 5L, 3L, 6L, 5L, 4L, 5L, 3L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 4L, 5L, 4L, 4L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 2L, 6L, 5L, 6L, 4L, 5L, 6L, 5L, 5L, 4L, 5L, 5L, 3L, 5L, 5L, 5L, 5L, 4L, 2L, 5L, 5L, 5L, 4L, 6L, 4L, 6L, 5L, 5L, 5L, 4L, 5L, 5L, 4L, 5L, 6L, 4L, 5L, 4L, 5L, 5L, 5L, 3L, 5L, 5L, 3L, 5L, 4L, 5L, 2L, 6L, 4L, 4L, 4L, 5L, 5L, 4L, 6L, 2L, 4L, 5L, 4L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 2L, 4L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 6L, 5L, 5L, 4L, 4L, 5L, 4L, 4L, 5L, 4L, 4L, 4L, 5L, 5L, 6L, 5L, 5L, 5L, 4L, 4L, 3L, 4L, 5L, 5L, 5L, 2L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 6L, 4L, 5L, 5L, 6L, 5L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 4L, 4L, 4L, 5L, 2L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 4L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 4L, 4L, 2L, 6L, 5L, 5L, 5L, 6L, 4L, 3L, 5L, 5L, 5L, 5L, 4L, 4L, 6L, 6L, 5L, 5L, 5L, 3L, 6L, 5L, 5L, 5L, 4L, 5L, 5L, 4L, 5L, 4L, 6L, 5L, 4L, 4L, 5L, 4L, 3L, 4L, 5L, 4L, 5L, 6L, 2L, 4L, 4L, 5L, 4L, 4L, 4L, 5L, 6L, 5L, 4L, 4L, 4L, 6L, 6L, 4L, 5L, 5L, 5L, 2L, 4L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 3L, 5L, 6L, 5L, 4L, 4L, 3L, 3L, 4L, 5L, 5L, 1L, 4L, 5L, 3L, 5L, 1L, 6L, 5L, 4L, 4L, 5L, 5L, 4L, 5L, 5L, 6L, 5L, 5L, 5L), V2 = c(5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 5L, 6L, 6L, 4L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 4L, 6L, 5L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 5L, 4L, 6L, 5L, 6L, 5L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 5L, 4L, 6L, 4L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 4L, 2L, 5L, 6L, 4L, 5L, 5L, 6L, 6L, 5L, 6L, 4L, 6L, 5L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 6L, 5L, 6L, 5L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 4L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 3L, 5L, 6L, 6L, 6L, 5L, 6L, 5L, 5L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 5L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 4L, 5L, 6L, 6L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 5L, 6L, 5L, 5L, 6L, 5L, 5L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 6L, 6L), V3 = c(5L, 5L, 3L, 6L, 5L, 2L, 4L, 4L, 4L, 4L, 6L, 5L, 5L, 4L, 4L, 4L, 4L, 5L, 4L, 4L, 1L, 3L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 1L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 4L, 3L, 2L, 5L, 6L, 6L, 5L, 5L, 5L, 6L, 5L, 6L, 5L, 4L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 3L, 5L, 4L, 5L, 5L, 5L, 3L, 5L, 5L, 5L, 3L, 5L, 4L, 5L, 4L, 5L, 4L, 3L, 5L, 3L, 5L, 3L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 2L, 5L, 4L, 5L, 4L, 6L, 4L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 6L, 5L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 4L, 4L, 5L, 5L, 4L, 6L, 4L, 6L, 5L, 4L, 4L, 5L, 5L, 5L, 4L, 4L, 1L, 5L, 4L, 4L, 5L, 5L, 4L, 6L, 3L, 4L, 4L, 5L, 4L, 4L, 5L, 3L, 5L, 6L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 3L, 5L, 4L, 6L, 5L, 4L, 6L, 3L, 4L, 2L, 4L, 4L, 5L, 4L, 4L, 5L, 3L, 4L, 5L, 5L, 4L, 3L, 3L, 5L, 5L, 5L, 5L, 4L, 3L, 4L, 2L, 5L, 5L, 6L, 4L, 5L, 4L, 4L, 5L, 4L, 6L, 6L, 4L, 4L, 4L, 5L, 5L, 6L, 5L, 4L, 6L, 5L, 5L, 5L, 4L, 6L, 6L, 3L, 2L, 3L, 6L, 4L, 5L, 3L, 6L, 3L, 4L, 4L, 5L, 4L, 6L, 4L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 4L, 5L, 5L, 6L, 5L, 4L,
[R] Why stacking rasters return NAs?
I have several rasters that I want to do some calculations ,basically calculating the moving average. dir2 - list.files(D:\\2010+2011, *.bin, full.names = TRUE) saf=stack(dir2) movi - overlay(stack(saf),fun=function(x) movingFun(x, fun=mean, n=3, na.rm=TRUE)) Error in .overlayList(x, fun = fun, filename = filename, ...) : cannot use this formula, probably because it is not vectorized I then checked the data but found that all values were returnd as NA and this may explain why i am getting the error. saf class : RasterStack dimensions : 720, 1440, 1036800, 601 (nrow, ncol, ncell, nlayers) resolution : 0.25, 0.25 (x, y) extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax) coord. ref. : +proj=longlat +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +no_defs names : Vol_025_H//00_1_wgs84, Vol_025_H//00_1_wgs84, Vol_025_H//00_1_wgs84,Vol_025_H//00_1_wgs84, Vol_025_H//00_1_wgs84, Vol_025_H//00_1_wgs84, , ... min values :NA,NA, NA, NA,NA,NA, NA, NA,NA,NA, NA,NA,NA,NA, NA, ... max values :NA,NA, NA,NA,NA,NA, NA,NA,NA,NA, NA, NA,NA,NA, NA, ... I wonder why this is happening, I checked the files separably(summary) and everything was right!as you can see bellow: ol_025_H14_2011092000_1_wgs84 Vol_025_H14_2011092100_1_wgs84 Vol_025_H14_2011092200_1_wgs84 Vol_025_H14_2011092300_1_wgs84 Vol_025_H14_2011092400_1_wgs84 Min. 0.0 0.000 0.000 0.000 0.000 1st Qu.0.31883 0.3163167 0.3146436 0.3113111 0.3064551 Median .0 .000 .000 .000 .000 3rd Qu. .0 .000 .000 .000 .000 Max..0 .000 .000 .000 .000 NA's 0.0 0.000 0.000 0.0 I am gratful to anyhelp -- View this message in context: http://r.789695.n4.nabble.com/Why-stacking-rasters-return-NAs-tp4661706.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] Loop or some other way to parse by data generated values when it is not linear
I'm sorry for the really vague subject line but I am not sure how to succinctly describe what I am doing and what the problem is. But, here goes: 1. I have data with two-way data with frequencies. Below is an example, though in reality I am looking at about 10 different variables that I am crossing so the values of X1 and X2 change. X1 and X2 are place holders. Here's the dataset (though using this first part does not happen in reality): X1 - matrix(c(0, 1, 2, 3, 4, 99), nrow=18, ncol=1, byrow=T) X2 - sort(matrix(c(0, 2, 4), nrow=18, ncol=1, byrow=T), decreasing=F) Y - matrix(c(83, 107, 47, 27, 38, 1, 12 ,25, 14, 4, 9, 0, 14, 27, 28, 13, 18, 0), nrow=18, ncol=1, byrow=T) tmp.n - data.frame(X1, X2, Y) The final data frame is what I actually get: X1 X2 Y 1 0 0 83 2 1 0 107 3 2 0 47 4 3 0 27 5 4 0 38 6 99 0 1 7 0 2 12 8 1 2 25 9 2 2 14 10 3 2 4 11 4 2 9 12 99 2 0 13 0 4 14 14 1 4 27 15 2 4 28 16 3 4 13 17 4 4 18 18 99 4 0 2. What I want is: 0 2 4 0 83 12 14 1 107 25 27 2 47 14 28 3 27 4 13 4 38 9 18 99 1 0 0 3. I've been trying to do it using this (which is inside a function so I can vary what variables X1 and X2 are): X1 - table(tmp.n[,1]) X2 - table(tmp.n[,2]) # Create the tmp.n.# datasets that contain the Y's. Do this in a loop to automate dta - NULL for (i in 0:length(X1)) { assign(tmp.n_, tmp.n[tmp.n[,1] == i, c(1,3)]) tmp.n_ - data.frame(tmp.n_[,2]) dta[i] - assign(paste(tmp.n., i, sep=), tmp.n_) dta } dta2 - (data.frame(matrix(unlist(dta), nrow=n2[1], byrow=T))) colnames(dta2) - names(X2) dta2 And that works so long as X1 and X2 are linear. In other words, if X1 - seq(0, 4, 1). But that 99 throws the whole thing off and it gives me this: X1 X2 1 107 25 2 27 47 3 14 28 4 27 4 5 13 38 6 9 18 It's basically breaks the whole thing. I've not been able to figure this out and I've been like a dog with a bone trying to make it work with modifications to the for loop. I know there is an easier way to do this, but my brain is no longer capable of thinking outside the box I've put it in. So, I am turning to you for help. Best, Jen __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] data.frame with NA
On 18-03-2013, at 16:49, Pete freeri...@gmail.com wrote: I have this little data.frame http://dl.dropbox.com/u/102669/nanotna.rdata Two column contains NA, so the best thing to do is use na.locf function (with fromLast = T) But locf function doesn't work because NA in my data.frame are not recognized as real NA. Is there a way to substitute fake NA with real NA? In this case na.locf function should work Your data are all characters. Do str(db) to see that. What is probably supposed to be numeric is also character, Somehow you have managed to read in data that R thinks is all chr. Your NA are NA in reality: a character string NA. You will have to review the method you used to get the data into R. And make sure that what you want to be numeric is indeed numeric. Then you can start to think about doing something about the NA's. 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] Loop or some other way to parse by data generated values when it is not linear
Hi, library(reshape2) dcast(tmp.n,X1~X2,value.var=Y) # X1 0 2 4 #1 0 83 12 14 #2 1 107 25 27 #3 2 47 14 28 #4 3 27 4 13 #5 4 38 9 18 #6 99 1 0 0 A.K. - Original Message - From: plessthanpointohf...@gmail.com plessthanpointohf...@gmail.com To: r-help@r-project.org Cc: Sent: Monday, March 18, 2013 1:44 PM Subject: [R] Loop or some other way to parse by data generated values when it is not linear I'm sorry for the really vague subject line but I am not sure how to succinctly describe what I am doing and what the problem is. But, here goes: 1. I have data with two-way data with frequencies. Below is an example, though in reality I am looking at about 10 different variables that I am crossing so the values of X1 and X2 change. X1 and X2 are place holders. Here's the dataset (though using this first part does not happen in reality): X1 - matrix(c(0, 1, 2, 3, 4, 99), nrow=18, ncol=1, byrow=T) X2 - sort(matrix(c(0, 2, 4), nrow=18, ncol=1, byrow=T), decreasing=F) Y - matrix(c(83, 107, 47, 27, 38, 1, 12 ,25, 14, 4, 9, 0, 14, 27, 28, 13, 18, 0), nrow=18, ncol=1, byrow=T) tmp.n - data.frame(X1, X2, Y) The final data frame is what I actually get: X1 X2 Y 1 0 0 83 2 1 0 107 3 2 0 47 4 3 0 27 5 4 0 38 6 99 0 1 7 0 2 12 8 1 2 25 9 2 2 14 10 3 2 4 11 4 2 9 12 99 2 0 13 0 4 14 14 1 4 27 15 2 4 28 16 3 4 13 17 4 4 18 18 99 4 0 2. What I want is: 0 2 4 0 83 12 14 1 107 25 27 2 47 14 28 3 27 4 13 4 38 9 18 99 1 0 0 3. I've been trying to do it using this (which is inside a function so I can vary what variables X1 and X2 are): X1 - table(tmp.n[,1]) X2 - table(tmp.n[,2]) # Create the tmp.n.# datasets that contain the Y's. Do this in a loop to automate dta - NULL for (i in 0:length(X1)) { assign(tmp.n_, tmp.n[tmp.n[,1] == i, c(1,3)]) tmp.n_ - data.frame(tmp.n_[,2]) dta[i] - assign(paste(tmp.n., i, sep=), tmp.n_) dta } dta2 - (data.frame(matrix(unlist(dta), nrow=n2[1], byrow=T))) colnames(dta2) - names(X2) dta2 And that works so long as X1 and X2 are linear. In other words, if X1 - seq(0, 4, 1). But that 99 throws the whole thing off and it gives me this: X1 X2 1 107 25 2 27 47 3 14 28 4 27 4 5 13 38 6 9 18 It's basically breaks the whole thing. I've not been able to figure this out and I've been like a dog with a bone trying to make it work with modifications to the for loop. I know there is an easier way to do this, but my brain is no longer capable of thinking outside the box I've put it in. So, I am turning to you for help. Best, Jen __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Loop or some other way to parse by data generated values when it is not linear
On Mar 18, 2013, at 12:44 PM, plessthanpointohf...@gmail.com wrote: I'm sorry for the really vague subject line but I am not sure how to succinctly describe what I am doing and what the problem is. But, here goes: 1. I have data with two-way data with frequencies. Below is an example, though in reality I am looking at about 10 different variables that I am crossing so the values of X1 and X2 change. X1 and X2 are place holders. Here's the dataset (though using this first part does not happen in reality): X1 - matrix(c(0, 1, 2, 3, 4, 99), nrow=18, ncol=1, byrow=T) X2 - sort(matrix(c(0, 2, 4), nrow=18, ncol=1, byrow=T), decreasing=F) Y - matrix(c(83, 107, 47, 27, 38, 1, 12 ,25, 14, 4, 9, 0, 14, 27, 28, 13, 18, 0), nrow=18, ncol=1, byrow=T) tmp.n - data.frame(X1, X2, Y) The final data frame is what I actually get: X1 X2 Y 1 0 0 83 2 1 0 107 3 2 0 47 4 3 0 27 5 4 0 38 6 99 0 1 7 0 2 12 8 1 2 25 9 2 2 14 10 3 2 4 11 4 2 9 12 99 2 0 13 0 4 14 14 1 4 27 15 2 4 28 16 3 4 13 17 4 4 18 18 99 4 0 2. What I want is: 0 2 4 0 83 12 14 1 107 25 27 2 47 14 28 3 27 4 13 4 38 9 18 99 1 0 0 3. I've been trying to do it using this (which is inside a function so I can vary what variables X1 and X2 are): X1 - table(tmp.n[,1]) X2 - table(tmp.n[,2]) # Create the tmp.n.# datasets that contain the Y's. Do this in a loop to automate dta - NULL for (i in 0:length(X1)) { assign(tmp.n_, tmp.n[tmp.n[,1] == i, c(1,3)]) tmp.n_ - data.frame(tmp.n_[,2]) dta[i] - assign(paste(tmp.n., i, sep=), tmp.n_) dta } dta2 - (data.frame(matrix(unlist(dta), nrow=n2[1], byrow=T))) colnames(dta2) - names(X2) dta2 And that works so long as X1 and X2 are linear. In other words, if X1 - seq(0, 4, 1). But that 99 throws the whole thing off and it gives me this: X1 X2 1 107 25 2 27 47 3 14 28 4 27 4 5 13 38 6 9 18 It's basically breaks the whole thing. I've not been able to figure this out and I've been like a dog with a bone trying to make it work with modifications to the for loop. I know there is an easier way to do this, but my brain is no longer capable of thinking outside the box I've put it in. So, I am turning to you for help. Best, Jen Something like this? xtabs(Y ~ X1 + X2, data = tmp.n) X2 X1 0 2 4 0 83 12 14 1 107 25 27 2 47 14 28 3 27 4 13 4 38 9 18 99 1 0 0 See ?xtabs Regards, Marc Schwartz __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Confirmatory factor analysis using the sem package. TLI CFI and RMSEA absent from model summary.
Dear Kevin, In earlier versions of the sem package, the summary() method for model objects produced many (and from version to version, an increasing number of) fit indices by default. For various reasons, I decided to provide only the LR chisquare test for the model, the BIC, and the AIC by default, and to leave it up the user to decide what fit indices to show. This can be controlled with the fit.indices option as well as the fit.indices argument to the summary() method. Best, John On Mon, 18 Mar 2013 16:46:03 + Kevin Cheung k.che...@derby.ac.uk wrote: Dear John, Thank you for taking the time to help me with this, I have been able to obtain fit indices using the information that you provided. Note to users searching archived R-help posts about this issue: The instructional video I looked at (http://vimeo.com/38941937) gave fit indices using the default summary() command without any additional arguments. This may have been due to a different version of R (I noticed that the instructor was using a mac based OS). With regards, Kevin Kevin Yet Fong Cheung, Bsc., MRes., MBPsS. Postgraduate Researcher Centre for Psychological Research University of Derby Kedleston Road Derby DE22 1GB k.che...@derby.ac.ukmailto:k.che...@derby.ac.uk 01332592081 http://derby.academia.edu/KevinCheung -Original Message- From: John Fox [mailto:j...@mcmaster.ca] Sent: 18 March 2013 15:55 To: Kevin Cheung Cc: r-help@r-project.org Subject: Re: [R] Confirmatory factor analysis using the sem package. TLI CFI and RMSEA absent from model summary. Dear Kevin, See ?summary.objectiveML, and in particular the description of the fit.indices argument. By default, the summary() method doesn't print many fit indices, but many are available optionally. I hope this helps, John John Fox Sen. William McMaster Prof. of Social Statistics Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox/ On Mon, 18 Mar 2013 15:00:06 + Kevin Cheung k.che...@derby.ac.uk wrote: Hi R-help, I am using the sem package to run confirmatory factor analysis (cfa) on some questionnaire data collected from 307 participants. I have been running R-2.15.3 in Windows in conjunction with R studio. The model I am using was developed from exploratory factor analysis of a separate dataset (n=439); it includes 18 items that load onto 3 factors. I have used the sem package documentation and this video (http://vimeo.com/38941937) to run the cfa and obtain a chi-square statistic for the model. However, when I use the summary() function, the model does not provide indices of fit; I need these as part of my analysis output. In particular, I am looking for the Tucker Lewis Index (TLI), Comparative Fit Index (CFI), the Root Mean Square of Approximation (RMSEA). I have checked the documentation and cannot seem to find any reason for this; none of the arguments listed with the sem command indicate that I have to specify these as part of the output. In addition, the analysis examp! le! f! rom the video includes these indices as part of the output, but my analysis does not provide them. I have included my code with comments below: library(sem) validation.data - structure(list(V1 = c(5L, 4L, 2L, 4L, 5L, 6L, 6L, 4L, 5L, 3L, 6L, 5L, 4L, 5L, 3L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 4L, 5L, 4L, 4L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 2L, 6L, 5L, 6L, 4L, 5L, 6L, 5L, 5L, 4L, 5L, 5L, 3L, 5L, 5L, 5L, 5L, 4L, 2L, 5L, 5L, 5L, 4L, 6L, 4L, 6L, 5L, 5L, 5L, 4L, 5L, 5L, 4L, 5L, 6L, 4L, 5L, 4L, 5L, 5L, 5L, 3L, 5L, 5L, 3L, 5L, 4L, 5L, 2L, 6L, 4L, 4L, 4L, 5L, 5L, 4L, 6L, 2L, 4L, 5L, 4L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 2L, 4L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 5L, 6L, 5L, 5L, 4L, 4L, 5L, 4L, 4L, 5L, 4L, 4L, 4L, 5L, 5L, 6L, 5L, 5L, 5L, 4L, 4L, 3L, 4L, 5L, 5L, 5L, 2L, 5L, 5L, 5L, 4L, 5L, 5L, 5L, 5L, 6L, 4L, 5L, 5L, 6L, 5L, 5L, 5L, 5L, 4L, 5L, 4L, 5L, 5L, 4L, 4L, 4L, 4L, 5L, 2L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 4L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 5L, 4L, 4L, 5L, 5L, 5L, 5L, 4L, 4L, 2L, 6L, 5L, 5L, 5L, 6L, 4L, 3L, 5L, 5L, 5L, 5L, 4L, 4L, 6L, 6L, 5L, 5L, 5L, 3L, 6L, 5L, 5L, 5L, 4L, 5L, 5L, 4L, 5L, 4L, 6L, 5L, 4L, 4L, 5L, 4L, 3L, 4L, 5L, 4L, 5L, 6L, 2L, 4L, 4L, 5L, 4L, 4L, 4L, 5L, 6L, 5L, 4L, 4L, 4L, 6L, 6L, 4L, 5L, 5L, 5L, 2L, 4L, 4L, 5L, 4L, 5L, 5L, 4L, 5L, 3L, 5L, 6L, 5L, 4L, 4L, 3L, 3L, 4L, 5L, 5L, 1L, 4L, 5L, 3L, 5L, 1L, 6L, 5L, 4L, 4L, 5L, 5L, 4L, 5L, 5L, 6L, 5L, 5L, 5L), V2 = c(5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 6L, 5L, 6L, 6L, 4L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 6L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 5L, 5L, 6L, 5L, 6L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 4L, 6L, 5L, 6L, 5L, 5L,
[R] Lattice strip with custom shingle.interval
Trying to make a strip in lattice xyplots to display shingle intervals with a sensible number of digits: First I make a 2x2 matrix with rounded values of the shingle itnervals: lec - round(matrix(unlist(levels(equal.count(sd$Frother))), ncol = 2, byrow = TRUE), 3) Next I make a custom strip function to handle the difference between a factor proper and a shingle: my.strip - function(which.given, which.panel, ..., shingle.intervals) { if (which.given == 1) { strip.default(which.given, which.panel, ...,strip.name = FALSE, strip.levels = TRUE, shingle.intervals = lec) } else { strip.default(which.given, which.panel, ..., strip.name = FALSE, strip.levels = TRUE) } } Then the call to xyplot. xyplot(Grade ~ Recovery| equal.count(Frother) + Row, groups = interaction(Row1.Mode, Row3.Mode, drop = TRUE), data = sd, pch = 1, layout = c(6,2), main = 'By Row and Frother dosing', strip = my.strip, sub = 'Colour indicates control mode', auto.key = list(columns = 5)) This however still gives full precision levels. Any help much appreciated. Alex van der Spek __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Loop or some other way to parse by data generated values when it is not linear
That works great! I didn't even know there was a reshape2 package. Thanks! Jen On 03/18/2013 01:49 PM, arun wrote: Hi, library(reshape2) dcast(tmp.n,X1~X2,value.var=Y) # X1 0 2 4 #1 0 83 12 14 #2 1 107 25 27 #3 2 47 14 28 #4 3 27 4 13 #5 4 38 9 18 #6 99 1 0 0 A.K. - Original Message - From: plessthanpointohf...@gmail.com plessthanpointohf...@gmail.com To: r-help@r-project.org Cc: Sent: Monday, March 18, 2013 1:44 PM Subject: [R] Loop or some other way to parse by data generated values when it is not linear I'm sorry for the really vague subject line but I am not sure how to succinctly describe what I am doing and what the problem is. But, here goes: 1. I have data with two-way data with frequencies. Below is an example, though in reality I am looking at about 10 different variables that I am crossing so the values of X1 and X2 change. X1 and X2 are place holders. Here's the dataset (though using this first part does not happen in reality): X1 - matrix(c(0, 1, 2, 3, 4, 99), nrow=18, ncol=1, byrow=T) X2 - sort(matrix(c(0, 2, 4), nrow=18, ncol=1, byrow=T), decreasing=F) Y - matrix(c(83, 107, 47, 27, 38, 1, 12 ,25, 14, 4, 9, 0, 14, 27, 28, 13, 18, 0), nrow=18, ncol=1, byrow=T) tmp.n - data.frame(X1, X2, Y) The final data frame is what I actually get: X1 X2 Y 1 0 0 83 2 1 0 107 3 2 0 47 4 3 0 27 5 4 0 38 6 99 0 1 7 0 2 12 8 1 2 25 9 2 2 14 10 3 2 4 11 4 2 9 12 99 2 0 13 0 4 14 14 1 4 27 15 2 4 28 16 3 4 13 17 4 4 18 18 99 4 0 2. What I want is: 0 2 4 0 83 12 14 1 107 25 27 2 47 14 28 3 27 4 13 4 38 9 18 99 1 0 0 3. I've been trying to do it using this (which is inside a function so I can vary what variables X1 and X2 are): X1 - table(tmp.n[,1]) X2 - table(tmp.n[,2]) # Create the tmp.n.# datasets that contain the Y's. Do this in a loop to automate dta - NULL for (i in 0:length(X1)) { assign(tmp.n_, tmp.n[tmp.n[,1] == i, c(1,3)]) tmp.n_ - data.frame(tmp.n_[,2]) dta[i] - assign(paste(tmp.n., i, sep=), tmp.n_) dta } dta2 - (data.frame(matrix(unlist(dta), nrow=n2[1], byrow=T))) colnames(dta2) - names(X2) dta2 And that works so long as X1 and X2 are linear. In other words, if X1 - seq(0, 4, 1). But that 99 throws the whole thing off and it gives me this: X1 X2 1 107 25 2 27 47 3 14 28 4 27 4 5 13 38 6 9 18 It's basically breaks the whole thing. I've not been able to figure this out and I've been like a dog with a bone trying to make it work with modifications to the for loop. I know there is an easier way to do this, but my brain is no longer capable of thinking outside the box I've put it in. So, I am turning to you for help. Best, Jen __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] OrgMassSpecR peak area issue
Hello! I am having an issue with the OrgMassSpecR package. I run my HPLC using a DAD detector. My raw data is exported form chemstation as a csv file. I then upload the csv into Rstudio no problem. Using the DrawChromatogram function, I get a nice chromatogram, and my retention time, peak area, and apex intensity values are given as well. The problem comes with the peak area value given. The peak area is much smaller than a value that would make sense. My peak area value is actually less than my apex intensity value. Is this because I am using a DAD detector rather than an MS? If so, is there a simply way to edit the peak area equation so that it will also work with absorbance values? Any help is greatly appreciated. Thanks for your time. Chris Beaver [[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] Loop or some other way to parse by data generated values when it is not linear
This works really great, too. Thanks for this option. I am glad I have two ways to accomplish this. I KNEW it was going to be that simple. I was convinced I was making more out of it than I needed to and I'm glad to see I was right. Best, Jen (It's fun to be new at something when I'm this old). On 03/18/2013 01:50 PM, Marc Schwartz wrote: On Mar 18, 2013, at 12:44 PM, plessthanpointohf...@gmail.com wrote: I'm sorry for the really vague subject line but I am not sure how to succinctly describe what I am doing and what the problem is. But, here goes: 1. I have data with two-way data with frequencies. Below is an example, though in reality I am looking at about 10 different variables that I am crossing so the values of X1 and X2 change. X1 and X2 are place holders. Here's the dataset (though using this first part does not happen in reality): X1 - matrix(c(0, 1, 2, 3, 4, 99), nrow=18, ncol=1, byrow=T) X2 - sort(matrix(c(0, 2, 4), nrow=18, ncol=1, byrow=T), decreasing=F) Y - matrix(c(83, 107, 47, 27, 38, 1, 12 ,25, 14, 4, 9, 0, 14, 27, 28, 13, 18, 0), nrow=18, ncol=1, byrow=T) tmp.n - data.frame(X1, X2, Y) The final data frame is what I actually get: X1 X2 Y 1 0 0 83 2 1 0 107 3 2 0 47 4 3 0 27 5 4 0 38 6 99 0 1 7 0 2 12 8 1 2 25 9 2 2 14 10 3 2 4 11 4 2 9 12 99 2 0 13 0 4 14 14 1 4 27 15 2 4 28 16 3 4 13 17 4 4 18 18 99 4 0 2. What I want is: 0 2 4 0 83 12 14 1 107 25 27 2 47 14 28 3 27 4 13 4 38 9 18 99 1 0 0 3. I've been trying to do it using this (which is inside a function so I can vary what variables X1 and X2 are): X1 - table(tmp.n[,1]) X2 - table(tmp.n[,2]) # Create the tmp.n.# datasets that contain the Y's. Do this in a loop to automate dta - NULL for (i in 0:length(X1)) { assign(tmp.n_, tmp.n[tmp.n[,1] == i, c(1,3)]) tmp.n_ - data.frame(tmp.n_[,2]) dta[i] - assign(paste(tmp.n., i, sep=), tmp.n_) dta } dta2 - (data.frame(matrix(unlist(dta), nrow=n2[1], byrow=T))) colnames(dta2) - names(X2) dta2 And that works so long as X1 and X2 are linear. In other words, if X1 - seq(0, 4, 1). But that 99 throws the whole thing off and it gives me this: X1 X2 1 107 25 2 27 47 3 14 28 4 27 4 5 13 38 6 9 18 It's basically breaks the whole thing. I've not been able to figure this out and I've been like a dog with a bone trying to make it work with modifications to the for loop. I know there is an easier way to do this, but my brain is no longer capable of thinking outside the box I've put it in. So, I am turning to you for help. Best, Jen Something like this? xtabs(Y ~ X1 + X2, data = tmp.n) X2 X1 0 2 4 0 83 12 14 1 107 25 27 2 47 14 28 3 27 4 13 4 38 9 18 99 1 0 0 See ?xtabs Regards, Marc Schwartz __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 stop set.seed() besides exiting out of R?
Hi list, I am curious how to stop the set.seed(), I don't want the same repeated random number. I know I can set it to a different seed, but I don't want to go through the trouble of setting different seed every time. Thanks, Mike [[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] new question
z.boxplot- function(lst){ new.list- lapply(lst,function(x) x[x$FDR0.01,]) print(new.list) par(mfrow=c(2,2)) b1-lapply(names(new.list),function(x) lapply(new.list[x],function(y) boxplot(FDR~z,data=y,xlab=Charge,ylab=FDR,main=x))) } z.boxplot(ListFacGroup) #prints new.list If you want to turn off that: z.boxplot- function(lst){ new.list- lapply(lst,function(x) x[x$FDR0.01,]) #print(new.list) par(mfrow=c(2,2)) b1-lapply(names(new.list),function(x) lapply(new.list[x],function(y) boxplot(FDR~z,data=y,xlab=Charge,ylab=FDR,main=x))) } z.boxplot(ListFacGroup) A.K. From: Vera Costa veracosta...@gmail.com To: arun smartpink...@yahoo.com Sent: Monday, March 18, 2013 1:59 PM Subject: Re: new question For example, if I run you code without pdf and dev.off I have what I want directory- C:/Users/Vera Costa/Desktop/dados.lixo #modified the function GetFileList - function(directory,number){ setwd(directory) filelist1-dir() lista-dir(directory,pattern = paste(MSMS_,number,PepInfo.txt,sep=), full.names = TRUE, recursive = TRUE) output- list(filelist1,lista) return(output) } file.list.names-GetFileList(directory,23)[[1]] lista-GetFileList(directory,23)[[2]] FacGroup-c(0,1,0,2,2,0,3) ReadDir-function(FacGroup){ list.new-lista[FacGroup!=0] read.list-lapply(list.new, function(x) read.table(x,header=TRUE, sep = \t)) names(read.list)-file.list.names[FacGroup!=0] return (read.list) } ListFacGroup-ReadDir(FacGroup) ListFacGroup z.boxplot- function(lst){ new.list- lapply(lst,function(x) x[x$FDR0.01,]) print(new.list) #pdf(VeraBP.pdf) par(mfrow=c(2,2)) lapply(names(new.list),function(x) lapply(new.list[x],function(y) boxplot(FDR~z,data=y,xlab=Charge,ylab=FDR,main=x))) #dev.off() } z.boxplot(ListFacGroup) But I have the results too (I don't need it) [[1]] [[1]]$a2 [[1]]$a2$stats [,1] [,2] [1,] 0.00 0.00 [2,] 0.00 0.00 [3,] 0.0001355197 0.0002175095 [4,] 0.0010588777 0.0004350190 [5,] 0.0016571381 0.0004350190 [[1]]$a2$n [1] 8 2 [[1]]$a2$conf [,1] [,2] [1,] -0.0004559846 -0.0002685062 [2,] 0.0007270240 0.0007035253 [[1]]$a2$out [1] 0.00494012 [[1]]$a2$group [1] 1 [[1]]$a2$names [1] 2 3 [[2]] [[2]]$c2 [[2]]$c2$stats [,1] [,2] [1,] 0.00 0.00 [2,] 0.00 0.00 [3,] 0.0001355197 0.0002175095 [4,] 0.0010588777 0.0004350190 [5,] 0.0016571381 0.0004350190 [[2]]$c2$n [1] 8 2 [[2]]$c2$conf [,1] [,2] [1,] -0.0004559846 -0.0002685062 [2,] 0.0007270240 0.0007035253 [[2]]$c2$out [1] 0.00494012 [[2]]$c2$group [1] 1 [[2]]$c2$names [1] 2 3 [[3]] [[3]]$c3 [[3]]$c3$stats [,1] [,2] [1,] 0.0 0.00 [2,] 0.0 0.00 [3,] 0.0 0.00 [4,] 0.002226409 0.0002086594 [5,] 0.002226409 0.0004173187 [[3]]$c3$n [1] 6 4 [[3]]$c3$conf [,1] [,2] [1,] -0.001436105 -0.0001648409 [2,] 0.001436105 0.0001648409 [[3]]$c3$out [1] 0.00560348 [[3]]$c3$group [1] 1 [[3]]$c3$names [1] 2 3 [[4]] [[4]]$t2 [[4]]$t2$stats [,1] [,2] [,3] [1,] 0 0.00 0 [2,] 0 0.00 0 [3,] 0 0.0002908668 0 [4,] 0 0.0025929827 0 [5,] 0 0.005251 0 [[4]]$t2$n [1] 1 10 5 [[4]]$t2$conf [,1] [,2] [,3] [1,] 0 -0.001004691 0 [2,] 0 0.001586424 0 [[4]]$t2$out [1] 0.0092051934 0.0007174888 [[4]]$t2$group [1] 2 3 [[4]]$t2$names [1] 1 2 3 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 stop set.seed() besides exiting out of R?
-Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of C W Sent: Monday, March 18, 2013 11:27 AM To: r-help Subject: [R] How to stop set.seed() besides exiting out of R? Hi list, I am curious how to stop the set.seed(), I don't want the same repeated random number. I know I can set it to a different seed, but I don't want to go through the trouble of setting different seed every time. Thanks, Mike Can you show us how you are using set.seed() that results in getting the same sequence repeatedly? If you are doing simulations in a loop, then set the seed once, outside the loop. Otherwise, I am not sure what you are doing that causes problems. A reproducible example would really help. 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] How to stop set.seed() besides exiting out of R?
set.seed(100) for (i in 1:100){ a - rnorm(1000, mean=0, sd=1) hist(a) } #Now say, I want to simulate without being confined to set.seed(100), I just want to get a simulation like how R normally does it. Mike On Mon, Mar 18, 2013 at 2:40 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 C W Sent: Monday, March 18, 2013 11:27 AM To: r-help Subject: [R] How to stop set.seed() besides exiting out of R? Hi list, I am curious how to stop the set.seed(), I don't want the same repeated random number. I know I can set it to a different seed, but I don't want to go through the trouble of setting different seed every time. Thanks, Mike Can you show us how you are using set.seed() that results in getting the same sequence repeatedly? If you are doing simulations in a loop, then set the seed once, outside the loop. Otherwise, I am not sure what you are doing that causes problems. A reproducible example would really help. 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. [[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] OrgMassSpecR peak area issue
If you type DrawChromatogram you can see the method used to calculate the peak area. Looks to me like you could easily hack it if you wanted. The relevant part about peak areas is this: for (j in 1:n) { k - (j%%n) + 1 x[j] - peakTime[j] * peakIntensity[k] - peakTime[k] * peakIntensity[j] } peakArea[i] - abs(sum(x)/2) which looks pretty standard to me, though I'm not clear right off the top of my head why they are dividing by 2. You can always contact the maintainer. Bryan On Mar 18, 2013, at 1:34 PM, Christopher Beaver christopher.bea...@gmail.com wrote: Hello! I am having an issue with the OrgMassSpecR package. I run my HPLC using a DAD detector. My raw data is exported form chemstation as a csv file. I then upload the csv into Rstudio no problem. Using the DrawChromatogram function, I get a nice chromatogram, and my retention time, peak area, and apex intensity values are given as well. The problem comes with the peak area value given. The peak area is much smaller than a value that would make sense. My peak area value is actually less than my apex intensity value. Is this because I am using a DAD detector rather than an MS? If so, is there a simply way to edit the peak area equation so that it will also work with absorbance values? Any help is greatly appreciated. Thanks for your time. Chris Beaver [[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] OrgMassSpecR peak area issue
Hi Chris, I am having an issue with the OrgMassSpecR package. I run my HPLC using a DAD detector. You are on a statistics IDE related mailing list. Have mercy with people from other fields and tell them that you are using a diode array to measure UV/VIS absorption. (And possibly let them know that you expect the absorbance A = lg I_0 - lg I ~ c.) My raw data is exported form chemstation as a csv file. I then upload the csv into Rstudio no problem. Using the DrawChromatogram function, I get a nice chromatogram, and my retention time, peak area, and apex intensity values are given as well. The problem comes with the peak area value given. The peak area is much smaller than a value that would make sense. How do you know that (see next comment)? My peak area value is actually less than my apex intensity value. This is not a good criterion to determine what area value would actually make sense: area and intensity have different units! Possible solution: a glance on the code in DrawChromatogram reveals that really the polygon area is calculated (as the manual specifies). Thus the area will be in counts*s or counts*min, and of course 1 count*min = 60 counts*s. How long does your analyte take to elute? Unless it is 2 min (if time is in min) or 2 s (for time scale in s), the numeric value of the area should be A_max (approximating the peak as triangule). Your apex (max) absorbance should ideally be a bit below 1, so a rough guesstimate for peak area would be 1/2 A_max * Δt which will be quite below 1 if you measure time in minutes. If you detect by mass spec, you get ion counts which are large numbers, so areas are likely to be 1 (regardless of min or s time scale). Is this because I am using a DAD detector rather than an MS? If so, is there a simply way to edit the peak area equation so that it will also work with absorbance values? Most probably you just want to get your units right! Hope that helps, Claudia PS: for future questions of this sort, you may want to consider asking on stackoverflow.com (or chemistry.stackexchange.com) where you can post nicely formatted code, calculation results and images with your question. -- Claudia Beleites Spectroscopy/Imaging Institute of Photonic Technology Albert-Einstein-Str. 9 07745 Jena Germany email: claudia.belei...@ipht-jena.de phone: +49 3641 206-133 fax: +49 2641 206-399 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 stop set.seed() besides exiting out of R?
As I understand it, how R âânormallyâ does it is to use the system clock to set the seed once per session, unless you use set.seed() to set a new seed. You chose to set the seed to a different value. But from that point on, the pseudo random number generation continues in the same way it ânormallyâ does. In your code below, each of your 100 histograms will be different. If you then execute the for loop again (but not the set.seed(100) statement), you will get a different set of histograms. The only way you would be âconfined to set.seed(100)â is if you keep resetting the seed to 100. 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 From: C W [mailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 11:50 AM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? set.seed(100) for (i in 1:100){ a - rnorm(1000, mean=0, sd=1) hist(a) } #Now say, I want to simulate without being confined to set.seed(100), I just want to get a simulation like how R normally does it. Mike On Mon, Mar 18, 2013 at 2:40 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.gov wrote: -Original Message- From: r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.org [mailto:r-help-bounces@r-mailto:r-help-bounces@r- project.orghttp://project.org] On Behalf Of C W Sent: Monday, March 18, 2013 11:27 AM To: r-help Subject: [R] How to stop set.seed() besides exiting out of R? Hi list, I am curious how to stop the set.seed(), I don't want the same repeated random number. I know I can set it to a different seed, but I don't want to go through the trouble of setting different seed every time. Thanks, Mike Can you show us how you are using set.seed() that results in getting the same sequence repeatedly? If you are doing simulations in a loop, then set the seed once, outside the loop. Otherwise, I am not sure what you are doing that causes problems. A reproducible example would really help. 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.orgmailto:R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] data.frame with NA
On 18-03-2013, at 19:57, Pietro freeri...@gmail.com wrote: Yes, it's true Berend! What i do is simply use read.xlsx function db - read.xlsx2(c:/mydb.xlsx,1,as.data.frame=T) This is excel file i use: http://dl.dropbox.com/u/102669/mydb.xlsx I don't have the required packages installed to read .xlsx files since I have no need. But I can load your xlsx in LibreOffice without any issues. I'm not sure what read.xslx2 actually does with columns. It seems that is regards the colClass of all columns as character. So maybe you should try read.xlsx or export the file as a .csv Or specify colClasses. Berend I can't find a way to import as numeric. My objective is to be able to work (in R) with my NA's At 18.46 18/03/2013, Berend Hasselman wrote: On 18-03-2013, at 16:49, Pete freeri...@gmail.com wrote: I have this little data.frame http://dl.dropbox.com/u/102669/nanotna.rdata Two column contains NA, so the best thing to do is use na.locf function (with fromLast = T) But locf function doesn't work because NA in my data.frame are not recognized as real NA. Is there a way to substitute fake NA with real NA? In this case na.locf function should work Your data are all characters. Do str(db) to see that. What is probably supposed to be numeric is also character, Somehow you have managed to read in data that R thinks is all chr. Your NA are NA in reality: a character string NA. You will have to review the method you used to get the data into R. And make sure that what you want to be numeric is indeed numeric. Then you can start to think about doing something about the NA's. 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.
[R] ksvm
Hi All, I'm developing a ksvm application. I need to use the results obtained from ksvm in R to construct or restore the ksvm formula and further use this formula to predict new data. The ksvm class contains xmatrix, ymatrix, coefficients, alpha, b and fitted value. Taking linear svm as an example ( vanilladot), the formula is like : F(x) = sign(sum(w*x*y*z)+b) Where, x = xmatrix, y = ymatrix, b is intercept and z is the input matrix? The w is alpha or coefficient from ksvm class? I did some calculation but the result is not the same as the fitted value. Can you give me some idea on the prediction of ksvm by the formula? Thanks a lot!!! Regards, Yan [[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] How to stop set.seed() besides exiting out of R?
Yes, I agree with you. I guess what I was really looking for is a function like UNset.seed()? By having set.seed(), I can have reproducible code. But what if I want to check my work against what's produced from set.seed(100)? I really want to escape from the shadow of set.seed(), can I unset it? On Mon, Mar 18, 2013 at 3:07 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.gov wrote: As I understand it, how R normally does it is to use the system clock to set the seed once per session, unless you use set.seed() to set a new seed. You chose to set the seed to a different value. But from that point on, the pseudo random number generation continues in the same way it normally does. In your code below, each of your 100 histograms will be different. If you then execute the for loop again (but not the set.seed(100) statement), you will get a different set of histograms. The only way you would be confined to set.seed(100) is if you keep resetting the seed to 100. 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 From: C W [mailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 11:50 AM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? set.seed(100) for (i in 1:100){ a - rnorm(1000, mean=0, sd=1) hist(a) } #Now say, I want to simulate without being confined to set.seed(100), I just want to get a simulation like how R normally does it. Mike On Mon, Mar 18, 2013 at 2:40 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.gov wrote: -Original Message- From: r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.org [mailto:r-help-bounces@r-mailto:r-help-bounces@r- project.orghttp://project.org] On Behalf Of C W Sent: Monday, March 18, 2013 11:27 AM To: r-help Subject: [R] How to stop set.seed() besides exiting out of R? Hi list, I am curious how to stop the set.seed(), I don't want the same repeated random number. I know I can set it to a different seed, but I don't want to go through the trouble of setting different seed every time. Thanks, Mike Can you show us how you are using set.seed() that results in getting the same sequence repeatedly? If you are doing simulations in a loop, then set the seed once, outside the loop. Otherwise, I am not sure what you are doing that causes problems. A reproducible example would really help. 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.orgmailto:R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] OrgMassSpecR peak area issue
Hi Bryan, Division by 2 is correct, comes from trapezoid calculation. The modulo line is a funny way of producing c (2 : n, 1) Best, Claudia Am Mon, 18 Mar 2013 15:00:06 -0400 schrieb Bryan Hanson han...@depauw.edu: If you type DrawChromatogram you can see the method used to calculate the peak area. Looks to me like you could easily hack it if you wanted. The relevant part about peak areas is this: for (j in 1:n) { k - (j%%n) + 1 x[j] - peakTime[j] * peakIntensity[k] - peakTime[k] * peakIntensity[j] } peakArea[i] - abs(sum(x)/2) which looks pretty standard to me, though I'm not clear right off the top of my head why they are dividing by 2. You can always contact the maintainer. Bryan On Mar 18, 2013, at 1:34 PM, Christopher Beaver christopher.bea...@gmail.com wrote: Hello! I am having an issue with the OrgMassSpecR package. I run my HPLC using a DAD detector. My raw data is exported form chemstation as a csv file. I then upload the csv into Rstudio no problem. Using the DrawChromatogram function, I get a nice chromatogram, and my retention time, peak area, and apex intensity values are given as well. The problem comes with the peak area value given. The peak area is much smaller than a value that would make sense. My peak area value is actually less than my apex intensity value. Is this because I am using a DAD detector rather than an MS? If so, is there a simply way to edit the peak area equation so that it will also work with absorbance values? Any help is greatly appreciated. Thanks for your time. Chris Beaver [[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. -- Claudia Beleites Spectroscopy/Imaging Institute of Photonic Technology Albert-Einstein-Str. 9 07745 Jena Germany email: claudia.belei...@ipht-jena.de phone: +49 3641 206-133 fax: +49 2641 206-399 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] data.frame with NA
Yes, it's true Berend! What i do is simply use read.xlsx function db - read.xlsx2(c:/mydb.xlsx,1,as.data.frame=T) This is excel file i use: http://dl.dropbox.com/u/102669/mydb.xlsx I can't find a way to import as numeric. My objective is to be able to work (in R) with my NA's At 18.46 18/03/2013, Berend Hasselman wrote: On 18-03-2013, at 16:49, Pete freeri...@gmail.com wrote: I have this little data.frame http://dl.dropbox.com/u/102669/nanotna.rdata Two column contains NA, so the best thing to do is use na.locf function (with fromLast = T) But locf function doesn't work because NA in my data.frame are not recognized as real NA. Is there a way to substitute fake NA with real NA? In this case na.locf function should work Your data are all characters. Do str(db) to see that. What is probably supposed to be numeric is also character, Somehow you have managed to read in data that R thinks is all chr. Your NA are NA in reality: a character string NA. You will have to review the method you used to get the data into R. And make sure that what you want to be numeric is indeed numeric. Then you can start to think about doing something about the NA's. 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] How to stop set.seed() besides exiting out of R?
No, you cannot unset the seed. You can set it to a different value, but a the random number generators always need a starting seed. If you donât set one, R will set one for you , you just wonât know what it is. And as a practical matter, given a sequence of random numbers you canât tell what the starting seed was. That is the point of good random number generators. Each sequence of random numbers for most intents and purposes can be considered independent from previous sets of numbers. Hope this is helpful, 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 From: C W [mailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 12:19 PM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? Yes, I agree with you. I guess what I was really looking for is a function like UNset.seed()? By having set.seed(), I can have reproducible code. But what if I want to check my work against what's produced from set.seed(100)? I really want to escape from the shadow of set.seed(), can I unset it? On Mon, Mar 18, 2013 at 3:07 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.gov wrote: As I understand it, how R âânormallyâ does it is to use the system clock to set the seed once per session, unless you use set.seed() to set a new seed. You chose to set the seed to a different value. But from that point on, the pseudo random number generation continues in the same way it ânormallyâ does. In your code below, each of your 100 histograms will be different. If you then execute the for loop again (but not the set.seed(100) statement), you will get a different set of histograms. The only way you would be âconfined to set.seed(100)â is if you keep resetting the seed to 100. 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 From: C W [mailto:tmrs...@gmail.commailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 11:50 AM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? set.seed(100) for (i in 1:100){ a - rnorm(1000, mean=0, sd=1) hist(a) } #Now say, I want to simulate without being confined to set.seed(100), I just want to get a simulation like how R normally does it. Mike On Mon, Mar 18, 2013 at 2:40 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.govmailto:nord...@dshs.wa.govmailto:nord...@dshs.wa.gov wrote: -Original Message- From: r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.org [mailto:r-help-bounces@r-mailto:r-help-bounces@r-mailto:r-help-bounces@r-mailto:r-help-bounces@r- project.orghttp://project.orghttp://project.org] On Behalf Of C W Sent: Monday, March 18, 2013 11:27 AM To: r-help Subject: [R] How to stop set.seed() besides exiting out of R? Hi list, I am curious how to stop the set.seed(), I don't want the same repeated random number. I know I can set it to a different seed, but I don't want to go through the trouble of setting different seed every time. Thanks, Mike Can you show us how you are using set.seed() that results in getting the same sequence repeatedly? If you are doing simulations in a loop, then set the seed once, outside the loop. Otherwise, I am not sure what you are doing that causes problems. A reproducible example would really help. 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.orgmailto:R-help@r-project.orgmailto:R-help@r-project.orgmailto:R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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.orgmailto:R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Creating SubScale Scores
Hello all, I have just completed a principal components analysis and now which to create sub scales out of my factors. Does anyone know how to do this in R? 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] Counting confidence intervals
Hello, There _is_ a function ?within. Maybe your function can be named 'between' Rui Barradas Em 18-03-2013 16:16, S Ellison escreveu: I want to cont how many times a number say 12 lies in the interval. Can anyone assist? Has anyone else ever wished there was a moderately general 'inside' or 'within' function in R for this problem? For example, something that behaves more or less like within - function(x, interval=NULL, closed=c(TRUE, TRUE), lower=min(interval), upper=max(interval)) { #interval must be a length 2 vector #closed is taken in the order (lower, upper) #lower and upper may be vectors and will be recycled (by etc) if not of length length(x) low.comp - if(closed[1]) = else high.comp - if(closed[2]) = else do.call(low.comp, list(lower, x)) do.call(high.comp, list(upper, x)) } #Examples within(1:5, c(2,4)) within(1:5, c(2,4), closed=c(FALSE, TRUE)) within(1:5, lower=5:1, upper=10:14) S Ellison LGC *** This email and any attachments are confidential. Any use...{{dropped:8}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 stop set.seed() besides exiting out of R?
I am not sure of this, but I think you can unset the seed by removing the dataset .Random.seed from the global environment. E.g., set.seed(1) runif(5) [1] 0.2655087 0.3721239 0.5728534 0.9082078 0.2016819 rm(list=.Random.seed, envir=globalenv()) runif(5) [1] 0.5952379 0.3355091 0.8820192 0.7633754 0.8064312 set.seed(1) runif(5) # same as before [1] 0.2655087 0.3721239 0.5728534 0.9082078 0.2016819 rm(list=.Random.seed, envir=globalenv()) runif(5) # different [1] 0.52023610 0.73407695 0.08824484 0.26977430 0.80089250 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 Nordlund, Dan (DSHS/RDA) Sent: Monday, March 18, 2013 12:44 PM To: C W Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? No, you cannot unset the seed. You can set it to a different value, but a the random number generators always need a starting seed. If you don’t set one, R will set one for you , you just won’t know what it is. And as a practical matter, given a sequence of random numbers you can’t tell what the starting seed was. That is the point of good random number generators. Each sequence of random numbers for most intents and purposes can be considered independent from previous sets of numbers. Hope this is helpful, 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 From: C W [mailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 12:19 PM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? Yes, I agree with you. I guess what I was really looking for is a function like UNset.seed()? By having set.seed(), I can have reproducible code. But what if I want to check my work against what's produced from set.seed(100)? I really want to escape from the shadow of set.seed(), can I unset it? On Mon, Mar 18, 2013 at 3:07 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.gov wrote: As I understand it, how R “‘normally” does it is to use the system clock to set the seed once per session, unless you use set.seed() to set a new seed. You chose to set the seed to a different value. But from that point on, the pseudo random number generation continues in the same way it “normally” does. In your code below, each of your 100 histograms will be different. If you then execute the for loop again (but not the set.seed(100) statement), you will get a different set of histograms. The only way you would be “confined to set.seed(100)” is if you keep resetting the seed to 100. 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 From: C W [mailto:tmrs...@gmail.commailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 11:50 AM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? set.seed(100) for (i in 1:100){ a - rnorm(1000, mean=0, sd=1) hist(a) } #Now say, I want to simulate without being confined to set.seed(100), I just want to get a simulation like how R normally does it. Mike On Mon, Mar 18, 2013 at 2:40 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.govmailto:nord...@dshs.wa.gov mailto:nord...@dshs.wa.gov wrote: -Original Message- From: r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.orgmailto:r- help-boun...@r-project.orgmailto:r-help-boun...@r-project.org [mailto:r-help- bounces@r-mailto:r-help-bounces@r-mailto:r-help-bounces@r-mailto:r-help- bounces@r- project.orghttp://project.orghttp://project.org] On Behalf Of C W Sent: Monday, March 18, 2013 11:27 AM To: r-help Subject: [R] How to stop set.seed() besides exiting out of R? Hi list, I am curious how to stop the set.seed(), I don't want the same repeated random number. I know I can set it to a different seed, but I don't want to go through the trouble of setting different seed every time. Thanks, Mike Can you show us how you are using set.seed() that results in getting the same sequence repeatedly? If you are doing simulations in a loop, then set the seed once, outside the loop. Otherwise, I am not sure what you are doing that causes problems. A reproducible example would really help. 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.orgmailto:R-help@r-project.orgmailto:R-help@r- project.orgmailto:R-help@r-project.org
Re: [R] How to stop set.seed() besides exiting out of R?
Your request is meaningless. The seed itself is effectively overwritten each time a random number is requested. It is only the repeatability of the sequence of random numbers following the set.seed call that is reproducible. You can set the seed to something else, but there is always a seed and it changes as numbers are requested. --- 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 --- Sent from my phone. Please excuse my brevity. C W tmrs...@gmail.com wrote: Yes, I agree with you. I guess what I was really looking for is a function like UNset.seed()? By having set.seed(), I can have reproducible code. But what if I want to check my work against what's produced from set.seed(100)? I really want to escape from the shadow of set.seed(), can I unset it? On Mon, Mar 18, 2013 at 3:07 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.gov wrote: As I understand it, how R ��normally� does it is to use the system clock to set the seed once per session, unless you use set.seed() to set a new seed. You chose to set the seed to a different value. But from that point on, the pseudo random number generation continues in the same way it �normally� does. In your code below, each of your 100 histograms will be different. If you then execute the for loop again (but not the set.seed(100) statement), you will get a different set of histograms. The only way you would be �confined to set.seed(100)� is if you keep resetting the seed to 100. 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 From: C W [mailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 11:50 AM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? set.seed(100) for (i in 1:100){ a - rnorm(1000, mean=0, sd=1) hist(a) } #Now say, I want to simulate without being confined to set.seed(100), I just want to get a simulation like how R normally does it. Mike On Mon, Mar 18, 2013 at 2:40 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.gov wrote: -Original Message- From: r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.org [mailto:r-help-bounces@r-mailto:r-help-bounces@r- project.orghttp://project.org] On Behalf Of C W Sent: Monday, March 18, 2013 11:27 AM To: r-help Subject: [R] How to stop set.seed() besides exiting out of R? Hi list, I am curious how to stop the set.seed(), I don't want the same repeated random number. I know I can set it to a different seed, but I don't want to go through the trouble of setting different seed every time. Thanks, Mike Can you show us how you are using set.seed() that results in getting the same sequence repeatedly? If you are doing simulations in a loop, then set the seed once, outside the loop. Otherwise, I am not sure what you are doing that causes problems. A reproducible example would really help. 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.orgmailto:R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible
Re: [R] How to stop set.seed() besides exiting out of R?
Thanks, William. You read my mind. Cheers, Mike On Mon, Mar 18, 2013 at 4:00 PM, Jeff Newmiller jdnew...@dcn.davis.ca.uswrote: Your request is meaningless. The seed itself is effectively overwritten each time a random number is requested. It is only the repeatability of the sequence of random numbers following the set.seed call that is reproducible. You can set the seed to something else, but there is always a seed and it changes as numbers are requested. --- 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 --- Sent from my phone. Please excuse my brevity. C W tmrs...@gmail.com wrote: Yes, I agree with you. I guess what I was really looking for is a function like UNset.seed()? By having set.seed(), I can have reproducible code. But what if I want to check my work against what's produced from set.seed(100)? I really want to escape from the shadow of set.seed(), can I unset it? On Mon, Mar 18, 2013 at 3:07 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.gov wrote: As I understand it, how R ��normally� does it is to use the system clock to set the seed once per session, unless you use set.seed() to set a new seed. You chose to set the seed to a different value. But from that point on, the pseudo random number generation continues in the same way it �normally� does. In your code below, each of your 100 histograms will be different. If you then execute the for loop again (but not the set.seed(100) statement), you will get a different set of histograms. The only way you would be �confined to set.seed(100)� is if you keep resetting the seed to 100. 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 From: C W [mailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 11:50 AM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? set.seed(100) for (i in 1:100){ a - rnorm(1000, mean=0, sd=1) hist(a) } #Now say, I want to simulate without being confined to set.seed(100), I just want to get a simulation like how R normally does it. Mike On Mon, Mar 18, 2013 at 2:40 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.gov wrote: -Original Message- From: r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.org [mailto:r-help-bounces@r-mailto:r-help-bounces@r- project.orghttp://project.org] On Behalf Of C W Sent: Monday, March 18, 2013 11:27 AM To: r-help Subject: [R] How to stop set.seed() besides exiting out of R? Hi list, I am curious how to stop the set.seed(), I don't want the same repeated random number. I know I can set it to a different seed, but I don't want to go through the trouble of setting different seed every time. Thanks, Mike Can you show us how you are using set.seed() that results in getting the same sequence repeatedly? If you are doing simulations in a loop, then set the seed once, outside the loop. Otherwise, I am not sure what you are doing that causes problems. A reproducible example would really help. 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.orgmailto:R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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,
Re: [R] data.frame with NA
Try this Open the spreadsheet in Excel. Select all of the data click Copy. Don't close Excel. Open R and type the following command: Foglio1 - read.table(clipboard-128, header=TRUE, sep=\t) Now take a look at the structure of the data.frame str(Foglio1) 'data.frame': 1489 obs. of 15 variables: $ Date: Factor w/ 1489 levels 1/10/2002,1/10/2003,..: 1275 1291 1295 1299 1304 1309 1321 1325 1329 1337 ... $ a : num 202 201 202 201 202 ... $ b : num 231 230 230 230 232 ... $ c : num 177 179 181 180 182 ... $ d : num 277 277 276 276 275 ... $ e : num NA NA NA NA NA NA NA NA NA NA ... $ f : num 275 277 279 279 279 ... $ g : num 91.7 90.7 90.8 91.1 91 ... $ h : num 11446 11258 11280 11396 11127 ... $ i : num 388 389 393 392 393 ... $ l : num 93.2 94 92.4 93.4 93.1 ... $ m : num 128 127 126 129 130 ... $ n : num NA NA NA NA NA NA NA NA NA NA ... $ o : num 133 133 133 133 133 ... $ p : num 107 107 107 107 107 ... -- David L Carlson Associate Professor of Anthropology Texas AM University College Station, TX 77843-4352 -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of Pietro Sent: Monday, March 18, 2013 1:57 PM To: Berend Hasselman Cc: r-h...@stat.math.ethz.ch Subject: Re: [R] data.frame with NA Yes, it's true Berend! What i do is simply use read.xlsx function db - read.xlsx2(c:/mydb.xlsx,1,as.data.frame=T) This is excel file i use: http://dl.dropbox.com/u/102669/mydb.xlsx I can't find a way to import as numeric. My objective is to be able to work (in R) with my NA's At 18.46 18/03/2013, Berend Hasselman wrote: On 18-03-2013, at 16:49, Pete freeri...@gmail.com wrote: I have this little data.frame http://dl.dropbox.com/u/102669/nanotna.rdata Two column contains NA, so the best thing to do is use na.locf function (with fromLast = T) But locf function doesn't work because NA in my data.frame are not recognized as real NA. Is there a way to substitute fake NA with real NA? In this case na.locf function should work Your data are all characters. Do str(db) to see that. What is probably supposed to be numeric is also character, Somehow you have managed to read in data that R thinks is all chr. Your NA are NA in reality: a character string NA. You will have to review the method you used to get the data into R. And make sure that what you want to be numeric is indeed numeric. Then you can start to think about doing something about the NA's. 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. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 stop set.seed() besides exiting out of R?
On 18/03/2013 19:59, William Dunlap wrote: I am not sure of this, but I think you can unset the seed by removing the dataset .Random.seed from the global environment. E.g., set.seed(1) runif(5) [1] 0.2655087 0.3721239 0.5728534 0.9082078 0.2016819 rm(list=.Random.seed, envir=globalenv()) runif(5) [1] 0.5952379 0.3355091 0.8820192 0.7633754 0.8064312 set.seed(1) runif(5) # same as before [1] 0.2655087 0.3721239 0.5728534 0.9082078 0.2016819 rm(list=.Random.seed, envir=globalenv()) runif(5) # different [1] 0.52023610 0.73407695 0.08824484 0.26977430 0.80089250 Yes, almost all the time. R does keep a copy in memory when it is using it and writes it back out when done with it: so assuming nothing asynchronous is going on (e.g. a GUI callback) removing .Random.seed will cause the RNG to be re-initialized at next use. 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 Nordlund, Dan (DSHS/RDA) Sent: Monday, March 18, 2013 12:44 PM To: C W Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? No, you cannot unset the seed. You can set it to a different value, but a the random number generators always need a starting seed. If you don’t set one, R will set one for you , you just won’t know what it is. And as a practical matter, given a sequence of random numbers you can’t tell what the starting seed was. That is the point of good random number generators. Each sequence of random numbers for most intents and purposes can be considered independent from previous sets of numbers. Hope this is helpful, 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 From: C W [mailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 12:19 PM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? Yes, I agree with you. I guess what I was really looking for is a function like UNset.seed()? By having set.seed(), I can have reproducible code. But what if I want to check my work against what's produced from set.seed(100)? I really want to escape from the shadow of set.seed(), can I unset it? On Mon, Mar 18, 2013 at 3:07 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.gov wrote: As I understand it, how R “‘normally” does it is to use the system clock to set the seed once per session, unless you use set.seed() to set a new seed. You chose to set the seed to a different value. But from that point on, the pseudo random number generation continues in the same way it “normally” does. In your code below, each of your 100 histograms will be different. If you then execute the for loop again (but not the set.seed(100) statement), you will get a different set of histograms. The only way you would be “confined to set.seed(100)” is if you keep resetting the seed to 100. 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 From: C W [mailto:tmrs...@gmail.commailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 11:50 AM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? set.seed(100) for (i in 1:100){ a - rnorm(1000, mean=0, sd=1) hist(a) } #Now say, I want to simulate without being confined to set.seed(100), I just want to get a simulation like how R normally does it. Mike On Mon, Mar 18, 2013 at 2:40 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.govmailto:nord...@dshs.wa.gov mailto:nord...@dshs.wa.gov wrote: -Original Message- From: r-help-boun...@r-project.orgmailto:r-help-boun...@r-project.orgmailto:r- help-boun...@r-project.orgmailto:r-help-boun...@r-project.org [mailto:r-help- bounces@r-mailto:r-help-bounces@r-mailto:r-help-bounces@r-mailto:r-help- bounces@r- project.orghttp://project.orghttp://project.org] On Behalf Of C W Sent: Monday, March 18, 2013 11:27 AM To: r-help Subject: [R] How to stop set.seed() besides exiting out of R? Hi list, I am curious how to stop the set.seed(), I don't want the same repeated random number. I know I can set it to a different seed, but I don't want to go through the trouble of setting different seed every time. Thanks, Mike Can you show us how you are using set.seed() that results in getting the same sequence repeatedly? If you are doing simulations in a loop, then set the seed once, outside the loop. Otherwise, I am not sure what you are doing that causes problems. A reproducible example would really help. Dan Daniel J. Nordlund Washington State Department of Social and Health Services Planning, Performance, and Accountability Research and Data Analysis Division
Re: [R] how to plot u-v wind by R?
The my.symbols and ms.arrows functions in the TeachingDemos package may help. On Mon, Mar 18, 2013 at 3:50 AM, Jie Tang totang...@gmail.com wrote: hi R users: I have a dataset including u wind in x-axis and v wind in y-axis. How can I plot the u,v wind data in vector or barb figure? which command ? thank you . -- TANG Jie [[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. -- Gregory (Greg) L. Snow Ph.D. 538...@gmail.com [[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] Superscript followed by number then superscript in text
Hi all, I'm having problems finding the correct format for a command. I would like to write some text on a plot. I'm using the following command: text(x,y,text here, srt=90) I would like the text to read: capacity 10^3 m^3 (with ^ denoting superscript (i.e. each '3' as superscript). I've tried fiddling around with expression(paste(etc to no avail. I receive errors such as: Error: unexpected numeric constant... Anyone had experience with this before? Any suggestions would be great. Many thanks in advance, Ben Gillespie Research Postgraduate School of Geography University of Leeds Leeds LS2 9JT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] data.frame with NA
It appears that you MUST use the colClasses= argument with read.xlsx2: Foglio1 - read.xlsx2(mydb.xlsx, 1, colClasses=c(Date, rep(numeric, 14))) However, e and n are converted to NaN not NA so you would need to convert those columns (at least, I didn't check for missing values in the other columns): Foglio1$e - ifelse(is.nan(Foglio1$e), NA, Foglio1$e) Foglio1$n - ifelse(is.nan(Foglio1$n), NA, Foglio1$n) str(Foglio1) 'data.frame': 1489 obs. of 15 variables: $ Date: Date, format: 2001-08-17 2001-08-20 ... $ a : num 202 201 202 201 202 ... $ b : num 231 230 230 230 232 ... $ c : num 177 179 181 180 182 ... $ d : num 277 277 276 276 275 ... $ e : num NA NA NA NA NA NA NA NA NA NA ... $ f : num 275 277 279 279 279 ... $ g : num 91.7 90.7 90.8 91.1 91 ... $ h : num 11446 11258 11280 11396 11127 ... $ i : num 388 389 393 392 393 ... $ l : num 93.2 94 92.4 93.4 93.1 ... $ m : num 128 127 126 129 130 ... $ n : num NA NA NA NA NA NA NA NA NA NA ... $ o : num 133 133 133 133 133 ... $ p : num 107 107 107 107 107 ... --- David -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of David L Carlson Sent: Monday, March 18, 2013 3:22 PM To: 'Pietro'; 'Berend Hasselman' Cc: r-h...@stat.math.ethz.ch Subject: Re: [R] data.frame with NA Try this Open the spreadsheet in Excel. Select all of the data click Copy. Don't close Excel. Open R and type the following command: Foglio1 - read.table(clipboard-128, header=TRUE, sep=\t) Now take a look at the structure of the data.frame str(Foglio1) 'data.frame': 1489 obs. of 15 variables: $ Date: Factor w/ 1489 levels 1/10/2002,1/10/2003,..: 1275 1291 1295 1299 1304 1309 1321 1325 1329 1337 ... $ a : num 202 201 202 201 202 ... $ b : num 231 230 230 230 232 ... $ c : num 177 179 181 180 182 ... $ d : num 277 277 276 276 275 ... $ e : num NA NA NA NA NA NA NA NA NA NA ... $ f : num 275 277 279 279 279 ... $ g : num 91.7 90.7 90.8 91.1 91 ... $ h : num 11446 11258 11280 11396 11127 ... $ i : num 388 389 393 392 393 ... $ l : num 93.2 94 92.4 93.4 93.1 ... $ m : num 128 127 126 129 130 ... $ n : num NA NA NA NA NA NA NA NA NA NA ... $ o : num 133 133 133 133 133 ... $ p : num 107 107 107 107 107 ... -- David L Carlson Associate Professor of Anthropology Texas AM University College Station, TX 77843-4352 -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of Pietro Sent: Monday, March 18, 2013 1:57 PM To: Berend Hasselman Cc: r-h...@stat.math.ethz.ch Subject: Re: [R] data.frame with NA Yes, it's true Berend! What i do is simply use read.xlsx function db - read.xlsx2(c:/mydb.xlsx,1,as.data.frame=T) This is excel file i use: http://dl.dropbox.com/u/102669/mydb.xlsx I can't find a way to import as numeric. My objective is to be able to work (in R) with my NA's At 18.46 18/03/2013, Berend Hasselman wrote: On 18-03-2013, at 16:49, Pete freeri...@gmail.com wrote: I have this little data.frame http://dl.dropbox.com/u/102669/nanotna.rdata Two column contains NA, so the best thing to do is use na.locf function (with fromLast = T) But locf function doesn't work because NA in my data.frame are not recognized as real NA. Is there a way to substitute fake NA with real NA? In this case na.locf function should work Your data are all characters. Do str(db) to see that. What is probably supposed to be numeric is also character, Somehow you have managed to read in data that R thinks is all chr. Your NA are NA in reality: a character string NA. You will have to review the method you used to get the data into R. And make sure that what you want to be numeric is indeed numeric. Then you can start to think about doing something about the NA's. 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. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 stop set.seed() besides exiting out of R?
-Original Message- From: Prof Brian Ripley [mailto:rip...@stats.ox.ac.uk] Sent: Monday, March 18, 2013 1:22 PM To: William Dunlap Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? On 18/03/2013 19:59, William Dunlap wrote: I am not sure of this, but I think you can unset the seed by removing the dataset .Random.seed from the global environment. E.g., set.seed(1) runif(5) [1] 0.2655087 0.3721239 0.5728534 0.9082078 0.2016819 rm(list=.Random.seed, envir=globalenv()) runif(5) [1] 0.5952379 0.3355091 0.8820192 0.7633754 0.8064312 set.seed(1) runif(5) # same as before [1] 0.2655087 0.3721239 0.5728534 0.9082078 0.2016819 rm(list=.Random.seed, envir=globalenv()) runif(5) # different [1] 0.52023610 0.73407695 0.08824484 0.26977430 0.80089250 Yes, almost all the time. R does keep a copy in memory when it is using it and writes it back out when done with it: so assuming nothing asynchronous is going on (e.g. a GUI callback) removing .Random.seed will cause the RNG to be re-initialized at next use. Thanks Brian, S+'s set.seed() generates a new seed if it is called with no arguments so as not to tempt people to directly fiddle around with magic datasets like .Random.seed. Bill Dunlap Spotfire, TIBCO Software wdunlap tibco.com 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 Nordlund, Dan (DSHS/RDA) Sent: Monday, March 18, 2013 12:44 PM To: C W Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? No, you cannot unset the seed. You can set it to a different value, but a the random number generators always need a starting seed. If you don’t set one, R will set one for you , you just won’t know what it is. And as a practical matter, given a sequence of random numbers you can’t tell what the starting seed was. That is the point of good random number generators. Each sequence of random numbers for most intents and purposes can be considered independent from previous sets of numbers. Hope this is helpful, 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 From: C W [mailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 12:19 PM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? Yes, I agree with you. I guess what I was really looking for is a function like UNset.seed()? By having set.seed(), I can have reproducible code. But what if I want to check my work against what's produced from set.seed(100)? I really want to escape from the shadow of set.seed(), can I unset it? On Mon, Mar 18, 2013 at 3:07 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.gov wrote: As I understand it, how R “‘normally” does it is to use the system clock to set the seed once per session, unless you use set.seed() to set a new seed. You chose to set the seed to a different value. But from that point on, the pseudo random number generation continues in the same way it “normally” does. In your code below, each of your 100 histograms will be different. If you then execute the for loop again (but not the set.seed(100) statement), you will get a different set of histograms. The only way you would be “confined to set.seed(100)” is if you keep resetting the seed to 100. 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 From: C W [mailto:tmrs...@gmail.commailto:tmrs...@gmail.com] Sent: Monday, March 18, 2013 11:50 AM To: Nordlund, Dan (DSHS/RDA) Cc: r-help Subject: Re: [R] How to stop set.seed() besides exiting out of R? set.seed(100) for (i in 1:100){ a - rnorm(1000, mean=0, sd=1) hist(a) } #Now say, I want to simulate without being confined to set.seed(100), I just want to get a simulation like how R normally does it. Mike On Mon, Mar 18, 2013 at 2:40 PM, Nordlund, Dan (DSHS/RDA) nord...@dshs.wa.govmailto:nord...@dshs.wa.govmailto:nord...@dshs.wa.gov mailto:nord...@dshs.wa.gov wrote: -Original Message- From: r-help-boun...@r-project.orgmailto:r-help-bounces@r- project.orgmailto:r- help-boun...@r-project.orgmailto:r-help-boun...@r-project.org [mailto:r-help- bounces@r-mailto:r-help-bounces@r-mailto:r-help-bounces@r-mailto:r-help- bounces@r- project.orghttp://project.orghttp://project.org] On Behalf Of C W Sent: Monday, March 18, 2013 11:27 AM To: r-help Subject: [R] How to stop set.seed() besides exiting out of R? Hi list, I
Re: [R] Superscript followed by number then superscript in text
On Tue, Mar 19, 2013 at 9:29 AM, Benjamin Gillespie gy...@leeds.ac.ukwrote: Hi all, I'm having problems finding the correct format for a command. I would like to write some text on a plot. I'm using the following command: text(x,y,text here, srt=90) I would like the text to read: capacity 10^3 m^3 (with ^ denoting superscript (i.e. each '3' as superscript). What did you try? Anyhow, this works text(1,1,expression(capacity~10^3~m^3)) -thomas -- Thomas Lumley Professor of Biostatistics University of Auckland [[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] Superscript followed by number then superscript in text
Hello, Something like this? plot(1, type = n) text(1,1, expression(capacity ~ 10^3 ~ m^3)) Hope this helps, Rui Barradas Em 18-03-2013 20:29, Benjamin Gillespie escreveu: Hi all, I'm having problems finding the correct format for a command. I would like to write some text on a plot. I'm using the following command: text(x,y,text here, srt=90) I would like the text to read: capacity 10^3 m^3 (with ^ denoting superscript (i.e. each '3' as superscript). I've tried fiddling around with expression(paste(etc to no avail. I receive errors such as: Error: unexpected numeric constant... Anyone had experience with this before? Any suggestions would be great. Many thanks in advance, Ben Gillespie Research Postgraduate School of Geography University of Leeds Leeds LS2 9JT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] melt with complications
HI, Not sure whether this is what you wanted. library(reshape2) result.dcast-dcast(meltTest,A+D~H,value.var=M) attr(result.reshape,reshapeWide)-NULL row.names(result.reshape)-1:nrow(result.reshape) identical(result.reshape,result.dcast) #[1] TRUE result.dcast # A D I J K L #1 B E 1 2 3 4 #2 B F 5 6 7 8 #3 B G 9 10 11 12 #4 C E 13 14 15 16 #5 C F 17 18 19 20 #6 C G 21 22 23 24 A.K. - Original Message - From: Richard M. Heiberger r...@temple.edu To: r-help r-help@r-project.org; Hadley Wickham h.wick...@gmail.com Cc: Sent: Monday, March 18, 2013 1:15 PM Subject: [R] melt with complications ## Can someone suggest a simpler expression than either of these, with the goal ## of taking a long matrix into a wide one with exactly one of the factors converted to ## columns and all the rest retained as factors. I want something that generalizes beyond ## the three factors illustrated here. ## Rich meltTest - data.frame(A=rep(c(B,C), each=12), D=rep(c(E,F,G), each=4, times=2), H=rep(c(I,J,K,L), times=6), M=1:24) meltTest result.melt - do.call(rbind, { tmp - cast(D ~ H | A, value=M, data=meltTest) lapply(names(tmp), function(x) cbind(A=x, tmp[[x]])) ## explicit use of name A }) result.melt result.reshape - reshape(meltTest, direction=wide, timevar=H, idvar=c(A,D)) names(result.reshape)[3:6] - unique(as.character(meltTest$H)) ## explicit use of name H result.reshape [[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] melt with complications
Richard M. Heiberger rmh at temple.edu writes: ## Can someone suggest a simpler expression than either of these, with the goal ## of taking a long matrix into a wide one with exactly one of the factors converted to ## columns and all the rest retained as factors. I want something that generalizes beyond ## the three factors illustrated here. ## Rich meltTest - data.frame(A=rep(c(B,C), each=12), D=rep(c(E,F,G), each=4, times=2), H=rep(c(I,J,K,L), times=6), M=1:24) amat - ftable( xtabs( M ~ A + D + H, meltTest ),row.vars=1:2 ) amat is such a matrix with a few attributes added. HTH, Chuck __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Fit a mixture of lognormal and normal distributions
Hello I am trying to find an automated way of fitting a mixture of normal and log-normal distributions to data which is clearly bimodal. Here's a simulated example: x.1-rnorm(6000, 2.4, 0.6) x.2-rlnorm(1, 1.3,0.1) X-c(x.1, x.2) hist(X,100,freq=FALSE, ylim=c(0,1.5))lines(density(x.1), lty=2, lwd=2) lines(density(x.2), lty=2, lwd=2) lines(density(X), lty=4) Currently i am using mixtools and just run: library(mixtools) mixmdl = normalmixEM(X, k=2, epsilon = 1e-08, maxit = 1000, maxrestarts=20, verb = TRUE, fast=FALSE, ECM = FALSE, arbmean = TRUE, arbvar = TRUE) plot(mixmdl,which=2)lines(density(X), lty=2, lwd=2) This is obviously not the best way of doing this. The estimates it gives are: mu 3.6595737(x.2 log()=1.29) 2.3113135(x.1) These are not too far off but I was wondering if someone knows of a better R package/way of doing this or has any hints that would help me coding it from scratch (EM+writing down the formulae for ML? but... would this really be better than using a more-advanced-already-available-package-made-by-pros?). The objective is to estimate threshold values at specific FDRs. (some help with that would also be most helpful.) Thanks to all in advance! To' __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Fit a mixture of lognormal and normal distributions
You posted this earlier. Do not repost, please. It was seen. -- Bert On Mon, Mar 18, 2013 at 3:28 PM, To . . kid...@hotmail.com wrote: Hello I am trying to find an automated way of fitting a mixture of normal and log-normal distributions to data which is clearly bimodal. Here's a simulated example: x.1-rnorm(6000, 2.4, 0.6) x.2-rlnorm(1, 1.3,0.1) X-c(x.1, x.2) hist(X,100,freq=FALSE, ylim=c(0,1.5))lines(density(x.1), lty=2, lwd=2) lines(density(x.2), lty=2, lwd=2) lines(density(X), lty=4) Currently i am using mixtools and just run: library(mixtools) mixmdl = normalmixEM(X, k=2, epsilon = 1e-08, maxit = 1000, maxrestarts=20, verb = TRUE, fast=FALSE, ECM = FALSE, arbmean = TRUE, arbvar = TRUE) plot(mixmdl,which=2)lines(density(X), lty=2, lwd=2) This is obviously not the best way of doing this. The estimates it gives are: mu 3.6595737(x.2 log()=1.29) 2.3113135(x.1) These are not too far off but I was wondering if someone knows of a better R package/way of doing this or has any hints that would help me coding it from scratch (EM+writing down the formulae for ML? but... would this really be better than using a more-advanced-already-available-package-made-by-pros?). The objective is to estimate threshold values at specific FDRs. (some help with that would also be most helpful.) Thanks to all in advance! To' __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Bert Gunter Genentech Nonclinical Biostatistics Internal Contact Info: Phone: 467-7374 Website: http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm [[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 web application development
Dear all, I am wondering if what would be the simple way to develop a simple web application that runs R. That is, the web application allows any user upload a dataframe as a variable to my web server, a linux-based apache, and then run a R package (my package) on the variable that should ideally be handled as a variable in memory instead of saving to the disk in my server , for security concern. After running, the result would be returned to the web. Any suggestion will be appreciated. Best, John [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] R web application development
I can't offer any advice, but I feel like you could probably get a good start on this by looking through the archives. On Mon 18 Mar 2013 06:55:33 PM CDT, John linux-user wrote: Dear all, I am wondering if what would be the simple way to develop a simple web application that runs R. That is, the web application allows any user upload a dataframe as a variable to my web server, a linux-based apache, and then run a R package (my package) on the variable that should ideally be handled as a variable in memory instead of saving to the disk in my server , for security concern. After running, the result would be returned to the web. Any suggestion will be appreciated. Best, John [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] How to automate this model selection algorithm?
I've got a complicated semi-parametric model that I'm fitting with mgcv. I start with a model based on theory. Its got lots of interaction terms. I want to winnow it down: removing each interaction term or un-interacted main effect one by one, checking the AIC, and retaining the model that gives me the lowest AIC. I then want to repeat the procedure on the retained model. Here is a simple example: set.seed(42) library(mgcv) N=220 fac = ceiling(runif(N)*2) a = rnorm(N); b = rnorm(N); c = rnorm(N); d = runif(N); e = rnorm(N); y = a^fac + b*e + d/(e+1) m1 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(a,b,c) +te(a,d,by=as.factor(fac)) ) m2 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(b,c) +te(a,d,by=as.factor(fac)) ) m3 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(a,c) +te(a,d,by=as.factor(fac)) ) m4 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(a,b) +te(a,d,by=as.factor(fac)) ) m5 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(a,b,c) +te(d,by=as.factor(fac)) ) m6 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(a,b,c) +te(a,by=as.factor(fac)) ) m7 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(a,b,c) +te(a,d) ) selection = AIC(m1,m2,m3,m4,m5,m6,m7) selection df AIC m1 14.53435 1626.611 m2 12.52501 1622.635 m3 12.54566 1622.615 m4 12.52652 1622.633 m5 13.14108 1620.759 m6 10.99684 1621.156 m7 10.98136 1622.229 m5 is the best m5 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(a,b,c) +te(d,by=as.factor(fac)) ) m5.1 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(b,c) +te(d,by=as.factor(fac)) ) m5.2 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(a,c) +te(d,by=as.factor(fac)) ) m5.3 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(a,b) +te(d,by=as.factor(fac)) ) m5.4 = gam(y~ as.factor(fac) +s(a) +s(b) +s(c) +s(d) +te(a,b,c) #+te(d,by=as.factor(fac)) ) selection.2 = AIC(m5,m5.1,m5.2,m5.3,m5.4) selection.2 df AIC m5 13.363029 1621.183 m5.1 9.671656 1617.641 m5.2 9.730047 1617.549 m5.3 9.706424 1617.569 m5.4 9.857504 1620.028 5.2 is the best ...and so on. I'd next try out each interaction term or un-interacted main effect in m5.2. Question is how I can automate this? To start with a model (m1, in my case), and have R give me the best model after running this algorithm to the point where there are no longer any better models according to this algorithm? Right now I can do it manually, but as I add data over time, model selection may change. Note the mgcv's select=TRUE functionality doesn't quite work for me (at least not directly), because I want to see if single parts of three-way interactions can/should be removed. Thanks in advance for any tips, and apologies for cross-posting (on stack overflow). Best, Andrew __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 web application development
Thanks for reply, but which archives? Thanks again. John From: Stephen Sefick ssef...@gmail.com Cc: r-help@r-project.org r-help@r-project.org Sent: Monday, March 18, 2013 8:02 PM Subject: Re: [R] R web application development I can't offer any advice, but I feel like you could probably get a good start on this by looking through the archives. On Mon 18 Mar 2013 06:55:33 PM CDT, John linux-user wrote: Dear all, I am wondering if what would be the simple way to develop a simple web application that runs R. That is, the web application allows any user upload a dataframe as a variable to my web server, a linux-based apache, and then run a R package (my package) on the variable that should ideally be handled as a variable in memory instead of saving to the disk in my server , for security concern. After running, the result would be returned to the web. Any suggestion will be appreciated. Best, John [[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] melt with complications
thank you all. the dcast formulation is what I am looking for. That version of the formula A+D ~ H also works with the original cast function. result.cast - cast(A+D ~ H, value=M, data=meltTest) It generalizes to a four-factor example, as requested. This ftable solution didn't do what I wanted because it didn't retain A and D as factors. On Mon, Mar 18, 2013 at 6:00 PM, Charles Berry ccbe...@ucsd.edu wrote: Richard M. Heiberger rmh at temple.edu writes: ## Can someone suggest a simpler expression than either of these, with the goal ## of taking a long matrix into a wide one with exactly one of the factors converted to ## columns and all the rest retained as factors. I want something that generalizes beyond ## the three factors illustrated here. ## Rich meltTest - data.frame(A=rep(c(B,C), each=12), D=rep(c(E,F,G), each=4, times=2), H=rep(c(I,J,K,L), times=6), M=1:24) amat - ftable( xtabs( M ~ A + D + H, meltTest ),row.vars=1:2 ) amat is such a matrix with a few attributes added. HTH, Chuck __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] How to automate this model selection algorithm?
Just to add, I've been playing around with select=TRUE in mgcv, and it does seem that it could work if I were to specify all of the nested two-way interactions in my three-way interactions (see the toy example below). But the problem is that I don't have enough degrees of freedom to feed such a model into GAM using my main dataset. N=200 a = rnorm(N) b = rnorm(N) c = rnorm(N) y = rnorm(N)+a+b+c+a*b m = gam(y~s(a)+s(b)+s(c)+te(a,b)+te(a,b,c)) msel = gam(y~s(a)+s(b)+s(c)+te(a,b)+te(a,b,c),select=TRUE) mdrop = gam(y~s(a)+s(b)+s(c)+te(a,b)) summary(m) summary(msel) summary(mdrop) plot(density(m$fitted.values)) lines(density(msel$fitted.values),col=red) lines(density(mdrop$fitted.values),col=blue) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 web application development
I generally use the Markmail search engine but others use the Newcastle resource. Some of other ones can be found with perusal of the Posting Guide and Mailing list info page. http://markmail.org/search/?q=list%3Aorg.r-project.r-help On Mar 18, 2013, at 5:11 PM, John linux-user wrote: Thanks for reply, but which archives? Thanks again. John From: Stephen Sefick ssef...@gmail.com Cc: r-help@r-project.org r-help@r-project.org Sent: Monday, March 18, 2013 8:02 PM Subject: Re: [R] R web application development I can't offer any advice, but I feel like you could probably get a good start on this by looking through the archives. On Mon 18 Mar 2013 06:55:33 PM CDT, John linux-user wrote: Dear all, I am wondering if what would be the simple way to develop a simple web application that runs R. That is, the web application allows any user upload a dataframe as a variable to my web server, a linux-based apache, and then run a R package (my package) on the variable that should ideally be handled as a variable in memory instead of saving to the disk in my server , for security concern. After running, the result would be returned to the web. Any suggestion will be appreciated. Best, John [[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. David Winsemius Alameda, CA, USA __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] RGL 3D plots are flat. Please Help
Hello, I have a matrix of simulated data that I want to plot as a 3D surface using RGL. Have followed the documentation carefully, but my plot only contains a weird single slice of the data. The actual RGL library and dependance seem fine as all of the demo code plots beautifully. The actual command I'm using is: rgl.surface(x,z,d, color=blue,alpha=0.5,shininess=128) The length of the x and z vectors match the dimensions of the y matrix. Can't figure out why I'm only getting a slice and not a full 3D surface. Without digging into too much detail, here is what my data looks like. x [1]0 100 200 300 400 500 600 700 800 900 1000 1100 1200 1300 1400 1500 [17] 1600 1700 1800 1900 2000 2100 2200 2300 2400 2500 z [1] 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 [21] 2.0 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3.0 y (Only first bit included here.) [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14] [1,]000000000 0 0 0 0 0 [2,]023568 10 11 131416181921 [3,]036 10 13 16 19 22 262932353842 [4,]05 10 14 19 24 29 34 384348535862 [5,]06 13 19 26 32 38 45 515864707783 [6,]08 16 24 32 40 48 56 6472808896 104 [7,]0 10 19 29 38 48 58 67 778696 106 115 125 [8,]0 11 22 34 45 56 67 78 90 101 112 123 134 146 [9,]0 13 26 38 51 64 77 90 102 115 128 141 154 166 [10,]0 14 29 43 58 72 86 101 115 130 144 158 173 187 [11,]0 16 32 48 64 80 96 112 128 144 160 176 192 208 [12,]0 18 35 53 70 88 106 123 141 158 176 194 211 229 [13,]0 19 38 58 77 96 115 134 154 173 192 211 230 250 [14,]0 21 42 62 83 104 125 146 166 187 208 229 250 270 [15,]0 22 45 67 90 112 134 157 179 202 224 246 269 291 [16,]0 24 48 72 96 120 144 168 192 216 240 264 288 312 [17,]0 26 51 77 102 128 154 179 205 230 256 282 307 333 [18,]0 27 54 82 109 136 163 190 218 245 272 299 326 354 [19,]0 29 58 86 115 144 173 202 230 259 288 317 346 374 [20,]0 30 61 91 122 152 182 213 243 274 304 334 365 395 [21,]0 32 64 96 128 160 192 224 256 288 320 352 384 416 [22,]0 34 67 101 134 168 202 235 269 302 336 370 403 437 [23,]0 35 70 106 141 176 211 246 282 317 352 387 422 458 [24,]0 37 74 110 147 184 221 258 294 331 368 405 442 478 [25,]0 38 77 115 154 192 230 269 307 346 384 422 461 499 [26,]0 40 80 120 160 200 240 280 320 360 400 440 480 520 -- Noah Silverman, M.S. UCLA Department of Statistics 8117 Math Sciences Building Los Angeles, CA 90095 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] RGL 3D plots are flat. Please Help
Hi, Your example is not reproducible. What is d in your command? Regards, Pascal On 19/03/13 10:21, Noah Silverman wrote: Hello, I have a matrix of simulated data that I want to plot as a 3D surface using RGL. Have followed the documentation carefully, but my plot only contains a weird single slice of the data. The actual RGL library and dependance seem fine as all of the demo code plots beautifully. The actual command I'm using is: rgl.surface(x,z,d, color=blue,alpha=0.5,shininess=128) The length of the x and z vectors match the dimensions of the y matrix. Can't figure out why I'm only getting a slice and not a full 3D surface. Without digging into too much detail, here is what my data looks like. x [1]0 100 200 300 400 500 600 700 800 900 1000 1100 1200 1300 1400 1500 [17] 1600 1700 1800 1900 2000 2100 2200 2300 2400 2500 z [1] 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 [21] 2.0 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3.0 y (Only first bit included here.) [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14] [1,]000000000 0 0 0 0 0 [2,]023568 10 11 1314161819 21 [3,]036 10 13 16 19 22 2629323538 42 [4,]05 10 14 19 24 29 34 3843485358 62 [5,]06 13 19 26 32 38 45 5158647077 83 [6,]08 16 24 32 40 48 56 6472808896 104 [7,]0 10 19 29 38 48 58 67 778696 106 115 125 [8,]0 11 22 34 45 56 67 78 90 101 112 123 134 146 [9,]0 13 26 38 51 64 77 90 102 115 128 141 154 166 [10,]0 14 29 43 58 72 86 101 115 130 144 158 173 187 [11,]0 16 32 48 64 80 96 112 128 144 160 176 192 208 [12,]0 18 35 53 70 88 106 123 141 158 176 194 211 229 [13,]0 19 38 58 77 96 115 134 154 173 192 211 230 250 [14,]0 21 42 62 83 104 125 146 166 187 208 229 250 270 [15,]0 22 45 67 90 112 134 157 179 202 224 246 269 291 [16,]0 24 48 72 96 120 144 168 192 216 240 264 288 312 [17,]0 26 51 77 102 128 154 179 205 230 256 282 307 333 [18,]0 27 54 82 109 136 163 190 218 245 272 299 326 354 [19,]0 29 58 86 115 144 173 202 230 259 288 317 346 374 [20,]0 30 61 91 122 152 182 213 243 274 304 334 365 395 [21,]0 32 64 96 128 160 192 224 256 288 320 352 384 416 [22,]0 34 67 101 134 168 202 235 269 302 336 370 403 437 [23,]0 35 70 106 141 176 211 246 282 317 352 387 422 458 [24,]0 37 74 110 147 184 221 258 294 331 368 405 442 478 [25,]0 38 77 115 154 192 230 269 307 346 384 422 461 499 [26,]0 40 80 120 160 200 240 280 320 360 400 440 480 520 -- Noah Silverman, M.S. UCLA Department of Statistics 8117 Math Sciences Building Los Angeles, CA 90095 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] RGL 3D plots are flat. Please Help
Oops, that was a type. d is just y. If you want to reproduce, I put the entire matrix, as a csv here: http://pastebin.com/gniyD4Rc Thanks, -- Noah Silverman, M.S. UCLA Department of Statistics 8117 Math Sciences Building Los Angeles, CA 90095 On Mar 18, 2013, at 6:58 PM, Pascal Oettli kri...@ymail.com wrote: Hi, Your example is not reproducible. What is d in your command? Regards, Pascal On 19/03/13 10:21, Noah Silverman wrote: Hello, I have a matrix of simulated data that I want to plot as a 3D surface using RGL. Have followed the documentation carefully, but my plot only contains a weird single slice of the data. The actual RGL library and dependance seem fine as all of the demo code plots beautifully. The actual command I'm using is: rgl.surface(x,z,d, color=blue,alpha=0.5,shininess=128) The length of the x and z vectors match the dimensions of the y matrix. Can't figure out why I'm only getting a slice and not a full 3D surface. Without digging into too much detail, here is what my data looks like. x [1]0 100 200 300 400 500 600 700 800 900 1000 1100 1200 1300 1400 1500 [17] 1600 1700 1800 1900 2000 2100 2200 2300 2400 2500 z [1] 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 [21] 2.0 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3.0 y (Only first bit included here.) [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14] [1,]000000000 0 0 0 0 0 [2,]023568 10 11 1314161819 21 [3,]036 10 13 16 19 22 2629323538 42 [4,]05 10 14 19 24 29 34 3843485358 62 [5,]06 13 19 26 32 38 45 5158647077 83 [6,]08 16 24 32 40 48 56 6472808896 104 [7,]0 10 19 29 38 48 58 67 778696 106 115 125 [8,]0 11 22 34 45 56 67 78 90 101 112 123 134 146 [9,]0 13 26 38 51 64 77 90 102 115 128 141 154 166 [10,]0 14 29 43 58 72 86 101 115 130 144 158 173 187 [11,]0 16 32 48 64 80 96 112 128 144 160 176 192 208 [12,]0 18 35 53 70 88 106 123 141 158 176 194 211 229 [13,]0 19 38 58 77 96 115 134 154 173 192 211 230 250 [14,]0 21 42 62 83 104 125 146 166 187 208 229 250 270 [15,]0 22 45 67 90 112 134 157 179 202 224 246 269 291 [16,]0 24 48 72 96 120 144 168 192 216 240 264 288 312 [17,]0 26 51 77 102 128 154 179 205 230 256 282 307 333 [18,]0 27 54 82 109 136 163 190 218 245 272 299 326 354 [19,]0 29 58 86 115 144 173 202 230 259 288 317 346 374 [20,]0 30 61 91 122 152 182 213 243 274 304 334 365 395 [21,]0 32 64 96 128 160 192 224 256 288 320 352 384 416 [22,]0 34 67 101 134 168 202 235 269 302 336 370 403 437 [23,]0 35 70 106 141 176 211 246 282 317 352 387 422 458 [24,]0 37 74 110 147 184 221 258 294 331 368 405 442 478 [25,]0 38 77 115 154 192 230 269 307 346 384 422 461 499 [26,]0 40 80 120 160 200 240 280 320 360 400 440 480 520 -- Noah Silverman, M.S. UCLA Department of Statistics 8117 Math Sciences Building Los Angeles, CA 90095 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Convert to date and time of the year
Dear R Users, I have data for more than 3 years. For each year I want to find the day corresponding to Jaunary 1 of that year. For example: x - c('5/5/2007','12/31/2007','1/2/2008') #Convert to day of year (julian date) - strptime(x,%m/%d/%Y)$yday+1 [1] 125 365 2 I want to know how to do the same thing but with time added. But I still get the day not time. Can anyone suggest what is the better way to find the julian date with date and time ? x1 - c('5/5/2007 02:00','12/31/2007 05:58','1/2/2008 16:25') #Convert to day of year (julian date) - strptime(x1,%m/%d/%Y %H:%M)$yday+1 [1] 125 365 2 Thank you so much. Janesh [[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] Convert to date and time of the year
Hi Janesh, On Mar 18, 2013, at 10:49 PM, Janesh Devkota wrote: Dear R Users, I have data for more than 3 years. For each year I want to find the day corresponding to Jaunary 1 of that year. For example: x - c('5/5/2007','12/31/2007','1/2/2008') #Convert to day of year (julian date) - strptime(x,%m/%d/%Y)$yday+1 [1] 125 365 2 I want to know how to do the same thing but with time added. But I still get the day not time. Can anyone suggest what is the better way to find the julian date with date and time ? x1 - c('5/5/2007 02:00','12/31/2007 05:58','1/2/2008 16:25') #Convert to day of year (julian date) - strptime(x1,%m/%d/%Y %H:%M)$yday+1 [1] 125 365 2 Would this do it for you? x1 - c('5/5/2007 02:00','12/31/2007 05:58','1/2/2008 16:25') x1 - as.POSIXct(x1, format = %m/%d/%Y %H:%M) day - as.numeric(format(x1, %j)) hour - as.numeric(format(x1, %H))/24 minute - as.numeric(format(x1, %M))/(60*24) second - as.numeric(format(x1, %S))/(60*60*24) day + hour + minute + second [1] 125.08 365.248611 2.684028 Cheers, Ben Thank you so much. Janesh [[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. Ben Tupper Bigelow Laboratory for Ocean Sciences 60 Bigelow Drive, P.O. Box 380 East Boothbay, Maine 04544 http://www.bigelow.org __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.