Re: [R] return one vector that is the sum of a list of vectors
On Apr 9, 2012, at 12:01 AM, Frank Tamborello wrote: Dear R Help, I am attempting to write a function that takes a list of variable groups and a vector of numbers (e.g., jtsv would indicate one group of variables, jtsv1, jtsv2, etc, each of which name a variable in a data frame), and returns a list of variable group sums. For example, given list(jtsv, ptsv) and c(1:10), where jtsv1, etc, are all numeric vectors of the same length, the function should return a list of jtsv and ptsv, where jtsv - jtsv1 + jtsv2 + jtsv3 + jtsv4 + jtsv5 + jtsv6 + jtsv7 + jtsv8 + jtsv9 + jtsv10 ptsv - ptsv1 + ptsv2 + ptsv3 + ptsv4 + ptsv5 + ptsv6 + ptsv7 + ptsv8 + ptsv9 + ptsv10 So far I've used a for loop to paste together the names of the variables as character vectors, so that I can at least name the variables that I want to get, but I'm stuck on how to add together the numeric vectors that those character vectors are ultimately supposed to name. It seems like I should be able to map the primitive + operator onto a list or vector that contains the variable names. If that's a good way to go, what would that look like? Or if that's not a good way to go, would someone please point me in a good direction? Perhaps ( depending on detail not provided): jtsv - rowSums( dfrm[ , grep(jtsv, names(dfrm))] ) -- David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Error in integrate(int.fn, lower = 0, upper = Inf) : evaluation of function gave a result of wrong length
Sorry but i didn't understand. Thanks for your answer. -- View this message in context: http://r.789695.n4.nabble.com/Error-in-integrate-int-fn-lower-0-upper-Inf-evaluation-of-function-gave-a-result-of-wrong-length-tp4541250p4542101.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] Creating a dummy variable
HI I want to create a dummy variable for a dataset which has various dates Eg I've one year data. The user can choose any date range startdate : 2012-01-01 enddate: 2012-02-01 startdate: 2012-06-01 enddate=2012-07-01 the columns of dataset= X ,Y ,date Thanks and re - Thanks in Advance Arun -- View this message in context: http://r.789695.n4.nabble.com/Creating-a-dummy-variable-tp4542150p4542150.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Summary Statistics Help
I've got this solved via Talks Stat mod.1-lm(Patents~FHouse, data=datpat) summary(mod.1) anova(mod.1) xtable(mod.1) -- View this message in context: http://r.789695.n4.nabble.com/Summary-Statistics-Help-tp4541923p4542103.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problems with CRAN package Ryacas
check your firewall? --- 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. Kjetil Halvorsen kjetilbrinchmannhalvor...@gmail.com wrote: First: sessionInfo() R version 2.15.0 (2012-03-30) Platform: x86_64-pc-linux-gnu (64-bit) locale: [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C [3] LC_TIME=en_US.UTF-8LC_COLLATE=en_US.UTF-8 [5] LC_MONETARY=en_US.UTF-8LC_MESSAGES=en_US.UTF-8 [7] LC_PAPER=C LC_NAME=C [9] LC_ADDRESS=C LC_TELEPHONE=C [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] Ryacas_0.2-11 XML_3.9-4 pnmath_0.0-3 MASS_7.3-17 loaded via a namespace (and not attached): [1] compiler_2.15.0 fortunes_1.5-0 tools_2.15.0 Operating system is debian (wheezy+some sid) I have installed yacas version 1.3.1 via synaptic. Then: library(Ryacas) yacas(expression(Factor(x^2-1))) # copied from yacas help page [1] Starting Yacas! Accepting requests from port 9734 Error in socketConnection(host = 127.0.0.1, port = 9734, server = FALSE, : cannot open the connection In addition: Warning message: In socketConnection(host = 127.0.0.1, port = 9734, server = FALSE, : 127.0.0.1:9734 cannot be opened ¿Any ideas? Kjetil __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] suggestions argument in rbga function in genalg package
On Thu, Sep 22, 2011 at 06:46:23PM +, Joseph Boyer wrote: Would someone be so kind as to provide example code where they use the suggestions argument in the rgba function In genalg? I can't get it to work. The following code works just fine: GenFit -rbga(Lower, Upper, evalFunc = evaluate) Lower and Upper are each numeric vectors with 7 elements. Evaluate is an objective function. However, when I want to use a suggested chromosome, I get an error message. My code is start - c(1,0.1,10, 100,1,100,1) suggestions - list(start) GenFit -rbga(Lower, Upper, suggestions = suggestions, evalFunc = evaluate) The error message is: Error in 1:suggestionCount : argument of length 0 Thanks. rbga( c( -1, -1), c( 1, 1), evalFunc= function( string){ sum( string**2)**0.5}, suggestions=t( as.data.frame( c( 0, 0 Regards __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Comparing 2 means. pls help
I am interested in the difference of 2 data: mat1= c(2.2, 2.3, 2.2,2.5) mat2= c(2.6, 2.8, 2.7,2.4) mat= mat2-mat1 I perform an action on both mat1 and mat2, and I get mat1prime and mat2prime: mat1prime= c(2.5, 2.5, 2.3,2.5) mat2prime= c(2.6, 2.8, 2.7,2.6) matprime= mat2prime-mat1prime I want to check wether mat and matprime have the same mean. How should I do this in R? -- View this message in context: http://r.789695.n4.nabble.com/Comparing-2-means-pls-help-tp4542237p4542237.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] slanted stacked bar graphs?
On Mon, Apr 9, 2012 at 7:29 AM, Susanna Makela susanna.m.mak...@gmail.com wrote: Hello R users, I would like to generate slanted stacked bar graphs like those on the bottom of pages 1 and 2 in this document: http://www.wssinfo.org/fileadmin/user_upload/resources/JMP-Snapshot-SWA-HLM.pdf . I've also attached the file to this email (pdf). Does anyone know if this is possible in R? I have tried googling and searching the R help archives, and it seems like ggplot2 might be able to make such graphs, but I'm not familiar enough with graphics in R to know for sure. (I personally don't feel that these slanted bar graphs - not sure if they have an actual name - convey the intended information very well, but I have to try and make them all the same. However, I am open to alternative suggestions for visualizing similar data if anyone has ideas.) These exact charts have been critiqued on the Junk Charts blog: http://junkcharts.typepad.com/junk_charts/2010/02/cousin-misfit.html and you'll even find some ggplot code in the comments for doing them. If you still want to... I just did a google image search for 'ggplot stacked' and there they were. Barry -- blog: http://geospaced.blogspot.com/ web: http://www.maths.lancs.ac.uk/~rowlings web: http://www.rowlingson.com/ twitter: http://twitter.com/geospacedman pics: http://www.flickr.com/photos/spacedman __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Drawing a line in xyplot
PS: The following shows possibilities that are available using latticeExtra layering: ## Best make type and attend factors xdat = data.frame(mortality =c(5, 8,7,5,8,10,11,6,4,5,20,25,27,30,35,32,28,21,20,34,11,15,18,12,15,12,10,15,19,20), type= factor(c(1, 1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2,3,3,3,3,3,3,3,3,3,3)), attend = factor(c(1, 0,1,1,1,1,0,0,0,1,0,1,0,0,0,0,1,0,0,0,1,1,1,0,1,0,1,0,0,0))) gph - xyplot (mortality ~ attend|type, pch=16, data = xdat, aspect=2:1,layout=c(3,1)) one4all - layer(avy - median(xdat$mortality), panel.segments(0.1, avy, 0.3, avy, col=red,lwd=4), panel.segments(0.7, avy, 1, avy, col=red,lwd=4)) medbytype - layer(avy - median(y), panel.segments(0.1, avy, 0.3, avy, col=red,lwd=4), panel.segments(0.7, avy, 1, avy, col=red,lwd=4)) interact - layer(panel.average(x, y, fun=median, col='red', lwd=4)) Compare gph + one4all (shows overall median lines in all 3 panels) gph + medbytype (shows separate median lines for the separate panels) gph+interact (gives a form of interaction plot) NB x (if its values are used) and y are local to the individual panel NB also that layer() accepts a data argument. The following is an alternative way to calculate one4all: one4all - layer(data=xdat, avy - median(mortality), panel.segments(0.1, avy, 0.3, avy, col=red,lwd=4), panel.segments(0.7, avy, 1, avy, col=red, lwd=4)) John Maindonald. This can be simplified by using the layering abilities that Felix Andrews made available in latticeExtra. These are too little known. These pretty much make it unnecessary to resort to trellis.focus(), at least in such cases as this. These layering abilities are too little known: library(latticeExtra) x11(height=8,width=11) xdat = data.frame(mortality =c(5, 8,7,5,8,10,11,6,4,5,20,25,27,30,35,32,28,21,20,34,11,15,18,12,15,12,10,15,19,20), type= c(1, 1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2,3,3,3,3,3,3,3,3,3,3), attend = c(1, 0,1,1,1,1,0,0,0,1,0,1,0,0,0,0,1,0,0,0,1,1,1,0,1,0,1,0,0,0)) gph - xyplot ( mortality ~ attend|type, panel=function(x,y) { panel.grid(lty=5) panel.xyplot(x,y,pch=16,jitter.x=TRUE,col=1) for(i in 1:3) { panel.segments(x0=c(0.7, 1.7), x1=c(1.3, 2.3), y0=with(xdat, c(tapply(mortality[type==i], attend[type==i], median))), y1=with(xdat, c(tapply(mortality[type==i], attend[type==i], median))),col=red,lwd=4) } }, data = xdat, aspect=2:1,layout=c(3,1)) ## Now create a layer object that will add the further segments. addlayer - layer(panel.segments(x0=1,y0=20,x1=1.2, y1=20), panel.segments(x0=2,y0=30,x1=2.2, y1=30)) gph+addlayer The code that produces the object gph would also be simpler and easier to follow if some relevant part was separated out into a separate layer. See also my notices on layering of lattice objects at: http://www.maths.anu.edu.au/%7Ejohnm/r-book/add-graphics.html John Maindonald email: john.maindon...@anu.edu.au phone : +61 2 (6125)3473fax : +61 2(6125)5549 Centre for Mathematics Its Applications, Room 1194, John Dedman Mathematical Sciences Building (Building 27) Australian National University, Canberra ACT 0200. http://www.maths.anu.edu.au/~johnm On 08/04/2012, at 8:00 PM, r-help-requ...@r-project.org wrote: From: David Winsemius dwinsem...@comcast.net Subject: Re: [R] Drawing a line in xyplot Date: 8 April 2012 2:17:22 PM AEST To: wcheckle wchec...@jhsph.edu Cc: r-help@r-project.org On Apr 7, 2012, at 10:29 PM, wcheckle wrote: Thank you David, the bwplot option does what I need: x11(height=8,width=11) bwplot ( mortality~ attend|type, pch=95,cex=5,col=2, par.settings=list( box.rectangle = list(col = transparent), box.umbrella = list(col = transparent), plot.symbol = list(col = transparent) ), panel=function(x,y,...){ panel.grid(lty=5) panel.xyplot(x,y,pch=16,jitter.x=TRUE,col=1) panel.bwplot(x,y,...) }, data = x,aspect=2:1,layout=c(3,1)) However, I am interested in also learning how to do it in xyplot as well. I wasnt able to follow the last two set of instructions (condition on packet.number and loop over segments), wondering if I can ask for your help for the correct code (my attempt resulted in all three mean lines within each panel): x11(height=8,width=11) xyplot ( mortality ~ attend|type, panel=function(x,y) { panel.grid(lty=5) panel.xyplot(x,y,pch=16,jitter.x=TRUE,col=1) for(i in 1:3) { panel.segments(x0=c(0.7, 1.7), x1=c(1.3, 2.3), y0=c(tapply(x$mortality[x$type==i], x$attend[x$type==i], median)), y1=c(tapply(x$mortality[x$type==i], x$attend[x$type==i], median)),col=red,lwd=4) } }, data = x,aspect=2:1,layout=c(3,1)) thank you again. I also found your info on data.frame useful. I haven't figured out how to do it with the interaction of 'attend' and ''type' from inside xyplot. I'm thinking I might be able to succeed using trellis.focus() to address separate columns within a
[R] Stepwise procedure with force.in command
Dear R-helpers, I am trying to do a stepwise procedure in which I want to force some variables in the model. I have searched around and it seems that only leaps package allows to force the variable in the stepwise procedure. I use the leaps package and use the regsubsets(lm1, force.in = 1, data) to force 1 variable in the model. However, the force.in command only allow me to force 1 variable in. I want to force several variables in the procedure. How could I force them in? Thanks a lot, Hien __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 get the confidence interval of area under the time dependent roc curve
On 04/07/2012 05:00 AM, r-help-requ...@r-project.org wrote: It is possible to calculate the c-index for time dependent outcomes (such as disease) using the survivalROC package in R. My question is : is it possible to produce a p-value for the c-index that is calculated (at a specific point in time)?How to get the confidence interval of area under then time dependent roc curve (or by bootstrap)? The current version of survival computes the c-index and it's standard error for time-dependent covariates. Use summary(cox model fit) or the survConcordance function. Terry T. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Order sapply
Good Afternoon, I get it, my question was really how the data were organized, I thought I was doing something wrong. Thanks -- View this message in context: http://r.789695.n4.nabble.com/Order-sapply-tp4537496p4542518.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] Question on harmonic (Fourier) analysis of sinusoidal time series
Hello, I will try to explain the problem, sorry if it will be a little long... I'm using R to analyze results of cyclic mechanical testing, like this: - apply quasi-sinusoidal load - measure quasi-sinusoidal vertical and horizontal deformations (quasi-sinusoidal load means that load should be sinusoidal, but testing machine puts in some noise...) I enclose a sample of data at the end of this message. I used harmonic regression to fit the data: force = constant + [1st order Fourier harmonic] + [2nd order Fourier harmonic] + higher order harmonics... + error = = v0 + [a1 sin (2*pi*f*t) + b1 cos (2*pi*f*t) ] + [a2 sin (2*2*pi*f*t) + b2 cos (2*2*pi*f*t)] + ... + error where f is the frequency of the sinusoidal load applied during the test. What I'm expecting to find is: - v0 (average value) - a1, b1 high values of coefficients of the first fundamental harmonic that the machine is applying - an, bn low values of coefficients of the higher harmonics - some noise ~ NID I used lm() to fit the linear model using the code suggested in Introductory Time Series with R by Paul S.P. Cowpertwait (2006 SPRINGER SCIENCE+BUSINESS MEDIA, LLC.). Started with 1st harmonic component and then increased up to 6. Here's the results (6 harmonic components, but results are similar with less harmonics): Call: lm(formula = force.fit ~ s + c) Residuals: Min 1Q Median 3QMax -0.0147651 -0.0052199 -0.0005713 0.0052815 0.0185080 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) -1.3264070 0.0002227 -5957.00 2e-16 *** s1 -1.3083551 0.0003149 -4154.91 2e-16 *** s2 -0.0538174 0.0003149 -170.91 2e-16 *** s3 -0.0346178 0.0003149 -109.94 2e-16 *** s4 -0.0070427 0.0003149 -22.36 2e-16 *** s5 -0.0069028 0.0003149 -21.92 2e-16 *** s6 0.0155151 0.000314949.27 2e-16 *** c1 0.3746976 0.0003149 1189.92 2e-16 *** c2 -0.0149771 0.0003149 -47.56 2e-16 *** c3 0.0272273 0.000314986.47 2e-16 *** c4 -0.0282915 0.0003149 -89.84 2e-16 *** c5 -0.0084926 0.0003149 -26.97 2e-16 *** c6 -0.0093216 0.0003149 -29.60 2e-16 *** Fundamental harmonic is there: a1 = -1.3083551, b1 = 0.3746976 and, as expected, higher order harmonics have much smaller coefficients. BUT residuals are PERIODIC: http://r.789695.n4.nabble.com/file/n4542389/residualsplot.jpg Since my error doesn't look like noise ( is strongly autocorrelated) did I make any mistakes? Thank you for reading down to this line... I hope you have suggestions! Thanks Andrea ### Here's a sample of data: data[1:10,] time forcevdef hdef 1 0. 0.000 0.000 0.000 2 0.0025 -0.007 -0.229 0.000 3 0.0050 -0.166 -1.829 0.453 4 0.0075 -0.544 -7.086 2.038 5 0.0100 -0.957 -13.714 3.849 6 0.0125 -1.304 -19.886 5.660 7 0.0150 -1.550 -25.143 7.245 8 0.0175 -1.719 -29.714 8.377 9 0.0200 -1.838 -33.372 9.510 10 0.0225 -1.933 -36.343 10.642 -- View this message in context: http://r.789695.n4.nabble.com/Question-on-harmonic-Fourier-analysis-of-sinusoidal-time-series-tp4542389p4542389.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] Creating Better Table in R
Could anyone please direct me on how to make a nicer table in R? THANKS FOR ALL THE HELP! I would like to make a table with the following in it: estimate, t value, significance, beta, standard errors, adjusted r squared, and residual standard error (3 decimal points if possible, but I can do it by hand). Also I'd love to include the source of my information on the bottom of the table if possible. file : http://r.789695.n4.nabble.com/file/n4542349/datpat.csv datpat.csv this is the code I am using mod.1-lm(Patents~FHouse, data=datpat) summary(mod.1) xtable(summary(mod.1)) *table I'm creating* http://r.789695.n4.nabble.com/file/n4542349/chart.jpg *Detailed code* *Summary Stats:* Residuals: Min 1Q Median 3Q Max -1.7540 -0.8833 -0.5123 0.1183 11.9858 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) 2.617760.52866 4.952 1.34e-06 *** FHouse -0.187920.08489 -2.214 0.0277 * --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 1.77 on 254 degrees of freedom Multiple R-squared: 0.01893,Adjusted R-squared: 0.01507 F-statistic: 4.901 on 1 and 254 DF, p-value: 0.02773 *X Table Code:* \begin{table}[ht] \begin{center} \begin{tabular}{r} \hline Estimate Std. Error t value Pr($$$|$t$|$) \\ \hline (Intercept) 2.6178 0.5287 4.95 0. \\ FHouse -0.1879 0.0849 -2.21 0.0277 \\ \hline \end{tabular} \end{center} \end{table} -- View this message in context: http://r.789695.n4.nabble.com/Creating-Better-Table-in-R-tp4542349p4542349.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] Shewhart Control Charts for Time Series data (qcc package)
Hello everybody! I am trying to apply qcc function to time series data. There is an example from documentation: ### library(qcc) data(pistonrings) attach(pistonrings) diameter - qcc.groups(diameter, sample) qcc(diameter[1:25,], type=xbar) detach(pistonrings) ### So, we have 5 iteams at each row - sample size =5 But what should I do if I have just a one dimensional time series data? ### out.m - aggregate(pistonrings$diameter, list(sample = pistonrings$sample), FUN=mean) qcc(out.m$x, type=xbar) - Error: group sizes must be larger than one ### Probably in this case samples are a consecutive subsets of the data but they are not grouped as in an example. I would like to apply qcc() like rollapply(width = sample_size, by = step) in such a way that 'width' would stands for a samples taken from an input data. Thank you in advance ! -- Kind regards, Alexandr [[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] i am not getting time series prediction results?
i have attached file, http://r.789695.n4.nabble.com/file/n4542543/1BHP_A.txt 1BHP_A.txt i just want to validate predicting results, i give less 10 digits for prediction compare results it with previous input data(input file) i-read.table(file.choose()) #attached file taking i-i$V1 length(i) [1] 44 j-i[1:34] pre-predict(ar(i),n.ahead=10) ts.plot(ts(j),pre$pred,col=c(1,5)) pre$pred Time Series: Start = 45 End = 54 Frequency = 1 [1] -0.01136364 -0.01136364 -0.01136364 -0.01136364 -0.01136364 [6] -0.01136364 -0.01136364 -0.01136364 -0.01136364 -0.01136364 prediction gives similar values,i.e. straight line. why this so? -- View this message in context: http://r.789695.n4.nabble.com/i-am-not-getting-time-series-prediction-results-tp4542543p4542543.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How do Sweave users collaborate with Word users?
You might want to consider SWord, which provides similar facilities for the Word and R user. Word-oriented co-authors can modify the Word part of the document without impacting the R part of the document. SWord is by Thomas Baier tho...@statconn.com, author of the statconnDCOM interface that is underneath RExcel. See rcom.univie.ac.at for information and download and to sign up on the rcom email list. Rich On Sat, Apr 7, 2012 at 4:54 PM, Alexander Shenkin ashen...@ufl.edu wrote: Hello All, I'm getting my workflow switched over to Sweave, which is very cool. However, I collaborate with folks (as many of you must as well) who use Word to Track Changes amongst a group while crafting a paper. In the simplest case, there will just be two people (one Sweave user and one Word user) editing a paper. I'm wondering, how do Sweave users go about this? I could convert a sweave file to a .docx easily enough via an intermediary pdf, rtf, html or otherwise. However, once the file has been marked up with changes, the challenge is to migrate those (accepted) changes back to the sweave document. Perhaps the most straightforward way is to manually back-propagate changes, but I imagine that could be a painstaking process. Ideally, I imagine a tool that puts invisible tags in the word document when it is originally produced from Sweave, and is then able to propagate changes back to that sweave file after markup. I'd be pleasantly surprised if such a tool existed. Perhaps there are other ways of making this work. Any thoughts are kindly appreciated! Thanks, Allie __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.htmlhttp://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 using R 2.14.2
Sir, I am a student in biostat and bioinformatics. I am interested to use R for my research work related to codon usage analysis.This software was used in a publication entitled Online synonymous codon usage analyses with the ade4 and seqinR packages and a weblink http://pbil.univ-lyon1.fr/datasets/charif04/ was given for readers to use the program. Now I would like to use the webpage for codon usage analysis for my selected genome.However, the program is not working for any other genome other than its default data set.I can filter the genome to keep only sequences more than 300bp in a text file in fasta formata and I would like to run the program/weblink for this file. Could you please let me know, how to use the program for my dataset on the weblink / a edited R script for running the program on R 2.14.2 for this set of sequences. Sincerely Kinshuk Chandra Nayak Institute of Life Sciences Department of Biotechnology,Govt India Bhubaneswar-751023, India [[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] Error in integrate(int.fn, lower = 0, upper = Inf) : evaluation of function gave a result of wrong length
On Apr 8, 2012, at 11:48 PM, Guaramy wrote: Sorry but i didn't understand. The path to understanding is through study and application of the advice in the Posting Guide. -- David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How do Sweave users collaborate with Word users?
If you're considering SWord, you should remember that the licence is not the normal R licence and in commercial use will require a commercial licence. While some academic disciplines use Word etc, the issue may be more common outside academia. For those of us where such requirements involve a procurement process, the need to purchase something when other users (and the administration) are happy with their Word/Excel solutions may be an insurmountable barrier. Good luck Paul Bivand Centre for Economic and Social Inclusion (a non-profit organisation) London On 9 April 2012 13:23, Richard M. Heiberger r...@temple.edu wrote: You might want to consider SWord, which provides similar facilities for the Word and R user. Word-oriented co-authors can modify the Word part of the document without impacting the R part of the document. SWord is by Thomas Baier tho...@statconn.com, author of the statconnDCOM interface that is underneath RExcel. See rcom.univie.ac.at for information and download and to sign up on the rcom email list. Rich On Sat, Apr 7, 2012 at 4:54 PM, Alexander Shenkin ashen...@ufl.edu wrote: Hello All, I'm getting my workflow switched over to Sweave, which is very cool. However, I collaborate with folks (as many of you must as well) who use Word to Track Changes amongst a group while crafting a paper. In the simplest case, there will just be two people (one Sweave user and one Word user) editing a paper. I'm wondering, how do Sweave users go about this? I could convert a sweave file to a .docx easily enough via an intermediary pdf, rtf, html or otherwise. However, once the file has been marked up with changes, the challenge is to migrate those (accepted) changes back to the sweave document. Perhaps the most straightforward way is to manually back-propagate changes, but I imagine that could be a painstaking process. Ideally, I imagine a tool that puts invisible tags in the word document when it is originally produced from Sweave, and is then able to propagate changes back to that sweave file after markup. I'd be pleasantly surprised if such a tool existed. Perhaps there are other ways of making this work. Any thoughts are kindly appreciated! Thanks, Allie __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.htmlhttp://www.r-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Stepwise procedure with force.in command
On 2012-04-09 00:44, Hien Nguyen wrote: Dear R-helpers, I am trying to do a stepwise procedure in which I want to force some variables in the model. I have searched around and it seems that only leaps package allows to force the variable in the stepwise procedure. I use the leaps package and use the regsubsets(lm1, force.in = 1, data) to force 1 variable in the model. However, the force.in command only allow me to force 1 variable in. I want to force several variables in the procedure. How could I force them in? Well, I don't know where you've searched, but the help page for step() or for stepAIC() in pkg MASS tells you that there's a scope= argument which should do what you want. Peter Ehlers Thanks a lot, Hien __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] R problem nls
Hi, I will try to explain what it is I need to do, how far I am in doing it yet and where my problem is: I have a lot of x,y values I need to fit a non linear function through. Subsequently, I need to find the intersection point of this fitted curve with y=1.01 The problem is I have a lot of values so I want to be able to do it all at once. I already imported my excel file of the points I have to use. The things I already have are following (not all of the data is visible because otherwise the file would be too long for this email, but i just showed part of my data so you could better understand my problem): __ X1242 X1242.5 X1243 X1243.5 X1244 X1244.5 X1245 X1245.5 X1246 X1246.5 X1247 X1247.5 X1248 X1248.5 X1249 = names of each data set ( corresponds to wavelengths, ranges from 400 to 2500 with 0.5 steps = a lot of points!) 18.14 0.9860316 0.9860272 0.9860203 0.9860121 0.9860044 0.9859994 0.9859971 0.9859976 0.985 0.9860035 0.9860069 0.9860103 0.9860128 0.9860151 0.9860178 15.8 0.9857134 0.9857106 0.9857063 0.9857011 0.9856958 0.9856917 0.9856887 0.9856874 0.9856880 0.9856893 0.9856906 0.9856919 0.9856922 0.9856916 0.9856912 13.77 0.9930109 0.9930015 0.9929921 0.9929833 0.9929765 0.9929714 0.9929674 0.9929637 0.9929603 0.9929569 0.9929533 0.9929501 0.9929469 0.9929440 0.9929416 9.03 0.9875374 0.9875321 0.9875242 0.9875140 0.9875024 0.9874921 0.9874840 0.9874793 0.9874780 0.9874802 0.9874835 0.9874869 0.9874877 0.9874857 0.9874801 6.14 0.9900554 0.9900465 0.9900376 0.9900286 0.9900204 0.9900122 0.9900032 0.9899924 0.9899802 0.9899669 0.9899544 0.9899445 0.9899372 0.9899333 0.9899317 4.27 1.0050327 1.0050242 1.0050175 1.0050129 1.0050107 1.0050101 1.0050094 1.0050070 1.0050025 1.0049959 1.0049885 1.0049812 1.0049746 1.0049691 1.0049647 2.77 0.9892697 0.9892585 0.9892454 0.9892311 0.9892164 0.9892030 0.9891906 0.9891796 0.9891704 0.9891624 0.9891550 0.9891480 0.9891401 0.9891320 0.9891235 1.52 0.9979284 0.9979430 0.9979548 0.9979644 0.9979739 0.9979850 0.9979984 0.9980137 0.9980312 0.9980498 0.9980691 0.9980897 0.9981105 0.9981323 0.9981542 These are my x-values and y-values X1249.5 X1250 X1250.5 X1251 X1251.5 X1252 X1252.5 X1253 X1253.5 X1254 X1254.5 X1255 X1255.5 X1256 X1256.5 18.14 0.9860214 0.9860261 0.9860320 0.9860377 0.9860425 0.9860456 0.9860462 0.9860449 0.9860433 0.9860422 0.9860417 0.9860428 0.9860444 0.9860456 0.9860456 15.8 0.9856911 0.9856918 0.9856934 0.9856958 0.9856984 0.9857014 0.9857040 0.9857069 0.9857099 0.9857132 0.9857153 0.9857156 0.9857132 0.9857080 0.9857016 13.77 0.9929406 0.9929405 0.9929425 0.9929445 0.9929464 0.9929476 0.9929476 0.9929455 0.9929431 0.9929409 0.9929377 0.9929356 0.9929331 0.9929304 0.9929279 9.03 0.9874721 0.9874641 0.9874578 0.9874542 0.9874529 0.9874541 0.9874556 0.9874563 0.9874565 0.9874551 0.9874527 0.9874498 0.9874461 0.9874415 0.9874371 6.14 0.9899318 0.9899319 0.9899304 0.9899263 0.9899192 0.9899095 0.9898986 0.9898873 0.9898777 0.9898703 0.9898641 0.9898591 0.9898546 0.9898495 0.9898439 4.27 1.0049605 1.0049564 1.0049528 1.0049495 1.0049466 1.0049435 1.0049404 1.0049357 1.0049298 1.0049230 1.0049155 1.0049093 1.0049049 1.0049019 1.0049002 2.77 0.9891158 0.9891100 0.9891058 0.9891036 0.9891013 0.9890986 0.9890943 0.9890873 0.9890789 0.9890691 0.9890583 0.9890485 0.9890395 0.9890323 0.9890266 1.52 0.9981760 0.9981970 0.9982176 0.9982372 0.9982565 0.9982753 0.9982935 0.9983102 0.9983259 0.9983408 0.9983533 0.9983647 0.9983749 0.9983841 0.9983929 X1257 X1257.5 X1258 X1258.5 X1259 X1259.5 X1260 X1260.5 X1261 X1261.5 X1262 X1262.5 X1263 X1263.5 X1264 18.14 0.9860445 0.9860437 0.9860438 0.9860467 0.9860527 0.9860613 0.9860705 0.9860782 0.9860830 0.9860844 0.9860832 0.9860806 0.9860774 0.9860746 0.9860716 15.8 0.9856955 0.9856925 0.9856930 0.9856967 0.9857032 0.9857098 0.9857162 0.9857213 0.9857248 0.9857268 0.9857283 0.9857298 0.9857314 0.9857340 0.9857374 13.77 0.9929253 0.9929242 0.9929234 0.9929239 0.9929254 0.9929276 0.9929302 0.9929321 0.9929329 0.9929328 0.9929315 0.9929294 0.9929275 0.9929257 0.9929243 9.03 0.9874330 0.9874313 0.9874312 0.9874335 0.9874375 0.9874418 0.9874468 0.9874510 0.9874547 0.9874577 0.9874602 0.9874623 0.9874640 0.9874659 0.9874680 6.14 0.9898392 0.9898366 0.9898360 0.9898381 0.9898417 0.9898455 0.9898484 0.9898489 0.9898469 0.9898432 0.9898392 0.9898363 0.9898354 0.9898371 0.9898402 4.27 1.0048982 1.0048966 1.0048942 1.0048923 1.0048914 1.0048917 1.0048929 1.0048939 1.0048933 1.0048910 1.0048865 1.0048803 1.0048742 1.0048693 1.0048658 2.77 0.9890227 0.9890210 0.9890197 0.9890184 0.9890164 0.9890143 0.9890120 0.9890102 0.9890088 0.9890087 0.9890089 0.9890086 0.9890073 0.9890056 0.9890024 1.52 0.9984034 0.9984168 0.9984328 0.9984514 0.9984707 0.9984899 0.9985069 0.9985223 0.9985357 0.9985479 0.9985593 0.9985694 0.9985778
Re: [R] Error in integrate(int.fn, lower = 0, upper = Inf) : evaluation of function gave a result of wrong length
I read it believe me . The reason that a i post this is because i am making a thesis and a i am having this problem for over 2 weeks. I can´t solve it and its causing me real problems. Sorry for any inconvenience Thanks -- View this message in context: http://r.789695.n4.nabble.com/Error-in-integrate-int-fn-lower-0-upper-Inf-evaluation-of-function-gave-a-result-of-wrong-length-tp4541250p4542677.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Creating Better Table in R
Hi, Yes, please check package tables. Regards, Carlos Ortega www.qualityexcellence.es 2012/4/9 bobo bleza...@gmail.com Could anyone please direct me on how to make a nicer table in R? THANKS FOR ALL THE HELP! I would like to make a table with the following in it: estimate, t value, significance, beta, standard errors, adjusted r squared, and residual standard error (3 decimal points if possible, but I can do it by hand). Also I'd love to include the source of my information on the bottom of the table if possible. file : http://r.789695.n4.nabble.com/file/n4542349/datpat.csv datpat.csv this is the code I am using mod.1-lm(Patents~FHouse, data=datpat) summary(mod.1) xtable(summary(mod.1)) *table I'm creating* http://r.789695.n4.nabble.com/file/n4542349/chart.jpg *Detailed code* *Summary Stats:* Residuals: Min 1Q Median 3Q Max -1.7540 -0.8833 -0.5123 0.1183 11.9858 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) 2.617760.52866 4.952 1.34e-06 *** FHouse -0.187920.08489 -2.214 0.0277 * --- Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1 Residual standard error: 1.77 on 254 degrees of freedom Multiple R-squared: 0.01893,Adjusted R-squared: 0.01507 F-statistic: 4.901 on 1 and 254 DF, p-value: 0.02773 *X Table Code:* \begin{table}[ht] \begin{center} \begin{tabular}{r} \hline Estimate Std. Error t value Pr($$$|$t$|$) \\ \hline (Intercept) 2.6178 0.5287 4.95 0. \\ FHouse -0.1879 0.0849 -2.21 0.0277 \\ \hline \end{tabular} \end{center} \end{table} -- View this message in context: http://r.789695.n4.nabble.com/Creating-Better-Table-in-R-tp4542349p4542349.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. -- Saludos, Carlos Ortega www.qualityexcellence.es [[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 problem nls
Your formula could be simplified to y ~ 1 + a * exp(-b * x) Solve this equation for x x = ln[(y - 1)/a]/b Use this equation to find the intersection point at a given value of y. For example, when y = 1.01 x = ln(0.01/a)/b Jean Karen Vandepoel wrote on 04/09/2012 07:33:19 AM: Hi, I will try to explain what it is I need to do, how far I am in doing it yet and where my problem is: I have a lot of x,y values I need to fit a non linear function through. Subsequently, I need to find the intersection point of this fitted curve with y=1.01 The problem is I have a lot of values so I want to be able to do it all at once. I already imported my excel file of the points I have to use. The things I already have are following (not all of the data is visible because otherwise the file would be too long for this email, but i just showed part of my data so you could better understand my problem): __ X1242 X1242.5 X1243 X1243.5 X1244 X1244.5 X1245 X1245.5 X1246 X1246.5 X1247 X1247.5 X1248 X1248.5 X1249 = names of each data set ( corresponds to wavelengths, ranges from 400 to 2500 with 0.5 steps = a lot of points!) 18.14 0.9860316 0.9860272 0.9860203 0.9860121 0.9860044 0.9859994 0.9859971 0.9859976 0.985 0.9860035 0.9860069 0.9860103 0.9860128 0.9860151 0.9860178 15.8 0.9857134 0.9857106 0.9857063 0.9857011 0.9856958 0.9856917 0.9856887 0.9856874 0.9856880 0.9856893 0.9856906 0.9856919 0.9856922 0.9856916 0.9856912 13.77 0.9930109 0.9930015 0.9929921 0.9929833 0.9929765 0.9929714 0.9929674 0.9929637 0.9929603 0.9929569 0.9929533 0.9929501 0.9929469 0.9929440 0.9929416 9.03 0.9875374 0.9875321 0.9875242 0.9875140 0.9875024 0.9874921 0.9874840 0.9874793 0.9874780 0.9874802 0.9874835 0.9874869 0.9874877 0.9874857 0.9874801 6.14 0.9900554 0.9900465 0.9900376 0.9900286 0.9900204 0.9900122 0.9900032 0.9899924 0.9899802 0.9899669 0.9899544 0.9899445 0.9899372 0.9899333 0.9899317 4.27 1.0050327 1.0050242 1.0050175 1.0050129 1.0050107 1.0050101 1.0050094 1.0050070 1.0050025 1.0049959 1.0049885 1.0049812 1.0049746 1.0049691 1.0049647 2.77 0.9892697 0.9892585 0.9892454 0.9892311 0.9892164 0.9892030 0.9891906 0.9891796 0.9891704 0.9891624 0.9891550 0.9891480 0.9891401 0.9891320 0.9891235 1.52 0.9979284 0.9979430 0.9979548 0.9979644 0.9979739 0.9979850 0.9979984 0.9980137 0.9980312 0.9980498 0.9980691 0.9980897 0.9981105 0.9981323 0.9981542 These are my x-values and y-values X1249.5 X1250 X1250.5 X1251 X1251.5 X1252 X1252.5 X1253 X1253.5 X1254 X1254.5 X1255 X1255.5 X1256 X1256.5 18.14 0.9860214 0.9860261 0.9860320 0.9860377 0.9860425 0.9860456 0.9860462 0.9860449 0.9860433 0.9860422 0.9860417 0.9860428 0.9860444 0.9860456 0.9860456 15.8 0.9856911 0.9856918 0.9856934 0.9856958 0.9856984 0.9857014 0.9857040 0.9857069 0.9857099 0.9857132 0.9857153 0.9857156 0.9857132 0.9857080 0.9857016 13.77 0.9929406 0.9929405 0.9929425 0.9929445 0.9929464 0.9929476 0.9929476 0.9929455 0.9929431 0.9929409 0.9929377 0.9929356 0.9929331 0.9929304 0.9929279 9.03 0.9874721 0.9874641 0.9874578 0.9874542 0.9874529 0.9874541 0.9874556 0.9874563 0.9874565 0.9874551 0.9874527 0.9874498 0.9874461 0.9874415 0.9874371 6.14 0.9899318 0.9899319 0.9899304 0.9899263 0.9899192 0.9899095 0.9898986 0.9898873 0.9898777 0.9898703 0.9898641 0.9898591 0.9898546 0.9898495 0.9898439 4.27 1.0049605 1.0049564 1.0049528 1.0049495 1.0049466 1.0049435 1.0049404 1.0049357 1.0049298 1.0049230 1.0049155 1.0049093 1.0049049 1.0049019 1.0049002 2.77 0.9891158 0.9891100 0.9891058 0.9891036 0.9891013 0.9890986 0.9890943 0.9890873 0.9890789 0.9890691 0.9890583 0.9890485 0.9890395 0.9890323 0.9890266 1.52 0.9981760 0.9981970 0.9982176 0.9982372 0.9982565 0.9982753 0.9982935 0.9983102 0.9983259 0.9983408 0.9983533 0.9983647 0.9983749 0.9983841 0.9983929 X1257 X1257.5 X1258 X1258.5 X1259 X1259.5 X1260 X1260.5 X1261 X1261.5 X1262 X1262.5 X1263 X1263.5 X1264 18.14 0.9860445 0.9860437 0.9860438 0.9860467 0.9860527 0.9860613 0.9860705 0.9860782 0.9860830 0.9860844 0.9860832 0.9860806 0.9860774 0.9860746 0.9860716 15.8 0.9856955 0.9856925 0.9856930 0.9856967 0.9857032 0.9857098 0.9857162 0.9857213 0.9857248 0.9857268 0.9857283 0.9857298 0.9857314 0.9857340 0.9857374 13.77 0.9929253 0.9929242 0.9929234 0.9929239 0.9929254 0.9929276 0.9929302 0.9929321 0.9929329 0.9929328 0.9929315 0.9929294 0.9929275 0.9929257 0.9929243 9.03 0.9874330 0.9874313 0.9874312 0.9874335 0.9874375 0.9874418 0.9874468 0.9874510 0.9874547 0.9874577 0.9874602 0.9874623 0.9874640 0.9874659 0.9874680 6.14 0.9898392 0.9898366 0.9898360 0.9898381 0.9898417 0.9898455 0.9898484 0.9898489 0.9898469 0.9898432 0.9898392 0.9898363 0.9898354 0.9898371 0.9898402 4.27 1.0048982
Re: [R] Error in integrate(int.fn, lower = 0, upper = Inf) : evaluation of function gave a result of wrong length
On 09-04-2012, at 14:54, Guaramy wrote: I read it believe me . The reason that a i post this is because i am making a thesis and a i am having this problem for over 2 weeks. I can´t solve it and its causing me real problems. You didn't give us a reproducible example. I'm assuming your original function is GNL.pdf.fn - function(x,mu,sigma,alpha,beta,rho) # missing in your original post { y - x-rho*mu cf.fn - function(s){ cplex = complex(1,0,1) temp1 = alpha*beta*exp(-sigma*s^2/2) temp2 = (alpha-cplex*s)*(beta+cplex*s) out = (temp1/temp2)^rho out } temp.fn - function(s){ Mod(cf.fn(s))*cos(Arg(cf.fn(s))-s*y) } int.fn - function(t){sapply(t,FUN=temp.fn)} te - integrate(int.fn,lower=0,upper=Inf,rel.tol=1e-10,subdivisions=100) temp3 - ifelse(te$message == OK,te$value/pi,NA) temp3 } Question: why are you using ifelse when te$message is a scalar? if( te$message == OK ) te$value/pi else NA would do what you want. David meant that you should have a look at Vectorize. And then do some experimenting. So I tried this GNL.pdf.fn(1, .2, .3, 1, 1, .6) GNL.vec - Vectorize(GNL.pdf.fn, x) GNL.vec(c(.9,1,1.1), .2, .3, 1, 1, .6) And now it's your turn for seq(-4,4,0.1). 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] Error in integrate(int.fn, lower = 0, upper = Inf) : evaluation of function gave a result of wrong length
On Apr 9, 2012, at 14:54 , Guaramy wrote: I read it believe me . The reason that a i post this is because i am making a thesis and a i am having this problem for over 2 weeks. I can´t solve it and its causing me real problems. Well, if you have a function that is not vectorized, i.e. it works for a single argument but not for a vector of arguments, and you want it to be vectorized, the tool is Vectorize(), the help page of which David pointed you to. (This essentially works around the problem by creating a new function which calls the function one vector element at the time.) The direct cause of your problem seems to be that if x is a vector, so is y, and your temp.fn returns a vector, so int.fn returns a matrix, and integrate() is unhappy with that. -- Peter Dalgaard, Professor, Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Email: pd@cbs.dk Priv: pda...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Specifying the ordering of a vector
Hi, I'm trying to create a vector (or matrix row) with a specific ordering. For example, I have the following vector: x-c(0.1,0.2,0.3,0.4,0.5,0.6) that has order order(x) [1] 1 2 3 4 5 6 I want another vector that has the same values as x, but with a different ordering. For example, I want y to have values 0.1, 0.2, etc. but in the order 1-2-5-6-3-4. The answer would be y [1] 0.1 0.2 0.5 0.6 0.3 0.4 Any help would be greatly appreciated. Thanks. -- View this message in context: http://r.789695.n4.nabble.com/Specifying-the-ordering-of-a-vector-tp4542766p4542766.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Comparing 2 means. pls help
You should look up what a t-test is. Michael On Mon, Apr 9, 2012 at 2:58 AM, ali_protocol mohammadianalimohammad...@gmail.com wrote: I am interested in the difference of 2 data: mat1= c(2.2, 2.3, 2.2,2.5) mat2= c(2.6, 2.8, 2.7,2.4) mat= mat2-mat1 I perform an action on both mat1 and mat2, and I get mat1prime and mat2prime: mat1prime= c(2.5, 2.5, 2.3,2.5) mat2prime= c(2.6, 2.8, 2.7,2.6) matprime= mat2prime-mat1prime I want to check wether mat and matprime have the same mean. How should I do this in R? -- View this message in context: http://r.789695.n4.nabble.com/Comparing-2-means-pls-help-tp4542237p4542237.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Avoid loop with the integrate function
Hi Berend, I now understand what you were explaining. It works great. Thanks a lot for your time and help. I really appreciate it. Navin Goyal On Mon, Apr 9, 2012 at 1:26 AM, Berend Hasselman b...@xs4all.nl wrote: On 09-04-2012, at 01:26, Navin Goyal wrote: Hi, I am not sure I follow the scalar part of the suggestion. Could you please elaborate it a bit more? Each row has a different value of score for each subject at each time point. That last remark has nothing to do with it. The point is that in your function func1 the returned value doesn't depend on the argument t. So func1 is returning a vector of length==length(t) of identical values except for the first entry which is zero. Do this after the loop for comb3 # No integral since expression in func1 doesn't depend on t comb4-comb2 comb4$cmhz=0 comb4-comb4[order(comb4$ID, comb4$TIME), ] for (q in 1:length(comb4$ID)) { comb4$cmhz[q]-comb4$TIME[q]*func1(comb4$TIME[q], cov1=0.001,beta1=0.02, change=comb4$score[q], other=comb4$frac[q]) } head(comb4) all.equal(comb3$cmhz, comb4$cmhz) You don't need integrate to get what you seem to be wanting. Berend Also I simplified the example but there are other terms that change over each row for each subject. Does this mean that loops is the only way to go ? Thanks again for your help and time. Navin Goyal On Sun, Apr 8, 2012 at 4:03 AM, Berend Hasselman b...@xs4all.nl wrote: On 08-04-2012, at 08:28, Navin Goyal wrote: Dear R users, I am running a loop with the integrate function. I have pasted the code below. I am integrating a function from time=0 to the time value in every row. I have to perform this integration over thousands of rows with different parameters in each row. Could someone please suggest if there is an efficient/faster/easier way to do this by avoiding the loops ? Thank you so much for your help and time. -- Navin Goyal # dose-10 time-0:5 id-1:5 data1-expand.grid(id,time,dose) names(data1)-c(ID,TIME, DOSE) data1-data1[order(data1$ID,data1$TIME),] ed-data1[!duplicated(data1$ID) , c(ID,DOSE)] set.seed(5324123) for (k in 1:length(ed$ID)) { ed$base[k]-100*exp(rnorm(1,0,0.05)) ed$drop[k]-0.2*exp(rnorm(1,0,0.01)) ed$frac[k]-0.5*exp(rnorm(1,0,0.1)) } Why not ed$base - 100*exp(rnorm(length(ed$ID), 0, 0.05)) etc. comb1-merge(data1[, c(ID,TIME)], ed) comb2-comb1 comb2$score-comb2$base*exp(-comb2$drop*comb2$TIME) func1-function(t,cov1,beta1, change,other) { ifelse(t==0,cov1, cov1*exp(beta1*change+other)) } AFAICS, cov1, beta1, change and other are scalars. So when an item of t is == 0 then the function value is cov1 and in all other cases it is the same scalar and does not depend on t. So func1 is a step function. You could calculate the integral directly. Berend -- Navin Goyal -- Navin Goyal [[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] Specifying the ordering of a vector
On 09/04/2012 9:38 AM, crmnaw wrote: Hi, I'm trying to create a vector (or matrix row) with a specific ordering. For example, I have the following vector: x-c(0.1,0.2,0.3,0.4,0.5,0.6) that has order order(x) [1] 1 2 3 4 5 6 I want another vector that has the same values as x, but with a different ordering. For example, I want y to have values 0.1, 0.2, etc. but in the order 1-2-5-6-3-4. The answer would be y [1] 0.1 0.2 0.5 0.6 0.3 0.4 Any help would be greatly appreciated. Thanks. x[c(1,2,5,6,3,4)] will do it. Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Specifying the ordering of a vector
Hi, Thanks for your reply. In the example I gave in the original post, your code works. But for others, it doesn't and I'm not sure why it works for some cases and not for others. For example: x-c(0.04,0.07,0.20,0.35,0.55,0.70) order(x) [1] 1 2 3 4 5 6 Let's say I want y to be in the order 1-2-5-3-6-4. So I do: y-x[c(1,2,5,3,6,4)] y [1] 0.04 0.07 0.55 0.20 0.70 0.35 order(y) [1] 1 2 4 6 3 5 So, y isn't in the order I want it to be in. Any help that can be provided would be greatly appreciated. -- View this message in context: http://r.789695.n4.nabble.com/Specifying-the-ordering-of-a-vector-tp4542766p4542815.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] R problem nls
Hi, And regarding how to extend the nls algorithm to a larger dataset it is a question of indicating in the nls() functions which data.frame to use. In you example you were using a small set of x's ad y's, so if you want to use a large set put it in a new data.frame and pass it to nls(). And if you want to repeat that calculation for many different sets of x,y (each one of them corresponding to a different wavelength) you could do that in a loop. Regards, Carlos Ortega www.qualityexcellence.es 2012/4/9 Karen Vandepoel karen.vandep...@gmail.com Hi, I will try to explain what it is I need to do, how far I am in doing it yet and where my problem is: I have a lot of x,y values I need to fit a non linear function through. Subsequently, I need to find the intersection point of this fitted curve with y=1.01 The problem is I have a lot of values so I want to be able to do it all at once. I already imported my excel file of the points I have to use. The things I already have are following (not all of the data is visible because otherwise the file would be too long for this email, but i just showed part of my data so you could better understand my problem): __ X1242 X1242.5 X1243 X1243.5 X1244 X1244.5 X1245 X1245.5 X1246 X1246.5 X1247 X1247.5 X1248 X1248.5 X1249 = names of each data set ( corresponds to wavelengths, ranges from 400 to 2500 with 0.5 steps = a lot of points!) 18.14 0.9860316 0.9860272 0.9860203 0.9860121 0.9860044 0.9859994 0.9859971 0.9859976 0.985 0.9860035 0.9860069 0.9860103 0.9860128 0.9860151 0.9860178 15.8 0.9857134 0.9857106 0.9857063 0.9857011 0.9856958 0.9856917 0.9856887 0.9856874 0.9856880 0.9856893 0.9856906 0.9856919 0.9856922 0.9856916 0.9856912 13.77 0.9930109 0.9930015 0.9929921 0.9929833 0.9929765 0.9929714 0.9929674 0.9929637 0.9929603 0.9929569 0.9929533 0.9929501 0.9929469 0.9929440 0.9929416 9.03 0.9875374 0.9875321 0.9875242 0.9875140 0.9875024 0.9874921 0.9874840 0.9874793 0.9874780 0.9874802 0.9874835 0.9874869 0.9874877 0.9874857 0.9874801 6.14 0.9900554 0.9900465 0.9900376 0.9900286 0.9900204 0.9900122 0.9900032 0.9899924 0.9899802 0.9899669 0.9899544 0.9899445 0.9899372 0.9899333 0.9899317 4.27 1.0050327 1.0050242 1.0050175 1.0050129 1.0050107 1.0050101 1.0050094 1.0050070 1.0050025 1.0049959 1.0049885 1.0049812 1.0049746 1.0049691 1.0049647 2.77 0.9892697 0.9892585 0.9892454 0.9892311 0.9892164 0.9892030 0.9891906 0.9891796 0.9891704 0.9891624 0.9891550 0.9891480 0.9891401 0.9891320 0.9891235 1.52 0.9979284 0.9979430 0.9979548 0.9979644 0.9979739 0.9979850 0.9979984 0.9980137 0.9980312 0.9980498 0.9980691 0.9980897 0.9981105 0.9981323 0.9981542 These are my x-values and y-values X1249.5 X1250 X1250.5 X1251 X1251.5 X1252 X1252.5 X1253 X1253.5 X1254 X1254.5 X1255 X1255.5 X1256 X1256.5 18.14 0.9860214 0.9860261 0.9860320 0.9860377 0.9860425 0.9860456 0.9860462 0.9860449 0.9860433 0.9860422 0.9860417 0.9860428 0.9860444 0.9860456 0.9860456 15.8 0.9856911 0.9856918 0.9856934 0.9856958 0.9856984 0.9857014 0.9857040 0.9857069 0.9857099 0.9857132 0.9857153 0.9857156 0.9857132 0.9857080 0.9857016 13.77 0.9929406 0.9929405 0.9929425 0.9929445 0.9929464 0.9929476 0.9929476 0.9929455 0.9929431 0.9929409 0.9929377 0.9929356 0.9929331 0.9929304 0.9929279 9.03 0.9874721 0.9874641 0.9874578 0.9874542 0.9874529 0.9874541 0.9874556 0.9874563 0.9874565 0.9874551 0.9874527 0.9874498 0.9874461 0.9874415 0.9874371 6.14 0.9899318 0.9899319 0.9899304 0.9899263 0.9899192 0.9899095 0.9898986 0.9898873 0.9898777 0.9898703 0.9898641 0.9898591 0.9898546 0.9898495 0.9898439 4.27 1.0049605 1.0049564 1.0049528 1.0049495 1.0049466 1.0049435 1.0049404 1.0049357 1.0049298 1.0049230 1.0049155 1.0049093 1.0049049 1.0049019 1.0049002 2.77 0.9891158 0.9891100 0.9891058 0.9891036 0.9891013 0.9890986 0.9890943 0.9890873 0.9890789 0.9890691 0.9890583 0.9890485 0.9890395 0.9890323 0.9890266 1.52 0.9981760 0.9981970 0.9982176 0.9982372 0.9982565 0.9982753 0.9982935 0.9983102 0.9983259 0.9983408 0.9983533 0.9983647 0.9983749 0.9983841 0.9983929 X1257 X1257.5 X1258 X1258.5 X1259 X1259.5 X1260 X1260.5 X1261 X1261.5 X1262 X1262.5 X1263 X1263.5 X1264 18.14 0.9860445 0.9860437 0.9860438 0.9860467 0.9860527 0.9860613 0.9860705 0.9860782 0.9860830 0.9860844 0.9860832 0.9860806 0.9860774 0.9860746 0.9860716 15.8 0.9856955 0.9856925 0.9856930 0.9856967 0.9857032 0.9857098 0.9857162 0.9857213 0.9857248 0.9857268 0.9857283 0.9857298 0.9857314 0.9857340 0.9857374 13.77 0.9929253 0.9929242 0.9929234 0.9929239 0.9929254 0.9929276 0.9929302 0.9929321 0.9929329 0.9929328 0.9929315 0.9929294 0.9929275 0.9929257 0.9929243 9.03 0.9874330 0.9874313 0.9874312 0.9874335 0.9874375 0.9874418 0.9874468 0.9874510 0.9874547 0.9874577 0.9874602 0.9874623
[R] Calculating the overlapping area of ellipsoides
I'm generating several scatter3d visualizations of vegetation community distribution against a number parameters using the scatter3d function car package, The visualizations are fantastic and they provide great support for the analysis. As far as I can determine, scatter3d does not save the output of the routine as an object so I can not x - scatter3d(MADep2010 ~ Day7Max2010 + S7DryFreq2010 | as.factor(VegClass), data = Veg , axis.scale = TRUE, sphere.size=1.5, surface=FALSE, parallel=FALSE, ellipsoid=TRUE, revolutions=1, surface.col =c(red, green, blue, gold, magenta)) str(x) NULL I would like to extend the work by calculating the percent of overlapping area among the ellipsoides on a community by community basis. Can anyone point me to a package or suggest an approach to accomplishing this. Thanks Steve Steve Friedman Ph. D. Ecologist / Spatial Statistical Analyst Everglades and Dry Tortugas National Park 950 N Krome Ave (3rd Floor) Homestead, Florida 33034 steve_fried...@nps.gov Office (305) 224 - 4282 Fax (305) 224 - 4147 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Error in integrate(int.fn, lower = 0, upper = Inf) : evaluation of function gave a result of wrong length
Tanks a lot to all of you that take the time to help me. This is a really useful and helpful forum and i will try to help as much as i can Best Regards Jorge -- View this message in context: http://r.789695.n4.nabble.com/Error-in-integrate-int-fn-lower-0-upper-Inf-evaluation-of-function-gave-a-result-of-wrong-length-tp4541250p4542842.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Error in integrate(int.fn, lower = 0, upper = Inf) : evaluation of function gave a result of wrong length
On 09.04.2012 16:08, Guaramy wrote: Tanks a lot to all of you that take the time to help me. This is a really useful and helpful forum and i will try to help as much as i can Thanks for offering help to the R-help *mailing list*. Please try to use the mailing list interface (rather than Nabble) and quote the context. Also, please send your message not only to the list, but also to those who need helped (or in this case, who helped you). Best, Uwe Ligges Best Regards Jorge -- View this message in context: http://r.789695.n4.nabble.com/Error-in-integrate-int-fn-lower-0-upper-Inf-evaluation-of-function-gave-a-result-of-wrong-length-tp4541250p4542842.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Specifying the ordering of a vector
On 09/04/2012 9:58 AM, crmnaw wrote: Hi, Thanks for your reply. In the example I gave in the original post, your code works. But for others, it doesn't and I'm not sure why it works for some cases and not for others. For example: x-c(0.04,0.07,0.20,0.35,0.55,0.70) order(x) [1] 1 2 3 4 5 6 Let's say I want y to be in the order 1-2-5-3-6-4. So I do: y-x[c(1,2,5,3,6,4)] y [1] 0.04 0.07 0.55 0.20 0.70 0.35 order(y) [1] 1 2 4 6 3 5 So, y isn't in the order I want it to be in. Any help that can be provided would be greatly appreciated. Okay, I didn't understand what you meant. You should think of the result of order(y) as a permutation to apply to y so that it is in sorted order. That is, y[order(y)] will always be in sorted order. So what you want is a permutation of x that will be put into sorted order by y[c(1,2,5,3,6,4)], i.e. you want the inverse of the permutation c(1,2,5,3,6,4). Turns out that if you apply order to a permutation you get its inverse. So what you want is y - x[order(c(1,2,5,3,6,4))] Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How do Sweave users collaborate with Word users?
Thanks Rich, While it doesn't tickle me the way sweave/knitr does, SWord sounds more or less like the thing I'm looking for. However, poking around http://rcom.univie.ac.at/ and http://www.statconn.com doesn't reveal any download links. It seems as though SWord has been pulled from their lineup? Not sure - I've shot them an email to inquire. Thanks, Allie On 4/9/2012 7:23 AM, Richard M. Heiberger wrote: You might want to consider SWord, which provides similar facilities for the Word and R user. Word-oriented co-authors can modify the Word part of the document without impacting the R part of the document. SWord is by Thomas Baier tho...@statconn.com mailto:tho...@statconn.com, author of the statconnDCOM interface that is underneath RExcel. See rcom.univie.ac.at http://rcom.univie.ac.at for information and download and to sign up on the rcom email list. Rich On Sat, Apr 7, 2012 at 4:54 PM, Alexander Shenkin ashen...@ufl.edu mailto:ashen...@ufl.edu wrote: Hello All, I'm getting my workflow switched over to Sweave, which is very cool. However, I collaborate with folks (as many of you must as well) who use Word to Track Changes amongst a group while crafting a paper. In the simplest case, there will just be two people (one Sweave user and one Word user) editing a paper. I'm wondering, how do Sweave users go about this? I could convert a sweave file to a .docx easily enough via an intermediary pdf, rtf, html or otherwise. However, once the file has been marked up with changes, the challenge is to migrate those (accepted) changes back to the sweave document. Perhaps the most straightforward way is to manually back-propagate changes, but I imagine that could be a painstaking process. Ideally, I imagine a tool that puts invisible tags in the word document when it is originally produced from Sweave, and is then able to propagate changes back to that sweave file after markup. I'd be pleasantly surprised if such a tool existed. Perhaps there are other ways of making this work. Any thoughts are kindly appreciated! Thanks, Allie __ R-help@r-project.org mailto:R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html 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 do Sweave users collaborate with Word users?
Thanks for the heads up, Paul. That's good to know. I happen to be in academics, as are my collaborators, but I could imagine running into problems down the road. The good thing seems to be that it's just the users who want to interact with R who need the software. If collaborators are just touching the text, or are inserting graphics of their own, then sounds like things are still good. But, that's just my impression from the brief description I've seen online. Thanks, Allie On 4/9/2012 7:42 AM, Paul Bivand wrote: If you're considering SWord, you should remember that the licence is not the normal R licence and in commercial use will require a commercial licence. While some academic disciplines use Word etc, the issue may be more common outside academia. For those of us where such requirements involve a procurement process, the need to purchase something when other users (and the administration) are happy with their Word/Excel solutions may be an insurmountable barrier. Good luck Paul Bivand Centre for Economic and Social Inclusion (a non-profit organisation) London On 9 April 2012 13:23, Richard M. Heiberger r...@temple.edu wrote: You might want to consider SWord, which provides similar facilities for the Word and R user. Word-oriented co-authors can modify the Word part of the document without impacting the R part of the document. SWord is by Thomas Baier tho...@statconn.com, author of the statconnDCOM interface that is underneath RExcel. See rcom.univie.ac.at for information and download and to sign up on the rcom email list. Rich On Sat, Apr 7, 2012 at 4:54 PM, Alexander Shenkin ashen...@ufl.edu wrote: Hello All, I'm getting my workflow switched over to Sweave, which is very cool. However, I collaborate with folks (as many of you must as well) who use Word to Track Changes amongst a group while crafting a paper. In the simplest case, there will just be two people (one Sweave user and one Word user) editing a paper. I'm wondering, how do Sweave users go about this? I could convert a sweave file to a .docx easily enough via an intermediary pdf, rtf, html or otherwise. However, once the file has been marked up with changes, the challenge is to migrate those (accepted) changes back to the sweave document. Perhaps the most straightforward way is to manually back-propagate changes, but I imagine that could be a painstaking process. Ideally, I imagine a tool that puts invisible tags in the word document when it is originally produced from Sweave, and is then able to propagate changes back to that sweave file after markup. I'd be pleasantly surprised if such a tool existed. Perhaps there are other ways of making this work. Any thoughts are kindly appreciated! Thanks, Allie __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.htmlhttp://www.r-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] i am not getting time series prediction results?
On 09.04.2012 13:27, sagarnikam123 wrote: i have attached file, http://r.789695.n4.nabble.com/file/n4542543/1BHP_A.txt 1BHP_A.txt i just want to validate predicting results, i give less 10 digits for prediction compare results it with previous input data(input file) i-read.table(file.choose()) #attached file taking i-i$V1 length(i) [1] 44 j-i[1:34] pre-predict(ar(i),n.ahead=10) ts.plot(ts(j),pre$pred,col=c(1,5)) pre$pred Time Series: Start = 45 End = 54 Frequency = 1 [1] -0.01136364 -0.01136364 -0.01136364 -0.01136364 -0.01136364 [6] -0.01136364 -0.01136364 -0.01136364 -0.01136364 -0.01136364 prediction gives similar values,i.e. straight line. why this so? ar(i) gives: Call: ar(x = i) Order selected 0 sigma^2 estimated as 0.0 hence an AR 0 process was fitted, i.e. no ar component, just th mean is used for prediction (and is the only estimated parameter). mean(i) -0.01136364 Uwe Ligges -- View this message in context: http://r.789695.n4.nabble.com/i-am-not-getting-time-series-prediction-results-tp4542543p4542543.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Difference between spec.pgram spec.ar
On 08.04.2012 20:39, Bazman76 wrote: Hi there, Can someone explain what the difference between spec.pgram and spec.ar is? I understand that they attempt to do the same thing one using an AR estimation of the underlying series to estimate teh sensity the other using the FFT. However when applied to teh same data set they seem to be giving quite different results? http://r.789695.n4.nabble.com/file/n4541358/R_example.jpg Clearly the spec.ar() result seems to be smoothed but the results are also generally very different only really sharing the peak as the frequencies go towards zero. Can someone please explain why these two functions produce such different results on the same data set? Because they really measure different things? Why do you expect to get the same output in time as well as in frequency domain? Uwe Ligges code shown below: library(waveslim) vols=read.csv(file=C:/Users/ocuk/My Documents/Abs Vol.csv, header=TRUE, sep=,) x-ts(vols[,1]) #x ## LA(8) vol.la8- mra(x, la8, 4, modwt) names(vol.la8)- c(d1, d2, d3, d4, s4) ## plot multiresolution analysis of IBM data #par(mfcol=c(6,1), pty=m, mar=c(5-2,4,4-2,2)) par(mfcol=c(3,1), pty=m, mar=c(5-2,4,4-2,2)) plot.ts(x, axes=F, ylab=, main=(abs rtns)) #for(i in 1:5) # plot.ts(vol.la8[[i]], axes=F, ylab=names(vol.la8)[i]) #axis(side=1, at=seq(0,1600,by=100), # labels=c(0,,200,,400,,600,,800,,1000,,1200,,1400,,1600)) spectrum(vol.la8[[1]]) specx- spec.ar(x) -- View this message in context: http://r.789695.n4.nabble.com/Difference-between-spec-pgram-spec-ar-tp4541358p4541358.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Panel.abline would not show beyond a certain slope value
On 09.04.2012 06:31, Aziz, Muhammad Fayez wrote: Hi, Please find the code and data following. Problem appears in lines 37 - 42 of the code. I am using R.2.13.0 on WinXP. Works for me with a recent version of R (i.e. 2.15.0) / grid and lattice (i.e. 0.20-6). Please always try with recent versions of R and packages before asking. Best, Uwe Ligges Regards, Fayez PowerLawGraphsAblineQn.Râ library('reshape') library(lattice) library(igraph) # to use power.law.fit function library(latticeExtra) # to use panel.lmlineq in loglog xyplot File2Open = C:\\Documents and Settings\\All Users\\Documents\\Desktop\\Fayez\\RPractice\\PowerLawGraphsAblineQnData.txt DataTable = read.table(File2Open, header = TRUE, sep = \t) md- melt(DataTable, id.vars = c('Range', 'TheValue', 'TheNo')) # removed 'DegType' column 030312 md$variable- as.numeric(substr(md$variable, 9, nchar(as.character(md$variable mypanel4loglog- function(x, # x is the variable column in melted data, equals the Domain Nos y, # y is the value column in melted data, the degrees ... # Rest of the arguments ) { kfreq- table(y); # compute frquency hash table of y, the values k- 1:max(y) for (i in k) { ichar- as.character(i) # convert to match the names(freq), the character-based hash key of freq, which is the value if (!(ichar %in% names(kfreq))) kfreq[ichar]- 0 } sortedkeys- as.character(k) kfreq- kfreq[sortedkeys] pk- kfreq / length(y) logk- log(k) logpk- log(pk) logpk[logpk == -Inf] = # remove the -Inf or log(p(k)) = 0 values for lm function, NULL is 0-length so use instead that has length of one null character logpk- as.numeric(logpk) # converts all values to character, lm needs numeric print(rbind(logk, logpk))#write.table(rbind(k, kfreq, pk, logk, logpk), paste(FilePath, \\data, sep=), sep = \t, append = TRUE) panel.xyplot(col=blue, logk, logpk, type = c('p', 'r')) lm.r = lm(logpk ~ logk) panel.abline(coef=c(-4.847634, -1.037480)) # --- This gets drawn panel.abline(coef=c(-4.847634, -1.037481)) # --- This doesn't get drawn print(coef(lm.r)) # -4.847634 -1.349699 ; -3.377894 -1.498693 } # end mypanel4loglog myprepanel4loglog- function(x, # x is the variable column in melted data, equals the Network's ages y, # y is the value column in melted data, the degrees ... # Rest of the arguments ) { FreqTable- as.data.frame(table(y)) FreqsVector- sort(FreqTable$Freq) Min- FreqsVector[1] # first element - the lowest value frequency #print(c(Max2, length(y))) list(ylim = c(log(Min / length(y)), 0), xlim = c(0, log(max(y # log(p(k)) is always -ve as p(k) is decimal, so max(log(p(k)) is 0 } # end myprepanel4loglog print(xyplot(value ~ variable | Range, data = md, xlab = log(k); Panel = Range, ylab = log(p(k)), main = log(k) vs. log(p(k)), groups = Range, pch = 20, # dots instead of circles panel = mypanel4loglog, prepanel = myprepanel4loglog, # to set the scale of k and pk scales = list(x = (relation = free), y = (relation = free)), # necessary to make x-axis in each panel adjustable according to k )) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Specifying the ordering of a vector
Thank you very much! This worked. -- View this message in context: http://r.789695.n4.nabble.com/Specifying-the-ordering-of-a-vector-tp4542766p4542907.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Creating a dummy variable
On 09.04.2012 07:16, arunkumar wrote: HI I want to create a dummy variable for a dataset which has various dates Eg I've one year data. The user can choose any date range startdate : 2012-01-01 enddate: 2012-02-01 dummy - date = startdate date = enddate Uwe Ligges startdate: 2012-06-01 enddate=2012-07-01 the columns of dataset= X ,Y ,date Thanks and re - Thanks in Advance Arun -- View this message in context: http://r.789695.n4.nabble.com/Creating-a-dummy-variable-tp4542150p4542150.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Difference between spec.pgram spec.ar
OK so I neeed to understan better what it it they are trying to measure. I understood (incorrectly it seems) that they were simply different methods to get the same result? -- View this message in context: http://r.789695.n4.nabble.com/Difference-between-spec-pgram-spec-ar-tp4541358p4542985.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Difference between spec.pgram spec.ar
On 09.04.2012 17:01, Bazman76 wrote: OK so I neeed to understan better what it it they are trying to measure. I understood (incorrectly it seems) that they were simply different methods to get the same result? Yes. Also note this is a mailing list and you are lucky I was able to remember the context. Please always quote the prior thread and do not only send to the list. Uwe Ligges -- View this message in context: http://r.789695.n4.nabble.com/Difference-between-spec-pgram-spec-ar-tp4541358p4542985.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to find best parameter values using deSolve n optim() ?
On 4/5/2012 6:32 PM, mhimanshu wrote: Hi Thomas, Thank you so much for your suggestion. I tried your code and it is working fine. Now when I change the values of Y in yobs I am getting so many warnings. say, yobs- data.frame( time = 0:7, Y = c(0.00, 3.40, 4.60 ,5.80, 5.80, 6.00, 6.00 ), Z = c(0.1, 0.11, 0.119, 0.128, 0.136, 0.145, 0.153, 0.16) ) This data set is wrong because time and Z have 8 elements but Y has only 7. Let's assume the missing Y is 6.0. So when i fit the model with the same code that you have written, i got the following warnings: DLSODA- Warning..Internal T (=R1) and H (=R2) are such that in the machine, T + H = T on the next step (H = step size). Solver will continue anyway. In above, R1 = 0.1484502806322D+01 R2 = 0.2264549048113D-16 and I have got so many such warnings. This means that the step adjustment algorithm was unable to determine an appropriate step size (h), possibly because of wrongly negative parameter values. Can you explain me why this is happening?? and Secondly, I dont understand why i am getting parameters values in negatives after fitting?? You can restrict the parameter ranges by using the optional arguments upper and lower, see ?modFit for details. Can you please help me out in this... :) The reason for this can be manifold, e.g. wrong model specification, unrealistic data, inappropriate accuracy settings, inadequate start parameters or just unidentifiable parameters, e.g. if some of them depend strongly on each other. The example below showed clearly that the parameters are strongly correlated (0.998 or even 1.000 !!!). This means that not all parameters can be identified simultaneously because of high interdependency (see my comment in my first post), so you should try to reduce the number of fitted parameters. Note however, that Ymax needs to be fitted (or adjusted) as well. Fitting nonlinear models can be an art and requires mathematical understanding of the applied techniques (differential equations, identifiability analysis, numerical optimization), a good understanding of the properties of your (differential equation) model -- and some feeling and experience. Thomas ##- library(deSolve) library(FME) ## function derivs derivs - function(time, y, pars) { with (as.list(c(y, pars)), { dy = a1 * Y * (1 - Y/Ymax) - a2 * (Y + 0.001) dz = a3 * Y - a4 * Z; return(list(c(dy, dz))) }) } ## parameters pars - c(a1 = 0.9, a2 = 0.7, a3 = 0.06, a4 = 0.02, Ymax=0.8) #Ymax - c(0.8) ## initial values y - c(Y = 0.2, Z = 0.1) ## time steps times - c(seq(0, 10, 0.1)) ## numerical solution of ODE out - ode(y = y, parms = pars, times = times, func = derivs) yobs - data.frame( time = 0:7, Y = c(0.00, 3.40, 4.60 ,5.80, 5.80, 6.00, 6.00, 6.00), Z = c(0.1, 0.11, 0.119, 0.128, 0.136, 0.145, 0.153, 0.16) ) modelCost - function(p) { out - ode(y = y, func = derivs, parms = p, times = yobs$time) return(modCost(out, yobs, weight=mean)) } ## start values for the parameters startpars - c(a1 = 5, a2 = 0.002, a3 = 0.002, a4 = 0.02, Ymax=7) plower - c(a1 = 1, a2 = 0.0001, a3 = 0.0001, a4 = 0.0001, Ymax=0.001) pupper - c(a1 = 10, a2 = 2, a3 = 1, a4 = 1, Ymax=10) ## fit the model; nprint = 1 shows intermediate results fit - modFit(f = modelCost, p = startpars, upper=pupper, lower=plower, control = list(nprint = 1)) summary(fit) ## graphical result out2 - ode(y = y, parms = startpars, times = times, func = derivs) out3 - ode(y = y, parms = fit$par, times = times, func = derivs) plot(out, out2, out3, obs = yobs) legend(topleft, legend=c(original, startpars, fitted), col = 1:3, lty = 1:3) summary(fit) #Parameters: # Estimate Std. Error t value Pr(|t|) #a1 5.610e+00 2.913e+03 0.0020.998 #a2 2.000e+00 2.908e+03 0.0010.999 #a3 1.872e-03 3.162e-03 0.5920.566 #a4 1.188e-04 1.224e-01 0.0010.999 #Ymax 8.908e+00 4.615e+03 0.0020.998 # #Residual standard error: 0.08001 on 11 degrees of freedom # #Parameter correlation: # a1 a2 a3 a4 Ymax #a11.0 1.0 -0.09916 -0.1031 1.0 #a21.0 1.0 -0.09917 -0.1032 1.0 #a3 -0.09916 -0.09917 1.0 0.9981 -0.09917 #a4 -0.10315 -0.10315 0.99807 1. -0.10316 #Ymax 1.0 1.0 -0.09917 -0.1032 1.0 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Difference between spec.pgram spec.ar
On Apr 9, 2012, at 16:55 , Uwe Ligges wrote: On 08.04.2012 20:39, Bazman76 wrote: Hi there, Can someone explain what the difference between spec.pgram and spec.ar is? I understand that they attempt to do the same thing one using an AR estimation of the underlying series to estimate teh sensity the other using the FFT. However when applied to teh same data set they seem to be giving quite different results? http://r.789695.n4.nabble.com/file/n4541358/R_example.jpg Clearly the spec.ar() result seems to be smoothed but the results are also generally very different only really sharing the peak as the frequencies go towards zero. Can someone please explain why these two functions produce such different results on the same data set? Because they really measure different things? Why do you expect to get the same output in time as well as in frequency domain? That wasn't the question There is a raw series and two spectra, and it is the difference between the latter that is remarkable. Offhand, this is of course a bit like comparing a model fit to a nonparametric smoother of ordinary measurements: Sometimes the data just do something that the model says that they can't do and you get huge discrepancies. So the first reaction must be that an autoregressive model is just wrong for these data. My second reaction would be that the original series look a bit like they might have an asymmetric distribution, so perhaps a log or a square root transformation is warranted. That being said, I'm still quite stumped trying to explain the pattern in the raw periodogram. Notice that what you are seeing is not so much an initial peak as a region in which the spectrum drops effectively to zero. Have the data by any chance been subject to some sort of preprocessing that removes low-frequency components? At any rate, this isn't really about R, is it? -- Peter Dalgaard, Professor, Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Email: pd@cbs.dk Priv: pda...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Difference between spec.pgram spec.ar
oops sorry n 08.04.2012 20:39, Bazman76 wrote: Hi there, Can someone explain what the difference between spec.pgram and spec.ar is? I understand that they attempt to do the same thing one using an AR estimation of the underlying series to estimate teh sensity the other using the FFT. However when applied to teh same data set they seem to be giving quite different results? http://r.789695.n4.nabble.com/file/n4541358/R_example.jpg Clearly the spec.ar() result seems to be smoothed but the results are also generally very different only really sharing the peak as the frequencies go towards zero. Can someone please explain why these two functions produce such different results on the same data set? « [hide part of quote] *Because they really measure different things? Why do you expect to get the same output in time as well as in frequency domain? * Uwe Ligges code shown below: library(waveslim) vols=read.csv(file=C:/Users/ocuk/My Documents/Abs Vol.csv, header=TRUE, sep=,) x-ts(vols[,1]) #x ## LA(8) vol.la8- mra(x, la8, 4, modwt) names(vol.la8)- c(d1, d2, d3, d4, s4) ## plot multiresolution analysis of IBM data #par(mfcol=c(6,1), pty=m, mar=c(5-2,4,4-2,2)) par(mfcol=c(3,1), pty=m, mar=c(5-2,4,4-2,2)) plot.ts(x, axes=F, ylab=, main=(abs rtns)) #for(i in 1:5) # plot.ts(vol.la8[[i]], axes=F, ylab=names(vol.la8)[i]) #axis(side=1, at=seq(0,1600,by=100), # labels=c(0,,200,,400,,600,,800,,1000,,1200,,1400,,1600)) spectrum(vol.la8[[1]]) specx- spec.ar(x) -- [hidden email] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. « [hide part of quote] *OK so I neeed to understan better what it it they are trying to measure. I understood (incorrectly it seems) that they were simply different methods to get the same result? * -- View this message in context: http://r.789695.n4.nabble.com/Difference-between-spec-pgram-spec-ar-tp4541358p4543140.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Help with Book example of Matrix application
I found this example in an Introductory R book in the chapter on Matrices and Arrays The array is m [,1] [,2] [,3] [,4] [,5] [1,]0 12 138 20 [2,] 120 15 28 88 [3,] 13 15069 [4,]8 2860 33 [5,] 20 889 330 The code is #returns the minimum value of d[i,j], i !=j and the row attaining #the minimum, for square symmetric d; no special policy on ties mind - function(d) { n - nrow(d) # add a column to identify row number for apply() dd - cbind(d,1:n) wmins - apply(dd[-n,],1,imin) # wins will be 2 X n, 1st row being indices and 2nd being values i - which.min(wmins[2,]) j - wmins[1,i] return(c(d[i,j],i,j)) } #this finds the location, value of the minimum in a row x imin - function(x) { lx - length(x) i - x[lx] # original row number j -which.min(x[(x+1):(lx-1)]) k - i+j return(c(k,x[k])) } The result s/b [1] 6 3 4 I am getting: Warning messages: 1: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used 2: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used 3: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used 4: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used I have check my typing a number of times. Does anyone see an error? Thanks, Jim L. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Help with Book example of Matrix application
On 09/04/2012 12:32 PM, James Lenihan wrote: I found this example in an Introductory R book in the chapter on Matrices and Arrays You should probably check with the author of the book (there might be an errata page posted somewhere), but it looks like a typo, in your code or the original: The array is m [,1] [,2] [,3] [,4] [,5] [1,]0 12 138 20 [2,] 120 15 28 88 [3,] 13 15069 [4,]8 2860 33 [5,] 20 889 330 The code is #returns the minimum value of d[i,j], i !=j and the row attaining #the minimum, for square symmetric d; no special policy on ties mind- function(d) { n- nrow(d) # add a column to identify row number for apply() dd- cbind(d,1:n) wmins- apply(dd[-n,],1,imin) # wins will be 2 X n, 1st row being indices and 2nd being values i- which.min(wmins[2,]) j- wmins[1,i] return(c(d[i,j],i,j)) } #this finds the location, value of the minimum in a row x imin- function(x) { lx- length(x) i- x[lx] # original row number j-which.min(x[(x+1):(lx-1)]) That line doesn't make sense. Putting i in place of the innermost x would make more sense, but wouldn't work in all cases: if i == lx, you'll get a garbage answer. Duncan Murdoch k- i+j return(c(k,x[k])) } The result s/b [1] 6 3 4 I am getting: Warning messages: 1: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used 2: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used 3: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used 4: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used I have check my typing a number of times. Does anyone see an error? Thanks, Jim L. [[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] Help with Book example of Matrix application
On Apr 9, 2012, at 12:32 PM, James Lenihan wrote: I found this example in an Introductory R book in the chapter on Matrices and Arrays The array is m [,1] [,2] [,3] [,4] [,5] [1,]0 12 138 20 [2,] 120 15 28 88 [3,] 13 15069 [4,]8 2860 33 [5,] 20 889 330 The code is #returns the minimum value of d[i,j], i !=j and the row attaining #the minimum, for square symmetric d; no special policy on ties Would you consider a more compact operation? That code seems very unR- ish. min( m[row(m) != col(m)] ) #[1] 6 which(m == min( m[row(m) != col(m)] ), arr.ind=TRUE) # row col [1,] 4 3 [2,] 3 4 --- c(min= min( m[row(m) != col(m)] ), which(m == min( m[row(m) != col(m)] ), arr.ind=TRUE)[1,]) min row col 6 4 3 I suppose the question could be how to debug that code, in which case you should be manually stepping through it to see what the values are during the implicit loop. mind - function(d) { n - nrow(d) # add a column to identify row number for apply() dd - cbind(d,1:n) wmins - apply(dd[-n,],1,imin) # wins will be 2 X n, 1st row being indices and 2nd being values i - which.min(wmins[2,]) j - wmins[1,i] return(c(d[i,j],i,j)) } #this finds the location, value of the minimum in a row x imin - function(x) { lx - length(x) i - x[lx] # original row number j -which.min(x[(x+1):(lx-1)]) k - i+j return(c(k,x[k])) } The result s/b [1] 6 3 4 I am getting: Warning messages: 1: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used 2: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used 3: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used 4: In (x + 1):(lx - 1) : numerical expression has 6 elements: only the first used I have check my typing a number of times. Does anyone see an error? Thanks, Jim L. [[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, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Appropriate method for sharing data across functions
Make OPCON an environment and pass it into the functions that may read it or alter it. There is no real need to pass it out, since environments are changed in-place (unlike lists). E.g., x - list2env(list(one=1, two=ii, three=3)) x environment: 0x03110890 objects(x) [1] one three two x[[two]] [1] ii with(x, three+one) [1] 4 f - function(z, env) { env[[newZ]] - z ; sqrt(z) } f(10, x) [1] 3.162278 x[[newZ]] # put there by f() [1] 10 The advantage of using a reference class is you get an object with defined behaviour, and something that you can document more easily. Hadley -- Assistant Professor / Dobelman Family Junior Chair Department of Statistics / Rice University http://had.co.nz/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Appropriate method for sharing data across functions
On 05/04/2012 4:20 PM, John C Nash wrote: In trying to streamline various optimization functions, I would like to have a scratch pad of working data that is shared across a number of functions. These can be called from different levels within some wrapper functions for maximum likelihood and other such computations. I'm sure there are other applications that could benefit from this. Below are two approaches. One uses the- assignment to a structure I call OPCON. The other attempts to create an environment with this name, but fails. Though I have looked at a number of references, I have so far not found an adequate description of how to specify where the OPCON environment is located. (Both the green and blue books do not cover this topic, at least not under environment in the index.) Is there a recommended approach to this? The one I would use is to create all of the functions that need access to the scratch pad as local functions within another. For example, makethem - function() { MAXIMIZE - TRUE PARSCALE - rep(1, npar) KFN - 0 ... etc ... add1 - function(){ KFN - 1 + KFN } list(add1 = add1, ... other functions here...) } A single call to makethem() will return a list of functions which all have access to MAXIMIZE, PARSCALE, etc. If you like, you can pull them out of the list to be called independently, e.g. fns - makethem() add1 - fns$add1 and they'll still have access. (If this is a one-off need, you can use local( ) around the definitions, but I prefer the wrapper-function style. I realize I could use argument lists, but they get long and tedious with the number of items I may need to pass, though passing the OPCON structure in and out might be the proper way. An onAttach() approach was suggested by Paul Gilbert and tried, but it has so far not succeeded and, unfortunately, does not seem to be usable from source() i.e., cannot be interpreted but must be built first. JN Example using- rm(list=ls()) optstart-function(npar){ # create structure for optimization computations # npar is number of parameters ?? test?? OPCON-list(MAXIMIZE=TRUE, PARSCALE=rep(1,npar), FNSCALE=1, KFN=0, KGR=0, KHESS=0) # may be other stuff ls(OPCON) } The code above creates OPCON as a global variable; that's dangerous, if you have multiple different optimizers potentially being used. add1-function(){ OPCON$KFN-1+OPCON$KFN test-OPCON$KFN } OPCON-list(MAXIMIZE=TRUE, PARSCALE=rep(1,4), FNSCALE=1, KFN=0, KGR=0, KHESS=0) ls(OPCON) print(add1()) print(add1()) print(ls.str()) rm(OPCON) # Try to remove the scratchpad print(ls()) tmp-readline(Now try from within a function) setup-optstart(4) # Need to sort out how to set this up appropriately cat(setup =) print(setup) print(add1()) print(add1()) rm(OPCON) # Try to remove the scratchpad == Example (failing) using new.env: rm(list=ls()) optstart-function(npar){ # create structure for optimization computations # npar is number of parameters ?? test?? OPCON-new.env(parent=globalenv()) OPCON-list(MAXIMIZE=TRUE, PARSCALE=rep(1,npar), FNSCALE=1, KFN=0, KGR=0, KHESS=0) This failed because the second assignment wiped out the first: you created a new environment, then threw it away when you assigned a list to OPCON. The other problem with this approach is that you never saved OPCON anywhere, so when optstart() was done, it would disappear. Hadley suggested using a reference class, and Bill suggested using an environment explicitly; those are both possible alternatives. My version creates an environment but never names it; it's the environment of the add1() function and other functions defined within makethem(). I'd say the more formal versions they suggested would be better if you had a longer (or open-ended) list of possible functions you wanted to work with, or if you wanted to pass your OPCON to the user to manipulate. Duncan Murdoch # may be other stuff ls(OPCON) } add1-function(){ OPCON$KFN-1+OPCON$KFN test-OPCON$KFN } OPCON-new.env(parent=globalenv()) OPCON-list(MAXIMIZE=TRUE, PARSCALE=rep(1,4), FNSCALE=1, KFN=0, KGR=0, KHESS=0) ls(OPCON) print(add1()) print(add1()) print(ls.str()) rm(OPCON) # Try to remove the scratchpad print(ls()) tmp-readline(Now try from within a function) setup-optstart(4) # Need to sort out how to set this up appropriately cat(setup =) print(setup) print(add1()) print(add1()) rm(OPCON) # Try to remove the scratchpad __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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
[R] Overall model significance for poisson GLM
Greetings, I am running glm models for species counts using a poisson link function. Normal summary functions for this provide summary statistics in the form of the deviance, AIC, and p-values for individual predictors. I would like to obtain the p-value for the overall model. So far, I have been using an analysis of deviance table to check a model against the null model with the intercept as the only predictor. Any advice on other methods to obtain the proper p-value would be appreciated. Thanks, Pat [[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] Difference between spec.pgram spec.ar
Yes I agree, there may be something pathalogical in the way at least one of the models handles the data. That's why I was trying to get a better handle on how the two functions spec.prgm() and spec.ar() work. The data has been processed by a wavelet analysis, so what you are seeing as the raw data is in fact the level1 details from the wavelet anlaysis. In theory this should only have high frequency components, that was why I am so surpirsed to see such strong components at low frequencies. That is not a R quesiton per se, but surely how spec.prgm() and spec.ar() work is? -- View this message in context: http://r.789695.n4.nabble.com/Difference-between-spec-pgram-spec-ar-tp4541358p4543191.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] binned tabulation
Hi, I am attempting to tabulate binned data. The '1' represents the appearance of the focal mouse pup, and '2' represents the disappearance of the focal mouse pup. The code written below is intended to calculate the total time spent appeared out of 3600s. For Sample 1, both the hand calculation and R code yield the same result, 50. A problem seems to occur when '1' is the last entry. For Sample 2, the total time appeared is 53 (hand calculation), however, using the R code below yields 55. If you have any suggestions for solving the problem, please let me know. Thank you in advance for any assistance you may provide. Delia Sample 1 0.0 1 40 2 45 1 55 2 Sample 2 0.0 1 40 2 45 1 55 2 57 1 name = read.table(file.choose(),header=F) # opening a data file colnames(name)-c(Time, Behavior) name = data.frame(name$Behavior, name$Time) colnames(name)-c(Behavior, Time) name-name[name$Time 3600, ]; ## run file x-seq(0,3600, by = 60) # total time partition by time which is 60 if (tail(name$Behavior, 1) == 1) {name-rbind(name, c(4, 3600))} else {name-rbind(name, c(1, 3600))} if (nrow(name) %% 2 != 0) {name -name[-c(nrow(name)), -c(nrow(name))]} q-NULL for (y in (1: nrow(name))) { if (y %% 2 != 0) q-c(q, (c(name$Time[y]:name$Time[y +1]))) } b-table(cut(q,x)) sum(b) [[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] Overall model significance for poisson GLM
Not really an R question. Post on a statistics site like http://stats.stackexchange.com/ -- Bert On Mon, Apr 9, 2012 at 9:51 AM, Pat Wilkins pwilk...@illinois.edu wrote: Greetings, I am running glm models for species counts using a poisson link function. Normal summary functions for this provide summary statistics in the form of the deviance, AIC, and p-values for individual predictors. I would like to obtain the p-value for the overall model. So far, I have been using an analysis of deviance table to check a model against the null model with the intercept as the only predictor. Any advice on other methods to obtain the proper p-value would be appreciated. Thanks, Pat [[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. -- 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 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Difference between spec.pgram spec.ar
On Mon, Apr 9, 2012 at 9:27 AM, Bazman76 h_a_patie...@hotmail.com wrote: Yes I agree, there may be something pathalogical in the way at least one of the models handles the data. That's why I was trying to get a better handle on how the two functions spec.prgm() and spec.ar() work. The data has been processed by a wavelet analysis, so what you are seeing as the raw data is in fact the level1 details from the wavelet anlaysis. In theory this should only have high frequency components, that was why I am so surpirsed to see such strong components at low frequencies. That is not a R quesiton per se, but surely how spec.prgm() and spec.ar() work is? Not necessarily, if they are e.g. C or Fortran programs merely called by R. Indeed, even if written in R, if the algorithms are the issue, then that is essentially a statistics/numerical analysis matter, not an R programming one. -- Bert -- View this message in context: http://r.789695.n4.nabble.com/Difference-between-spec-pgram-spec-ar-tp4541358p4543191.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. -- 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 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Listing the contents of an FTP directory via R?
R-helpers: I'd like to be able to store all the file information from an ftp site (e.g. file and foldernames) through an R command. Any ideas how to do this? Here's an example site to use: ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005 --j -- Jonathan A. Greenberg, PhD Assistant Professor Department of Geography and Geographic Information Science University of Illinois at Urbana-Champaign 607 South Mathews Avenue, MC 150 Urbana, IL 61801 Phone: 415-763-5476 AIM: jgrn307, MSN: jgrn...@hotmail.com, Gchat: jgrn307, Skype: jgrn3007 http://www.geog.illinois.edu/people/JonathanGreenberg.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] Listing the contents of an FTP directory via R?
A couple of ways. using Rcurl you can use the curlOption of dirlistonly. otherwise you can read the page and parse. I've got some code around here to do that. Steve On Mon, Apr 9, 2012 at 11:27 AM, Jonathan Greenberg j...@illinois.eduwrote: R-helpers: I'd like to be able to store all the file information from an ftp site (e.g. file and foldernames) through an R command. Any ideas how to do this? Here's an example site to use: ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005 --j -- Jonathan A. Greenberg, PhD Assistant Professor Department of Geography and Geographic Information Science University of Illinois at Urbana-Champaign 607 South Mathews Avenue, MC 150 Urbana, IL 61801 Phone: 415-763-5476 AIM: jgrn307, MSN: jgrn...@hotmail.com, Gchat: jgrn307, Skype: jgrn3007 http://www.geog.illinois.edu/people/JonathanGreenberg.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. [[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] Difference between spec.pgram spec.ar
On 09/04/2012 18:52, Bert Gunter wrote: On Mon, Apr 9, 2012 at 9:27 AM, Bazman76h_a_patie...@hotmail.com wrote: Yes I agree, there may be something pathalogical in the way at least one of the models handles the data. That's why I was trying to get a better handle on how the two functions spec.prgm() and spec.ar() work. The data has been processed by a wavelet analysis, so what you are seeing as the raw data is in fact the level1 details from the wavelet anlaysis. In theory this should only have high frequency components, that was why I am so surpirsed to see such strong components at low frequencies. That is not a R quesiton per se, but surely how spec.prgm() and spec.ar() work is? Not necessarily, if they are e.g. C or Fortran programs merely called by R. Indeed, even if written in R, if the algorithms are the issue, then that is essentially a statistics/numerical analysis matter, not an R programming one. Which is why the help pages often (and do here) have definitive references. So the best course of action is to start following up the references. -- 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] Appropriate method for sharing data across functions
Thanks to Hadley, William and Duncan for suggestions. I'm currently implementing a solution that is close to that of William and Duncan (and learning more about environments in the process). I suspect the reference classes are possibly a more reliable long term solution. I'll plead laziness until I find the other approach won't work. For anyone interested in these issues, over the next few weeks the outcome of my investigations and trials should be on the SVN repository of R-forge under Optimization and Solving Tools, https://r-forge.r-project.org/R/?group_id=395, in particular within the packages optfntools and optimx. However, it will be late April before things stabilize. John Nash __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Listing the contents of an FTP directory via R?
Steven: Thanks -- I seem to be running into the problem with the link I sent along: getURL(ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005/,verbose=TRUE,dirlistonly = TRUE) Error in function (type, msg, asError = TRUE) : RETR response: 550 I'm wondering if it might be a passive ftp issue, but: getURL(ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005/,verbose=TRUE,ftp.use.epsv=TRUE, dirlistonly = TRUE) Error in function (type, msg, asError = TRUE) : FTP response reading failed Does not seem to work... Thoughts? --j On Mon, Apr 9, 2012 at 1:32 PM, steven mosher mosherste...@gmail.com wrote: A couple of ways. using Rcurl you can use the curlOption of dirlistonly. otherwise you can read the page and parse. I've got some code around here to do that. Steve On Mon, Apr 9, 2012 at 11:27 AM, Jonathan Greenberg j...@illinois.edu wrote: R-helpers: I'd like to be able to store all the file information from an ftp site (e.g. file and foldernames) through an R command. Any ideas how to do this? Here's an example site to use: ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005 --j -- Jonathan A. Greenberg, PhD Assistant Professor Department of Geography and Geographic Information Science University of Illinois at Urbana-Champaign 607 South Mathews Avenue, MC 150 Urbana, IL 61801 Phone: 415-763-5476 AIM: jgrn307, MSN: jgrn...@hotmail.com, Gchat: jgrn307, Skype: jgrn3007 http://www.geog.illinois.edu/people/JonathanGreenberg.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. -- Jonathan A. Greenberg, PhD Assistant Professor Department of Geography and Geographic Information Science University of Illinois at Urbana-Champaign 607 South Mathews Avenue, MC 150 Urbana, IL 61801 Phone: 415-763-5476 AIM: jgrn307, MSN: jgrn...@hotmail.com, Gchat: jgrn307, Skype: jgrn3007 http://www.geog.illinois.edu/people/JonathanGreenberg.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] Listing the contents of an FTP directory via R?
Ya I hit the same error with my code that reads directories. I'll try some other stuff . I think I hit this error before I with usgs. On Mon, Apr 9, 2012 at 11:40 AM, Jonathan Greenberg j...@illinois.eduwrote: Steven: Thanks -- I seem to be running into the problem with the link I sent along: getURL(ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005/,verbose=TRUE,dirlistonly = TRUE) Error in function (type, msg, asError = TRUE) : RETR response: 550 I'm wondering if it might be a passive ftp issue, but: getURL(ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005/,verbose=TRUE,ftp.use.epsv=TRUE, dirlistonly = TRUE) Error in function (type, msg, asError = TRUE) : FTP response reading failed Does not seem to work... Thoughts? --j On Mon, Apr 9, 2012 at 1:32 PM, steven mosher mosherste...@gmail.com wrote: A couple of ways. using Rcurl you can use the curlOption of dirlistonly. otherwise you can read the page and parse. I've got some code around here to do that. Steve On Mon, Apr 9, 2012 at 11:27 AM, Jonathan Greenberg j...@illinois.edu wrote: R-helpers: I'd like to be able to store all the file information from an ftp site (e.g. file and foldernames) through an R command. Any ideas how to do this? Here's an example site to use: ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005 --j -- Jonathan A. Greenberg, PhD Assistant Professor Department of Geography and Geographic Information Science University of Illinois at Urbana-Champaign 607 South Mathews Avenue, MC 150 Urbana, IL 61801 Phone: 415-763-5476 AIM: jgrn307, MSN: jgrn...@hotmail.com, Gchat: jgrn307, Skype: jgrn3007 http://www.geog.illinois.edu/people/JonathanGreenberg.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. -- Jonathan A. Greenberg, PhD Assistant Professor Department of Geography and Geographic Information Science University of Illinois at Urbana-Champaign 607 South Mathews Avenue, MC 150 Urbana, IL 61801 Phone: 415-763-5476 AIM: jgrn307, MSN: jgrn...@hotmail.com, Gchat: jgrn307, Skype: jgrn3007 http://www.geog.illinois.edu/people/JonathanGreenberg.html [[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] Problems with CRAN package Ryacas
Apropos: I don't have the problems, the OP had, but on my ubuntu notebook, Ryacas does not return expressions (just the strings), and hence as.expression( yacas-result ) always gives NULL and e.g. the demo(Ryacas-Function) also fails: yacas(expression(deriv(BurrCDF(x,c,k k*c*x^(c-1)*(x^c+1)^(-(k+1)); yy - yacas(expression(deriv(BurrCDF(x,c,k yy k*c*x^(c-1)*(x^c+1)^(-(k+1)); str(yy) List of 2 $ : NULL $ YacasForm: chr k*c*x^(c-1)*(x^c+1)^(-(k+1)); - attr(*, class)= chr yacas I have - yacas 1.2.2 (standard ubuntu package) - re-installed Ryacas (manually from source, just now, under R 2.15.0 patched) sessionInfo() R version 2.15.0 Patched (2012-04-09 r58947) Platform: x86_64-unknown-linux-gnu (64-bit) locale: [1] LC_CTYPE=de_CH.UTF-8 LC_NUMERIC=C LC_TIME=en_US.UTF-8 [4] LC_COLLATE=de_CH.UTF-8 LC_MONETARY=en_US.UTF-8 LC_MESSAGES=de_CH.UTF-8 [7] LC_PAPER=C LC_NAME=C LC_ADDRESS=C [10] LC_TELEPHONE=C LC_MEASUREMENT=de_CH.UTF-8 LC_IDENTIFICATION=C attached base packages: [1] graphics grDevices datasets stats utils methods base other attached packages: [1] Ryacas_0.2-11 XML_3.9-4 fortunes_1.5-0 sfsmisc_1.0-20 loaded via a namespace (and not attached): [1] compiler_2.15.0 tools_2.15.0 -- On our Fedora computers at work (with an *older* version of yacas, installed manually from sources there), things work fine : str(yacas(expression(Factor(x^2-1 List of 2 $ text : expression((x + 1) * (x - 1)) $ OMForm: chr [1:15] OMOBJ OMA OMS cd=\arith1\ name=\times\/ OMA ... - attr(*, class)= chr yacas as.expression(yacas(expression(Factor(x^2-1 expression((x + 1) * (x - 1)) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Listing the contents of an FTP directory via R?
Its kinda gross, but you can just download.file() and then parse the result !DOCTYPE HTML PUBLIC -//IETF//DTD HTML//EN HTML HEAD TITLEFTP directory /MOTA/MCD15A3.005 at e4ftl01.cr.usgs.gov/TITLE /HEAD BODY H2 ID=WinINetFtpDirectoryFTP directory /MOTA/MCD15A3.005 at e4ftl01.cr.usgs.gov/H2 HR H4PRE ## ## *** WARNING TO USERS OF THIS SYSTEM *** THIS COMPUTER SYSTEM, INCLUDING ALL RELATED EQUIPMENT, NETWORKS, AND NETWORK DEVICES (INCLUDING INTERNET ACCESS), IS PROVIDED BY THE DEPARTMENT OF THE INTERIOR (DOI) IN ACCORDANCE WITH THE AGENCY POLICY FOR OFFICIAL USE AND LIMITED PERSONAL USE. ALL AGENCY COMPUTER SYSTEMS MAY BE MONITORED FOR ALL LAWFUL PURPOSES, INCLUDING BUT NOT LIMITED TO, ENSURING THAT USE IS AUTHORIZED, FOR MANAGEMENT OF THE SYSTEM, TO FACILITATE PROTECTION AGAINST UNAUTHORIZED ACCESS, AND TO VERIFY SECURITY PROCEDURES, SURVIVABILITY AND OPERATIONAL SECURITY. ANY INFORMATION ON THIS COMPUTER SYSTEM MAY BE EXAMINED, RECORDED, COPIED AND USED FOR AUTHORIZED PURPOSES AT ANY TIME. ALL INFORMATION, INCLUDING PERSONAL INFORMATION, PLACED OR SENT OVER THIS SYSTEM MAY BE MONITORED, AND USERS OF THIS SYSTEM ARE REMINDED THAT SUCH MONITORING DOES OCCUR. THEREFORE, THERE SHOULD BE NO EXPECTATION OF PRIVACY WITH RESPECT TO USE OF THIS SYSTEM. BY LOGGING INTO THIS AGENCY COMPUTER SYSTEM, YOU ACKNOWLEDGE AND CONSENT TO THE MONITORING OF THIS SYSTEM. EVIDENCE OF YOUR USE, AUTHORIZED OR UNAUTHORIZED, COLLECTED DURING MONITORING MAY BE USED FOR CIVIL, CRIMINAL, ADMINISTRATIVE, OR OTHER ADVERSE ACTION. UNAUTHORIZED OR ILLEGAL USE MAY SUBJECT YOU TO PROSECUTION. ## ## /PRE/H4 HR A HREF=..Up to higher level directory/ABRPRE 10/18/2011 05:23PM Directory A HREF=/MOTA/MCD15A3.005/2002.07.04/B2002.07.04/B/A 10/18/2011 05:08PM Directory A HREF=/MOTA/MCD15A3.005/2002.07.08/B2002.07.08/B/A 10/18/2011 05:05PM Directory A HREF=/MOTA/MCD15A3.005/2002.07.12/B2002.07.12/B/A 10/18/2011 06:00PM Directory A HREF=/MOTA/MCD15A3.005/2002.07.16/B2002.07.16/B/A 10/18/2011 06:00PM Directory A HREF=/MOTA/MCD15A3.005/2002.07.20/B2002.07.20/B/A 10/18/2011 06:00PM Directory A HREF=/MOTA/MCD15A3.005/2002.07.24/B2002.07.24/B/A 10/18/2011 06:01PM Directory A HREF=/MOTA/MCD15A3.005/2002.07.28/B2002.07.28/B/A 10/18/2011 06:01PM Directory A HREF=/MOTA/MCD15A3.005/2002.08.05/B2002.08.05/B/A 10/18/2011 06:01PM Directory A HREF=/MOTA/MCD15A3.005/2002.08.09/B2002.08.09/B/A 10/18/2011 05:55PM Directory A HREF=/MOTA/MCD15A3.005/2002.08.13/B2002.08.13/B/A 10/18/2011 05:54PM Directory A HREF=/MOTA/MCD15A3.005/2002.08.17/B2002.08.17/B/A 10/18/2011 05:37PM Directory A HREF=/MOTA/MCD15A3.005/2002.08.21/B2002.08.21/B/A 10/18/2011 04:26PM Directory A HREF=/MOTA/MCD15A3.005/2002.08.25/B2002.08.25/B/A 10/18/2011 05:19PM Directory A HREF=/MOTA/MCD15A3.005/2002.08.29/B2002.08.29/B/A 10/18/2011 05:19PM Directory A HREF=/MOTA/MCD15A3.005/2002.09.02/B2002.09.02/B/A 10/18/2011 05:14PM Directory A HREF=/MOTA/MCD15A3.005/2002.09.06/B2002.09.06/B/A 10/18/2011 03:39PM Directory A HREF=/MOTA/MCD15A3.005/2002.09.10/B2002.09.10/B/A 10/18/2011 05:23PM Directory A HREF=/MOTA/MCD15A3.005/2002.09.14/B2002.09.14/B/A 10/18/2011 05:25PM Directory A HREF=/MOTA/MCD15A3.005/2002.09.18/B2002.09.18/B/A 10/18/2011 05:25PM Directory A HREF=/MOTA/MCD15A3.005/2002.09.22/B2002.09.22/B/A 10/18/2011 05:21PM Directory A HREF=/MOTA/MCD15A3.005/2002.09.26/B2002.09.26/B/A 10/18/2011 05:59PM Directory A HREF=/MOTA/MCD15A3.005/2002.09.30/B2002.09.30/B/A 10/18/2011 06:00PM Directory A HREF=/MOTA/MCD15A3.005/2002.10.04/B2002.10.04/B/A 10/18/2011 04:37PM Directory A HREF=/MOTA/MCD15A3.005/2002.10.08/B2002.10.08/B/A 10/18/2011 04:03PM Directory A HREF=/MOTA/MCD15A3.005/2002.10.12/B2002.10.12/B/A 10/18/2011 06:00PM Directory A HREF=/MOTA/MCD15A3.005/2002.10.16/B2002.10.16/B/A 10/18/2011 05:36PM Directory A HREF=/MOTA/MCD15A3.005/2002.10.20/B2002.10.20/B/A 10/18/2011 03:54PM Directory A HREF=/MOTA/MCD15A3.005/2002.10.24/B2002.10.24/B/A 10/18/2011 03:54PM Directory A HREF=/MOTA/MCD15A3.005/2002.10.28/B2002.10.28/B/A 10/18/2011 03:50PM Directory A HREF=/MOTA/MCD15A3.005/2002.11.01/B2002.11.01/B/A 10/18/2011 01:25PM Directory A HREF=/MOTA/MCD15A3.005/2002.11.05/B2002.11.05/B/A 10/18/2011 01:27PM Directory A HREF=/MOTA/MCD15A3.005/2002.11.09/B2002.11.09/B/A 10/18/2011 01:29PM Directory A HREF=/MOTA/MCD15A3.005/2002.11.13/B2002.11.13/B/A 10/18/2011 01:31PM Directory A
Re: [R] Listing the contents of an FTP directory via R?
It seems to be choking on NLST: require(RCurl) getURL(ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005/,verbose=TRUE,ftp.use.epsv=TRUE, dirlistonly = TRUE) ... 230 Guest login ok, access restrictions apply. PWD 257 / is current directory. * Entry path is '/' CWD MOTA 250 CWD command successful. CWD MCD15A3.005 250 CWD command successful. EPSV * Connect data stream passively 500 'EPSV': command not understood. * disabling EPSV usage PASV 227 Entering Passive Mode (152,61,133,5,177,17) * Trying 152.61.133.5... * connected * Connecting to 152.61.133.5 (152.61.133.5) port 45329 TYPE A 200 Type set to A. NLST 550 No files found. * RETR response: 550 * Remembering we are in dir MOTA/MCD15A3.005/ * Connection #0 to host e4ftl01.cr.usgs.gov left intact Error in function (type, msg, asError = TRUE) : RETR response: 550 Thanks for looking into this! --j On Mon, Apr 9, 2012 at 1:45 PM, steven mosher mosherste...@gmail.com wrote: Ya I hit the same error with my code that reads directories. I'll try some other stuff . I think I hit this error before I with usgs. On Mon, Apr 9, 2012 at 11:40 AM, Jonathan Greenberg j...@illinois.edu wrote: Steven: Thanks -- I seem to be running into the problem with the link I sent along: getURL(ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005/,verbose=TRUE,dirlistonly = TRUE) Error in function (type, msg, asError = TRUE) : RETR response: 550 I'm wondering if it might be a passive ftp issue, but: getURL(ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005/,verbose=TRUE,ftp.use.epsv=TRUE, dirlistonly = TRUE) Error in function (type, msg, asError = TRUE) : FTP response reading failed Does not seem to work... Thoughts? --j On Mon, Apr 9, 2012 at 1:32 PM, steven mosher mosherste...@gmail.com wrote: A couple of ways. using Rcurl you can use the curlOption of dirlistonly. otherwise you can read the page and parse. I've got some code around here to do that. Steve On Mon, Apr 9, 2012 at 11:27 AM, Jonathan Greenberg j...@illinois.edu wrote: R-helpers: I'd like to be able to store all the file information from an ftp site (e.g. file and foldernames) through an R command. Any ideas how to do this? Here's an example site to use: ftp://e4ftl01.cr.usgs.gov/MOTA/MCD15A3.005 --j -- Jonathan A. Greenberg, PhD Assistant Professor Department of Geography and Geographic Information Science University of Illinois at Urbana-Champaign 607 South Mathews Avenue, MC 150 Urbana, IL 61801 Phone: 415-763-5476 AIM: jgrn307, MSN: jgrn...@hotmail.com, Gchat: jgrn307, Skype: jgrn3007 http://www.geog.illinois.edu/people/JonathanGreenberg.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. -- Jonathan A. Greenberg, PhD Assistant Professor Department of Geography and Geographic Information Science University of Illinois at Urbana-Champaign 607 South Mathews Avenue, MC 150 Urbana, IL 61801 Phone: 415-763-5476 AIM: jgrn307, MSN: jgrn...@hotmail.com, Gchat: jgrn307, Skype: jgrn3007 http://www.geog.illinois.edu/people/JonathanGreenberg.html -- Jonathan A. Greenberg, PhD Assistant Professor Department of Geography and Geographic Information Science University of Illinois at Urbana-Champaign 607 South Mathews Avenue, MC 150 Urbana, IL 61801 Phone: 415-763-5476 AIM: jgrn307, MSN: jgrn...@hotmail.com, Gchat: jgrn307, Skype: jgrn3007 http://www.geog.illinois.edu/people/JonathanGreenberg.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] Problems with CRAN package Ryacas
On Mon, Apr 9, 2012 at 2:47 PM, Martin Maechler maech...@stat.math.ethz.ch wrote: Apropos: I don't have the problems, the OP had, but on my ubuntu notebook, Ryacas does not return expressions (just the strings), and hence as.expression( yacas-result ) always gives NULL and e.g. the demo(Ryacas-Function) also fails: yacas(expression(deriv(BurrCDF(x,c,k k*c*x^(c-1)*(x^c+1)^(-(k+1)); yy - yacas(expression(deriv(BurrCDF(x,c,k yy k*c*x^(c-1)*(x^c+1)^(-(k+1)); str(yy) List of 2 $ : NULL $ YacasForm: chr k*c*x^(c-1)*(x^c+1)^(-(k+1)); - attr(*, class)= chr yacas I have - yacas 1.2.2 (standard ubuntu package) The latest version of Ryacas on CRAN is 0.2-11 http://cran.r-project.org/package=Ryacas -- Statistics Software Consulting GKX Group, GKX Associates Inc. tel: 1-877-GKX-GROUP email: ggrothendieck at gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problems with CRAN package Ryacas
On Mon, Apr 9, 2012 at 4:02 PM, Gabor Grothendieck ggrothendi...@gmail.com wrote: On Mon, Apr 9, 2012 at 2:47 PM, Martin Maechler maech...@stat.math.ethz.ch wrote: Apropos: I don't have the problems, the OP had, but on my ubuntu notebook, Ryacas does not return expressions (just the strings), and hence as.expression( yacas-result ) always gives NULL and e.g. the demo(Ryacas-Function) also fails: yacas(expression(deriv(BurrCDF(x,c,k k*c*x^(c-1)*(x^c+1)^(-(k+1)); yy - yacas(expression(deriv(BurrCDF(x,c,k yy k*c*x^(c-1)*(x^c+1)^(-(k+1)); str(yy) List of 2 $ : NULL $ YacasForm: chr k*c*x^(c-1)*(x^c+1)^(-(k+1)); - attr(*, class)= chr yacas I have - yacas 1.2.2 (standard ubuntu package) The latest version of Ryacas on CRAN is 0.2-11 http://cran.r-project.org/package=Ryacas and the version of yacas supported by Ryacas is: yacas(Version()) expression(1.0.63) The troubleshooting section on the home page that I referred to does mention this. -- Statistics Software Consulting GKX Group, GKX Associates Inc. tel: 1-877-GKX-GROUP email: ggrothendieck at gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] sdev, variance in prcomp
Hello, It might be a trivial question but I just wanted to find out the relationship between sdev and proportion of variance generated by prcomp. I got the following result from my data set PC1 PC2 PC3 Standard deviation 104.89454 15.40910 9.012047 Proportion of Variance 0.52344 0.01130 0.003860 Cumulative Proportion 0.52344 0.53474 0.538600 first, I had thought that the variance must be standard deviation power to 2 but it seems that variance is sdev/200. Did I misunderstood some thing? Best, Carol __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] linp code for fuzzy regression
Hi! I study fuzzy regression. Linear Programming(LP) methods are commonly used for fuzzy linear regression (FLR) because they are simple and easy to apply. I should make a simulation study with R and have come a distance. I use linp (for instence a part from my codes) (L-linp(E=NULL,F=NULL,Cost=obj,G=sınır,H=sol,ispos=FALSE)) there is ispose argument which allows all the variables(x) to get negative values when its false. my problem is there. I want to allow some of the variables(x)to get negative values and the others must be nonnegative.(for instence x1,x2,x3 could be negative or positive and x4,x5,x6 must be nonnegative.) Is there anything to solve that problem? thanks for any help. Bekir. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] sdev, variance in prcomp
It is reporting _proportion_ of the total variance. Not the variance. The total variance is the sum of all the squared standard deviations. Only the first 3 standard deviations are shown. Kevin On Mon, Apr 9, 2012 at 3:31 PM, carol white wht_...@yahoo.com wrote: Hello, It might be a trivial question but I just wanted to find out the relationship between sdev and proportion of variance generated by prcomp. I got the following result from my data set PC1 PC2 PC3 Standard deviation 104.89454 15.40910 9.012047 Proportion of Variance 0.52344 0.01130 0.003860 Cumulative Proportion0.52344 0.53474 0.538600 first, I had thought that the variance must be standard deviation power to 2 but it seems that variance is sdev/200. Did I misunderstood some thing? Best, Carol __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Kevin Wright [[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] Pairwise comparison matrix elements
Hi!, I'm really hoping someone out there will be able to help me. I recently started my MSc dissertation on Population Projection Matrices, which has been going well until now. I am trying to set-up a general script that does a pairwise comparison of all elements in my matrices. So for example, given that I have the following matrix S: S [,1] [,2] [,3] [1,] 0.0 0.007361183 0.0001634936 [2,] 13.88458 0.353988378 0.00 [3,] 0.0 16.086367098 0.3572819557 I'd like to create a matrix SA for each element in my matrix compared to another element like so: SA S-[1,1] [,1][,2][,3] [1,] [1,1]-[1,1] [1,2]-[1,1] [1,3]-[1,1] [2,] [2,1]-[1,1] [2,2]-[1,1] [2,3]-[1,1] [3,] [3,1]-[1,1] [3,2]-[1,1] [3,3]-[1,1] S-[1,2] [...] S-[1,3] [...] ... and so on The aim is to be able to prove existing rules and trends across matrix dimensions. For example, I've been trying to test whether the first line of values decreases or remains the same from column 1 to the penultimate column: S[1,j] = S[1,lj] The last thing I tried was: M-Matlab2R([0 24.724 1377.48;0.029328 0.26 0;0 0.014 0.78]) w - abs(Re(eigen(M)$vectors[,1])) v - abs(Re(eigen(t(M))$vectors[,1])) w - w/sum(w) v - v/(t(v)%*%w) S - (v%*%t(w))/as.vector(t(v)%*%w) S n- length(S) M.comp - array(S,dim=c(n,n,n,n)) for (i in 1:n) { for (j in 1:n) { for (k in 1:n) { for (l in 1:n) { poscompa[1,3,1,j-1] - S[i,j]=S[k,l] # if=1, then TRUE and if=0, then FALSE M.comp sum(poscompa[1,3,1,1:n]) # should give me the value for n. I know this doesn't test it accurately, but first I need to get the loops to work right. I'm getting an array of errors such as can't find function poscompa (or compa), subscripts out of range or type of closure not valid regardless of what I try. Using COMPAR (poscompa and compa) was the last recommendation I got, but I'm starting to wonder if there might be other ways to go about this. All out-of-the-box ideas I've come up with and tried haven't gotten me much farther. I've now practically exhausted my creative thinking, and I'm becoming very frustrated. I'd really like to get this script going since my current one would make my life hell for large populations (60x60 population matrices). Any ideas on how I could move forward? Many, many thanks! Marta [[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] binned tabulation
Delia, name - data.frame(Behavior=c(1, 2, 1, 2, 1), Time=c(0, 40, 45, 55, 57)) appear - name$Time[name$Behavior==1] disappear - name$Time[name$Behavior==2] if(length(appear) length(disappear)) disappear - c(disappear, 60) sum(disappear - appear) Delia Shelton wrote on 04/09/2012 12:30:23 PM: Hi, I am attempting to tabulate binned data. The '1' represents the appearance of the focal mouse pup, and '2' represents the disappearance of the focal mouse pup. The code written below is intended to calculate the total time spent appeared out of 3600s. For Sample 1, both the hand calculation and R code yield the same result, 50. A problem seems to occur when '1' is the last entry. For Sample 2, the total time appeared is 53 (hand calculation), however, using the R code below yields 55. If you have any suggestions for solving the problem, please let me know. Thank you in advance for any assistance you may provide. Delia Sample 1 0.0 1 40 2 45 1 55 2 Sample 2 0.0 1 40 2 45 1 55 2 57 1 name = read.table(file.choose(),header=F) # opening a data file colnames(name)-c(Time, Behavior) name = data.frame(name$Behavior, name$Time) colnames(name)-c(Behavior, Time) name-name[name$Time 3600, ]; ## run file x-seq(0,3600, by = 60) # total time partition by time which is 60 if (tail(name$Behavior, 1) == 1) {name-rbind(name, c(4, 3600))} else {name-rbind(name, c(1, 3600))} if (nrow(name) %% 2 != 0) {name -name[-c(nrow(name)), -c(nrow(name))]} q-NULL for (y in (1: nrow(name))) { if (y %% 2 != 0) q-c(q, (c(name$Time[y]:name$Time[y +1]))) } b-table(cut(q,x)) sum(b) [[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] Pairwise comparison matrix elements
On Apr 9, 2012, at 3:17 PM, Kerapi wrote: Hi!, I'm really hoping someone out there will be able to help me. I recently started my MSc dissertation on Population Projection Matrices, which has been going well until now. I am trying to set-up a general script that does a pairwise comparison of all elements in my matrices. So for example, given that I have the following matrix S: S [,1] [,2] [,3] [1,] 0.0 0.007361183 0.0001634936 [2,] 13.88458 0.353988378 0.00 [3,] 0.0 16.086367098 0.3572819557 I'd like to create a matrix SA for each element in my matrix compared to another element like so: SA S-[1,1] [,1][,2][,3] [1,] [1,1]-[1,1] [1,2]-[1,1] [1,3]-[1,1] [2,] [2,1]-[1,1] [2,2]-[1,1] [2,3]-[1,1] [3,] [3,1]-[1,1] [3,2]-[1,1] [3,3]-[1,1] Try this: kS - -1*kronecker(c(S),S, -) dim(kS) - c(3, 3, 3, 3) The row number is the thrid dimension and the col nimber is the second dimension. If you want them to be rearrange to be in the ( r,c) order you can use aperm. kSr - aperm(kS, c(1, 4, 2, 3)) #The zero entires were used because I could easily tell what a correct answer should be. kSr[,,3,1] [,1] [,2] [,3] [1,] 0.0 0.007361183 0.0001634936 [2,] 13.88458 0.353988378 0.00 [3,] 0.0 16.086367098 0.3572819557 kSr[,,1,1] [,1] [,2] [,3] [1,] 0.0 0.007361183 0.0001634936 [2,] 13.88458 0.353988378 0.00 [3,] 0.0 16.086367098 0.3572819557 kSr[,,2,3] [,1] [,2] [,3] [1,] 0.0 0.007361183 0.0001634936 [2,] 13.88458 0.353988378 0.00 [3,] 0.0 16.086367098 0.3572819557 Test a non-zero entry kSr[,,2,2] [,1] [,2] [,3] [1,] -0.3539884 -0.3466272 -0.353824884 [2,] 13.5305916 0.000 -0.353988378 [3,] -0.3539884 15.7323787 0.003293578 S - S[2,2] [,1] [,2] [,3] [1,] -0.3539884 -0.3466272 -0.353824884 [2,] 13.5305916 0.000 -0.353988378 [3,] -0.3539884 15.7323787 0.003293578 sapply(1:3, function(rw) sapply(1:3, function(cl) identical(S-S[rw,cl], kSr[,,rw,cl] )) ) [,1] [,2] [,3] [1,] TRUE TRUE TRUE [2,] TRUE TRUE TRUE [3,] TRUE TRUE TRUE -- David. S-[1,2] [...] S-[1,3] [...] ... and so on The aim is to be able to prove existing rules and trends across matrix dimensions. For example, I've been trying to test whether the first line of values decreases or remains the same from column 1 to the penultimate column: S[1,j] = S[1,lj] The last thing I tried was: M-Matlab2R([0 24.724 1377.48;0.029328 0.26 0;0 0.014 0.78]) w - abs(Re(eigen(M)$vectors[,1])) v - abs(Re(eigen(t(M))$vectors[,1])) w - w/sum(w) v - v/(t(v)%*%w) S - (v%*%t(w))/as.vector(t(v)%*%w) S n- length(S) M.comp - array(S,dim=c(n,n,n,n)) for (i in 1:n) { for (j in 1:n) { for (k in 1:n) { for (l in 1:n) { poscompa[1,3,1,j-1] - S[i,j]=S[k,l] # if=1, then TRUE and if=0, then FALSE M.comp sum(poscompa[1,3,1,1:n]) # should give me the value for n. I know this doesn't test it accurately, but first I need to get the loops to work right. I'm getting an array of errors such as can't find function poscompa (or compa), subscripts out of range or type of closure not valid regardless of what I try. Using COMPAR (poscompa and compa) was the last recommendation I got, but I'm starting to wonder if there might be other ways to go about this. All out-of-the-box ideas I've come up with and tried haven't gotten me much farther. I've now practically exhausted my creative thinking, and I'm becoming very frustrated. I'd really like to get this script going since my current one would make my life hell for large populations (60x60 population matrices). Any ideas on how I could move forward? Many, many thanks! Marta [[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, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] For loops
Hi, I am having trouble with syntax for a for loop. Here is what I am trying to do. class=c(rep(1,3),rep(2,3),rep(3,3)) out1=rnorm(length(class)) out2=rnorm(length(class)) out3=rnorm(length(class)) data=data.frame(class,out1,out2,out3) dat.split=split(data,data$class) for(i in 1:3){ sub[i]=dat.split[i] } However, the for loop doesn't work. I want to assign each split to a different data object. Better yet, how I could assign each class to a separate object and skip the splitting? Thanks, Chris [[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] For loops
class1 - data[1==data$class,] gets you a subset of the data into a dedicated object, but if you want to handle arbitrarily large amounts of data or class values then the list output of split is really much better to stay with. Also, data is a predefined function, so it is not a good idea to define objects with that name due to confusing code review. --- 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. Christopher Desjardins cddesjard...@gmail.com wrote: Hi, I am having trouble with syntax for a for loop. Here is what I am trying to do. class=c(rep(1,3),rep(2,3),rep(3,3)) out1=rnorm(length(class)) out2=rnorm(length(class)) out3=rnorm(length(class)) data=data.frame(class,out1,out2,out3) dat.split=split(data,data$class) for(i in 1:3){ sub[i]=dat.split[i] } However, the for loop doesn't work. I want to assign each split to a different data object. Better yet, how I could assign each class to a separate object and skip the splitting? Thanks, Chris [[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] For loops
On Apr 9, 2012, at 5:33 PM, Christopher Desjardins wrote: Hi, I am having trouble with syntax for a for loop. Here is what I am trying to do. class=c(rep(1,3),rep(2,3),rep(3,3)) out1=rnorm(length(class)) out2=rnorm(length(class)) out3=rnorm(length(class)) data=data.frame(class,out1,out2,out3) dat.split=split(data,data$class) for(i in 1:3){ sub[i]=dat.split[i] } However, the for loop doesn't work. I want to assign each split to a different data object. Why? What's wrong with leaving them in a list that split() provides. In that form you can easily use lapply() and start your 12 step program away from for-loop addiction. Better yet, how I could assign each class to a separate object and skip the splitting? -- David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Overall model significance for poisson GLM
Pat Wilkins pwilkin2 at illinois.edu writes: Greetings, I am running glm models for species counts using a poisson link function. Normal summary functions for this provide summary statistics in the form of the deviance, AIC, and p-values for individual predictors. I would like to obtain the p-value for the overall model. So far, I have been using an analysis of deviance table to check a model against the null model with the intercept as the only predictor. Any advice on other methods to obtain the proper p-value would be appreciated. What you're doing seems reasonable, although you can also dig the necessary values out of the summary and compute the p-value yourself: counts - c(18,17,15,20,10,20,25,13,12) outcome - gl(3,1,9) treatment - gl(3,3) d.AD - data.frame(treatment, outcome, counts) glm.D93 - glm(counts ~ outcome + treatment, family=poisson(), data=d.AD) ## as you have been doing anova(update(glm.D93,.~1),glm.D93,test=Chisq) ## or sg1 - summary(glm.D93) devdiff - with(sg1,null.deviance-deviance) dfdiff - with(sg1,df.null-df.residual) pchisq(abs(devdiff),df=dfdiff,lower.tail=FALSE) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Most efficient way to do this...
I have time-series data looking like this: dataIn[sample(c(1:nrow(dataIn)), 25),] accelerometer_y id data_block_epoch_time 782 0.8424 201300 133179733 1868 0.3432 202386 1331797384000 1828 0.3510 202346 1331797382000 1026 0.2184 201544 1331797342000 1569 0.3432 202087 1331797369000 1453 0.3588 201971 1331797363000 1204 0.3666 201722 1331797351000 1821 0.3588 202339 1331797382000 860 0.8658 201378 1331797333000 910 0.8580 201428 1331797336000 1488 0.3432 202006 1331797365000 578 0.9126 201096 1331797319000 1478 0.3666 201996 1331797364000 1183 0.3588 201701 133179735 29 -0.1716 200547 1331797292000 1540 0.3588 202058 1331797367000 392 -0.1560 200910 133179731 1533 0.3744 202051 1331797367000 1016 0.6318 201534 1331797341000 314 -0.1560 200832 1331797306000 410 -0.1638 200928 1331797311000 769 0.8580 201287 1331797329000 1101 0.3588 201619 1331797346000 403 -0.1638 200921 1331797311000 1794 0.3666 202312 133179738 The id field represents the subsecond value in the timestamp (n ids/second). What I'd like to do is combine the fields, so that I'm left with an equally-spaced timeseries with readings. What's the most efficient way to do this using R? There are an arbitrary number of timestamps per second. Many thanks! -- H -- Sent from my mobile device Envoyait de mon portable [[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 add 3d-points to bplot {rms} figure?
Hello! I have created a bplot-figure using this code: *file - 2dcali_red.ttt ux-as.matrix(read.table(file, dec = ,)) mode(ux)-'numeric' vel-ux[,1] ang-ux[,2] x-ux[,3] y-ux[,4] dat- data.frame(ang=ang, x=x,y=y) require(rms) ddist2 - datadist(dat) options(datadist=ddist2) fitn - lrm(ang ~ rcs(x,4) + rcs(y,4), data=dat) predi - Predict(fitn, x, y) bplot(predi,lfun=wireframe, screen = list(z = -40, x = -80), drape=TRUE)* The file 2dcali_red.ttt consists of 4 columns can be found here : http://www.color-space.de/upload/dl/2dcali_red.ttt link The code gives me the bplot-figure, which looks like this: http://r.789695.n4.nabble.com/file/n4544050/Screen_shot_2012-04-10_at_12.04.24_AM.png Now I want to add some 3d data (x,y,z) to this plot (like scatterplot3d), but i don't know how. I also don't know how to extract the x,y,z mesh as coordinates from the bplot. Some help would be very appreciated. Thank you! -- View this message in context: http://r.789695.n4.nabble.com/how-to-add-3d-points-to-bplot-rms-figure-tp4544050p4544050.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problems with CRAN package Ryacas
see below! On Mon, Apr 9, 2012 at 3:04 PM, Gabor Grothendieck ggrothendi...@gmail.com wrote: On Mon, Apr 9, 2012 at 4:02 PM, Gabor Grothendieck ggrothendi...@gmail.com wrote: On Mon, Apr 9, 2012 at 2:47 PM, Martin Maechler maech...@stat.math.ethz.ch wrote: Apropos: I don't have the problems, the OP had, but on my ubuntu notebook, Ryacas does not return expressions (just the strings), and hence as.expression( yacas-result ) always gives NULL and e.g. the demo(Ryacas-Function) also fails: yacas(expression(deriv(BurrCDF(x,c,k k*c*x^(c-1)*(x^c+1)^(-(k+1)); yy - yacas(expression(deriv(BurrCDF(x,c,k yy k*c*x^(c-1)*(x^c+1)^(-(k+1)); str(yy) List of 2 $ : NULL $ YacasForm: chr k*c*x^(c-1)*(x^c+1)^(-(k+1)); - attr(*, class)= chr yacas I have - yacas 1.2.2 (standard ubuntu package) The latest version of Ryacas on CRAN is 0.2-11 http://cran.r-project.org/package=Ryacas and the version of yacas supported by Ryacas is: yacas(Version()) expression(1.0.63) The troubleshooting section on the home page that I referred to does mention this. Thanks! I deinstalled the synaptic-installed yacas, downloaded yacas_1.0.63.tgz unpacked it , run ./configure --enable-server then did make But make does not terminate, it ends in error! I give the last few lines of output: rm -fr .libs/filescanner.la .libs/filescanner.* .libs/filescanner.* gcc -shared filescanner.lo plugin.lo -Wl,-soname -Wl,filescanner.so -o .libs/filescanner.so ar cru .libs/filescanner.a filescanner.o plugin.o ranlib .libs/filescanner.a creating filescanner.la (cd .libs rm -f filescanner.la ln -s ../filescanner.la filescanner.la) make[3]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/plugins/filescanner' make[3]: Entering directory `/home/kjetil/yacas/yacas-1.0.63/plugins' make[3]: Nothing to be done for `all-am'. make[3]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/plugins' make[2]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/plugins' Making all in proteus make[2]: Entering directory `/home/kjetil/yacas/yacas-1.0.63/proteus' make[2]: Nothing to be done for `all'. make[2]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/proteus' Making all in manmake make[2]: Entering directory `/home/kjetil/yacas/yacas-1.0.63/manmake' source='manripper.cpp' object='manripper.o' libtool=no \ depfile='.deps/manripper.Po' tmpdepfile='.deps/manripper.TPo' \ depmode=gcc3 /bin/bash ../depcomp \ g++ -DHAVE_CONFIG_H -I. -I. -I.. -g -O2 -Wall -c -o manripper.o `test -f 'manripper.cpp' || echo './'`manripper.cpp manripper.cpp: In function ‘void GetBf(char*, int, FILE*)’: manripper.cpp:20:21: error: ‘strlen’ was not declared in this scope manripper.cpp: In function ‘void ProcessFile(char*)’: manripper.cpp:40:23: error: ‘strlen’ was not declared in this scope manripper.cpp:43:30: error: ‘strncmp’ was not declared in this scope manripper.cpp:57:25: error: ‘memset’ was not declared in this scope manripper.cpp:58:39: error: ‘memcpy’ was not declared in this scope manripper.cpp:78:27: error: ‘strlen’ was not declared in this scope manripper.cpp:81:25: error: ‘memset’ was not declared in this scope manripper.cpp:82:41: error: ‘memcpy’ was not declared in this scope make[2]: *** [manripper.o] Error 1 make[2]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/manmake' make[1]: *** [all-recursive] Error 1 make[1]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63' make: *** [all] Error 2 kjetil@kjetil:~/yacas/yacas-1.0.63$ Kjetil -- Statistics Software Consulting GKX Group, GKX Associates Inc. tel: 1-877-GKX-GROUP email: ggrothendieck at gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Gradients in bar charts XXXX
On Mon, Apr 9, 2012 at 12:40 PM, Jason Rodriguez jason.rodrig...@dca.ga.gov wrote: Hello, I have a graphics-related question: I was wondering if anyone knows of a way to create a bar chart that is colored with a three-part gradient that changes at fixed y-values. Each bar needs to fade green-to-yellow at Y=.10 and from yellow-to-red at Y=.20. Is there an option in a package somewhere that offers an easy way to do this? ?rainbow ?hsv In R an easy way is an ill-defined term. In the absence of actual data: bpd - matrix(c(1,seq(0,1,l=64),2,1,seq(0,1,l=64),5,1,seq(0,1,l=64),7),nc=3) mycols - c('green',rainbow(64,start=0,end=.4)[64:1],'red') barplot(bpd,col=mycols,border=NA) Easy enough ? Cheers Attached is a chart I macgyvered together in Excel using a combination of a simple bar chart, fit line, and some drawing tools. I want to avoid doing it this way in the future by finding a way to replicate it in R. Any ideas? Thanks, Jason Michael Rodriguez Data Analyst State Housing Trust Fund for the Homeless Georgia Department of Community Affairs Email: jason.rodrig...@dca.ga.gov __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 convert seconds to 12 hour time format
Hello everyone, I am wondering if there is any routine in R which can convert time given in 'seconds' unit to the 12 hour time format. For example, suppose the data set looks like x=c(36885,84000,20) #x in seconds I want to get the output as [1] 11:14:45 AM [2] 11:20:00 PM [3] 12:20:00 AM Does anyone have any idea? Thanks in advance. Cassie [[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] Problems with CRAN package Ryacas
We have had Linux users that successfully used Ryacas but if you can't get yacas to build on your particular system then you could try running the Windows version under wine. A Windows binary of yacas is available so you won't have to build yacas. I don't know if anyone has tried that yet but its worth a try. Another possibility is to try rsympy or rmathpiper (see introductory paragraphs on Ryacas home page for links). On Mon, Apr 9, 2012 at 8:00 PM, Kjetil Halvorsen kjetilbrinchmannhalvor...@gmail.com wrote: see below! On Mon, Apr 9, 2012 at 3:04 PM, Gabor Grothendieck ggrothendi...@gmail.com wrote: On Mon, Apr 9, 2012 at 4:02 PM, Gabor Grothendieck ggrothendi...@gmail.com wrote: On Mon, Apr 9, 2012 at 2:47 PM, Martin Maechler maech...@stat.math.ethz.ch wrote: Apropos: I don't have the problems, the OP had, but on my ubuntu notebook, Ryacas does not return expressions (just the strings), and hence as.expression( yacas-result ) always gives NULL and e.g. the demo(Ryacas-Function) also fails: yacas(expression(deriv(BurrCDF(x,c,k k*c*x^(c-1)*(x^c+1)^(-(k+1)); yy - yacas(expression(deriv(BurrCDF(x,c,k yy k*c*x^(c-1)*(x^c+1)^(-(k+1)); str(yy) List of 2 $ : NULL $ YacasForm: chr k*c*x^(c-1)*(x^c+1)^(-(k+1)); - attr(*, class)= chr yacas I have - yacas 1.2.2 (standard ubuntu package) The latest version of Ryacas on CRAN is 0.2-11 http://cran.r-project.org/package=Ryacas and the version of yacas supported by Ryacas is: yacas(Version()) expression(1.0.63) The troubleshooting section on the home page that I referred to does mention this. Thanks! I deinstalled the synaptic-installed yacas, downloaded yacas_1.0.63.tgz unpacked it , run ./configure --enable-server then did make But make does not terminate, it ends in error! I give the last few lines of output: rm -fr .libs/filescanner.la .libs/filescanner.* .libs/filescanner.* gcc -shared filescanner.lo plugin.lo -Wl,-soname -Wl,filescanner.so -o .libs/filescanner.so ar cru .libs/filescanner.a filescanner.o plugin.o ranlib .libs/filescanner.a creating filescanner.la (cd .libs rm -f filescanner.la ln -s ../filescanner.la filescanner.la) make[3]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/plugins/filescanner' make[3]: Entering directory `/home/kjetil/yacas/yacas-1.0.63/plugins' make[3]: Nothing to be done for `all-am'. make[3]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/plugins' make[2]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/plugins' Making all in proteus make[2]: Entering directory `/home/kjetil/yacas/yacas-1.0.63/proteus' make[2]: Nothing to be done for `all'. make[2]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/proteus' Making all in manmake make[2]: Entering directory `/home/kjetil/yacas/yacas-1.0.63/manmake' source='manripper.cpp' object='manripper.o' libtool=no \ depfile='.deps/manripper.Po' tmpdepfile='.deps/manripper.TPo' \ depmode=gcc3 /bin/bash ../depcomp \ g++ -DHAVE_CONFIG_H -I. -I. -I.. -g -O2 -Wall -c -o manripper.o `test -f 'manripper.cpp' || echo './'`manripper.cpp manripper.cpp: In function ‘void GetBf(char*, int, FILE*)’: manripper.cpp:20:21: error: ‘strlen’ was not declared in this scope manripper.cpp: In function ‘void ProcessFile(char*)’: manripper.cpp:40:23: error: ‘strlen’ was not declared in this scope manripper.cpp:43:30: error: ‘strncmp’ was not declared in this scope manripper.cpp:57:25: error: ‘memset’ was not declared in this scope manripper.cpp:58:39: error: ‘memcpy’ was not declared in this scope manripper.cpp:78:27: error: ‘strlen’ was not declared in this scope manripper.cpp:81:25: error: ‘memset’ was not declared in this scope manripper.cpp:82:41: error: ‘memcpy’ was not declared in this scope make[2]: *** [manripper.o] Error 1 make[2]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/manmake' make[1]: *** [all-recursive] Error 1 make[1]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63' make: *** [all] Error 2 kjetil@kjetil:~/yacas/yacas-1.0.63$ Kjetil -- Statistics Software Consulting GKX Group, GKX Associates Inc. tel: 1-877-GKX-GROUP email: ggrothendieck at gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Statistics Software Consulting GKX Group, GKX Associates Inc. tel: 1-877-GKX-GROUP email: ggrothendieck at gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] quantmod getOptionChain Not Work
On a Unix-alike system, if you have an svn client installed, it suffices to type svn checkout svn://svn.r-forge.r-project.org/svnroot/quantmod/ at a command prompt -- then use R CMD install on the pkg/ directory to install the new version. No idea how to do it on Windows Michael On Sun, Apr 8, 2012 at 7:11 PM, Sparks, John James jspa...@uic.edu wrote: Michael, I have not had time to look at this for a while but still wanted to say thanks for looking into it and sending this solution. By the way, Jeff mentioned that the version of quantmod on the SVN (0.3.18) works for this. I tried to figure out how to download that version, but found the documentation on SVN's quite confusing. Is there anyway that you could make that version available? Much appreciated. --John Sparks On Fri, March 23, 2012 5:55 pm, R. Michael Weylandt wrote: Sorry about that: two small mistakes and I imagine there are a few more I've missed. This should actually work: ### library(XML) readYahooOptions - function(Symbols, Exp, ...){ parse.expiry - function(x) { if(is.null(x)) return(NULL) if(inherits(x, Date) || inherits(x, POSIXt)) return(format(x, %Y-%m)) if (nchar(x) == 5L) { x - sprintf(substring(x, 4, 5), match(substring(x, 1, 3), month.abb), fmt = 20%s-%02i) } else if (nchar(x) == 6L) { x - paste(substring(x, 1, 4), substring(x, 5, 6), sep = -) } return(x) } clean.opt.table - function(tableIn){ tableOut - sapply(tableIn[,-2], function(x) as.numeric(gsub(,,,x))) rownames(tableOut) - tableIn[,2] tableOut } if(missing(Exp)) optURL - paste(paste(http://finance.yahoo.com/q/op?s,Symbols,sep==;),Options,sep=+) else optURL - paste(paste(http://finance.yahoo.com/q/op?s=,Symbols,m=,parse.expiry(Exp),sep=),Options,sep=+) if(!missing(Exp) is.null(Exp)) { optPage - readLines(optURL) optPage - optPage[grep(View By Expiration, optPage)] allExp - gregexpr(m=, optPage)[[1]][-1] + 2 allExp - substring(optPage, allExp, allExp + 6) allExp - allExp[seq_len(length(allExp)-1)] # Last one seems useless ? return(structure(lapply(allExp, readYahooOptions, Symbols=Symbols), .Names=format(as.yearmon(allExp } stopifnot(require(XML)) optURL - readHTMLTable(optURL) # Not smart to hard code these but it's a 'good-enough' hack for now # Also, what is table 9 on this page? list(calls = clean.opt.table(optURL[[10]]), puts = clean.opt.table(optURL[[14]]), symbol = Symbols) } On Fri, Mar 23, 2012 at 6:44 PM, R. Michael Weylandt michael.weyla...@gmail.com wrote: I just got around to taking a look at this, but below is a fix. It seems like yahoo finance redesigned the page and rather than reparsing all their HTML, I'll use Duncan TL's XML package to make life happier. (I loathe HTML parsing) This isn't thoroughly tested and it'll break if yahoo redesigns things again (I hardcode the table numbers for now) but it seems to work well enough. Let me know if you have any errors with it. If Jeff likes it, it should be a drop-in replacement for the getOptionChain.yahoo for quantmod with a name change. Feedback welcome, Michael # library(XML) readYahooOptions - function(Symbols, Exp, ...){ áparse.expiry - function(x) { á áif(is.null(x)) á á áreturn(NULL) á áif(inherits(x, Date) || inherits(x, POSIXt)) á á áreturn(format(x, %Y-%m)) á áif (nchar(x) == 5L) { á á áx - sprintf(substring(x, 4, 5), match(substring(x, á á á á á á á á á á á á á á á á á á á á á á á á á á á 1, 3), month.abb), fmt = 20%s-%02i) á á} á áelse if (nchar(x) == 6L) { á á áx - paste(substring(x, 1, 4), substring(x, 5, 6), á á á á á á á á sep = -) á á} á áreturn(x) á} áclean.opt.table - function(tableIn){ á átableOut - lapply(tableIn[,-2], function(x) as.numeric(gsub(,,,x))) á árownames(tableOut) - tableIn[,2] á} áif(missing(Exp)) á áoptURL - paste(paste(http://finance.yahoo.com/q/op?s,Symbols,sep==;),Options,sep=+) áelse á áoptURL - paste(paste(http://finance.yahoo.com/q/op?s=,Symbols,m=,parse.expiry(Exp),sep=),Options,sep=+) áif(!missing(Exp) is.null(Exp)) { á áoptPage - readLines(optURL) á áoptPage - optPage[grep(View By Expiration, optPage)] á áallExp - gregexpr(m=, optPage)[[1]][-1] + 2 á áallExp - substring(optPage, allExp, allExp + 6) á áallExp - allExp[seq_len(length(allExp)-1)] # Last one seems useless ? Always true? á áreturn(structure(lapply(allExp, readYahooOptions, Symbols=Symbols), .Names=format(as.yearmon(allExp á} ástopifnot(require(XML)) áoptURL - readHTMLTable(optURL) á# Not smart to hard code these but it's a 'good-enough' hack for now á# Also, what is table 9 on this page? áCALLS - optURL[[10]] áPUTS - optURL[[14]] álist(calls =
Re: [R] how to convert seconds to 12 hour time format
Try this: Sys.setlocale(category = LC_TIME, locale = US) format(as.POSIXct('0001-01-01 00:00:00') + x, %I:%M:%S %p) On Mon, Apr 9, 2012 at 9:38 PM, cassie jones cassiejone...@gmail.com wrote: Hello everyone, I am wondering if there is any routine in R which can convert time given in 'seconds' unit to the 12 hour time format. For example, suppose the data set looks like x=c(36885,84000,20) #x in seconds I want to get the output as [1] 11:14:45 AM [2] 11:20:00 PM [3] 12:20:00 AM Does anyone have any idea? Thanks in advance. Cassie [[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. -- Henrique Dallazuanna Curitiba-Paraná-Brasil 25° 25' 40 S 49° 16' 22 O __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 convert seconds to 12 hour time format
Thanks Henrique.It worked!!! I appreciate your help. Cassie On Mon, Apr 9, 2012 at 8:06 PM, Henrique Dallazuanna www...@gmail.comwrote: Try this: Sys.setlocale(category = LC_TIME, locale = US) format(as.POSIXct('0001-01-01 00:00:00') + x, %I:%M:%S %p) On Mon, Apr 9, 2012 at 9:38 PM, cassie jones cassiejone...@gmail.com wrote: Hello everyone, I am wondering if there is any routine in R which can convert time given in 'seconds' unit to the 12 hour time format. For example, suppose the data set looks like x=c(36885,84000,20) #x in seconds I want to get the output as [1] 11:14:45 AM [2] 11:20:00 PM [3] 12:20:00 AM Does anyone have any idea? Thanks in advance. Cassie [[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. -- Henrique Dallazuanna Curitiba-Paraná-Brasil 25° 25' 40 S 49° 16' 22 O [[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] Problems with CRAN package Ryacas
see below! On Mon, Apr 9, 2012 at 7:38 PM, Gabor Grothendieck ggrothendi...@gmail.com wrote: We have had Linux users that successfully used Ryacas but if you can't get yacas to build on your particular system then you could try running the Windows version under wine. A Windows binary of yacas is available so you won't have to build yacas. I don't know if anyone has tried that yet but its worth a try. Another possibility is to try rsympy or rmathpiper (see introductory paragraphs on Ryacas home page for links). On Mon, Apr 9, 2012 at 8:00 PM, Kjetil Halvorsen kjetilbrinchmannhalvor...@gmail.com wrote: see below! On Mon, Apr 9, 2012 at 3:04 PM, Gabor Grothendieck ggrothendi...@gmail.com wrote: On Mon, Apr 9, 2012 at 4:02 PM, Gabor Grothendieck ggrothendi...@gmail.com wrote: On Mon, Apr 9, 2012 at 2:47 PM, Martin Maechler maech...@stat.math.ethz.ch wrote: Apropos: I don't have the problems, the OP had, but on my ubuntu notebook, Ryacas does not return expressions (just the strings), and hence as.expression( yacas-result ) always gives NULL and e.g. the demo(Ryacas-Function) also fails: yacas(expression(deriv(BurrCDF(x,c,k k*c*x^(c-1)*(x^c+1)^(-(k+1)); yy - yacas(expression(deriv(BurrCDF(x,c,k yy k*c*x^(c-1)*(x^c+1)^(-(k+1)); str(yy) List of 2 $ : NULL $ YacasForm: chr k*c*x^(c-1)*(x^c+1)^(-(k+1)); - attr(*, class)= chr yacas I have - yacas 1.2.2 (standard ubuntu package) The latest version of Ryacas on CRAN is 0.2-11 http://cran.r-project.org/package=Ryacas and the version of yacas supported by Ryacas is: yacas(Version()) expression(1.0.63) The troubleshooting section on the home page that I referred to does mention this. I am trying out rsympy by now, and that is working just fine! Thanks. Kjetil Thanks! I deinstalled the synaptic-installed yacas, downloaded yacas_1.0.63.tgz unpacked it , run ./configure --enable-server then did make But make does not terminate, it ends in error! I give the last few lines of output: rm -fr .libs/filescanner.la .libs/filescanner.* .libs/filescanner.* gcc -shared filescanner.lo plugin.lo -Wl,-soname -Wl,filescanner.so -o .libs/filescanner.so ar cru .libs/filescanner.a filescanner.o plugin.o ranlib .libs/filescanner.a creating filescanner.la (cd .libs rm -f filescanner.la ln -s ../filescanner.la filescanner.la) make[3]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/plugins/filescanner' make[3]: Entering directory `/home/kjetil/yacas/yacas-1.0.63/plugins' make[3]: Nothing to be done for `all-am'. make[3]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/plugins' make[2]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/plugins' Making all in proteus make[2]: Entering directory `/home/kjetil/yacas/yacas-1.0.63/proteus' make[2]: Nothing to be done for `all'. make[2]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/proteus' Making all in manmake make[2]: Entering directory `/home/kjetil/yacas/yacas-1.0.63/manmake' source='manripper.cpp' object='manripper.o' libtool=no \ depfile='.deps/manripper.Po' tmpdepfile='.deps/manripper.TPo' \ depmode=gcc3 /bin/bash ../depcomp \ g++ -DHAVE_CONFIG_H -I. -I. -I.. -g -O2 -Wall -c -o manripper.o `test -f 'manripper.cpp' || echo './'`manripper.cpp manripper.cpp: In function ‘void GetBf(char*, int, FILE*)’: manripper.cpp:20:21: error: ‘strlen’ was not declared in this scope manripper.cpp: In function ‘void ProcessFile(char*)’: manripper.cpp:40:23: error: ‘strlen’ was not declared in this scope manripper.cpp:43:30: error: ‘strncmp’ was not declared in this scope manripper.cpp:57:25: error: ‘memset’ was not declared in this scope manripper.cpp:58:39: error: ‘memcpy’ was not declared in this scope manripper.cpp:78:27: error: ‘strlen’ was not declared in this scope manripper.cpp:81:25: error: ‘memset’ was not declared in this scope manripper.cpp:82:41: error: ‘memcpy’ was not declared in this scope make[2]: *** [manripper.o] Error 1 make[2]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63/manmake' make[1]: *** [all-recursive] Error 1 make[1]: Leaving directory `/home/kjetil/yacas/yacas-1.0.63' make: *** [all] Error 2 kjetil@kjetil:~/yacas/yacas-1.0.63$ Kjetil -- Statistics Software Consulting GKX Group, GKX Associates Inc. tel: 1-877-GKX-GROUP email: ggrothendieck at gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Statistics Software Consulting GKX Group, GKX Associates Inc. tel: 1-877-GKX-GROUP email: ggrothendieck at gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide
Re: [R] binned tabulation
Hello, Jean V Adams wrote Delia, name - data.frame(Behavior=c(1, 2, 1, 2, 1), Time=c(0, 40, 45, 55, 57)) appear - name$Time[name$Behavior==1] disappear - name$Time[name$Behavior==2] if(length(appear) length(disappear)) disappear - c(disappear, 60) sum(disappear - appear) Delia Shelton wrote on 04/09/2012 12:30:23 PM: Hi, I am attempting to tabulate binned data. The '1' represents the appearance of the focal mouse pup, and '2' represents the disappearance of the focal mouse pup. The code written below is intended to calculate the total time spent appeared out of 3600s. For Sample 1, both the hand calculation and R code yield the same result, 50. A problem seems to occur when '1' is the last entry. For Sample 2, the total time appeared is 53 (hand calculation), however, using the R code below yields 55. If you have any suggestions for solving the problem, please let me know. Thank you in advance for any assistance you may provide. Delia Sample 1 0.0 1 40 2 45 1 55 2 Sample 2 0.0 1 40 2 45 1 55 2 57 1 name = read.table(file.choose(),header=F) # opening a data file colnames(name)-c(Time, Behavior) name = data.frame(name$Behavior, name$Time) colnames(name)-c(Behavior, Time) name-name[name$Time 3600, ]; ## run file x-seq(0,3600, by = 60) # total time partition by time which is 60 if (tail(name$Behavior, 1) == 1) {name-rbind(name, c(4, 3600))} else {name-rbind(name, c(1, 3600))} if (nrow(name) %% 2 != 0) {name -name[-c(nrow(name)), -c(nrow(name))]} q-NULL for (y in (1: nrow(name))) { if (y %% 2 != 0) q-c(q, (c(name$Time[y]:name$Time[y +1]))) } b-table(cut(q,x)) sum(b) [[alternative HTML version deleted]] __ R-help@ mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Or use a logical index into the data.frame's rows. fun.count - function(x, increment=60){ if (x$Behavior[nrow(x)] == 1) x - rbind(x, c(4, increment*ceiling(x$Time[nrow(x)] / increment))) inx - x$Behavior == 1 sum(x$Time[!inx] - x$Time[inx]) } fun.count(name) Hope this helps, Rui Barradas -- View this message in context: http://r.789695.n4.nabble.com/binned-tabulation-tp4543405p4544257.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Panel.abline would not show beyond a certain slope value
Yup. Latest version 2.15.0 for windows solved the problem alright! Thank you Uwe. Regards, Fayez From: Uwe Ligges [lig...@statistik.tu-dortmund.de] Sent: Monday, April 09, 2012 9:56 AM To: Aziz, Muhammad Fayez Cc: r-help@r-project.org Subject: Re: [R] Panel.abline would not show beyond a certain slope value On 09.04.2012 06:31, Aziz, Muhammad Fayez wrote: Hi, Please find the code and data following. Problem appears in lines 37 - 42 of the code. I am using R.2.13.0 on WinXP. Works for me with a recent version of R (i.e. 2.15.0) / grid and lattice (i.e. 0.20-6). Please always try with recent versions of R and packages before asking. Best, Uwe Ligges Regards, Fayez PowerLawGraphsAblineQn.Râ library('reshape') library(lattice) library(igraph) # to use power.law.fit function library(latticeExtra) # to use panel.lmlineq in loglog xyplot File2Open = C:\\Documents and Settings\\All Users\\Documents\\Desktop\\Fayez\\RPractice\\PowerLawGraphsAblineQnData.txt DataTable = read.table(File2Open, header = TRUE, sep = \t) md- melt(DataTable, id.vars = c('Range', 'TheValue', 'TheNo')) # removed 'DegType' column 030312 md$variable- as.numeric(substr(md$variable, 9, nchar(as.character(md$variable mypanel4loglog- function(x, # x is the variable column in melted data, equals the Domain Nos y, # y is the value column in melted data, the degrees ... # Rest of the arguments ) { kfreq- table(y); # compute frquency hash table of y, the values k- 1:max(y) for (i in k) { ichar- as.character(i) # convert to match the names(freq), the character-based hash key of freq, which is the value if (!(ichar %in% names(kfreq))) kfreq[ichar]- 0 } sortedkeys- as.character(k) kfreq- kfreq[sortedkeys] pk- kfreq / length(y) logk- log(k) logpk- log(pk) logpk[logpk == -Inf] = # remove the -Inf or log(p(k)) = 0 values for lm function, NULL is 0-length so use instead that has length of one null character logpk- as.numeric(logpk) # converts all values to character, lm needs numeric print(rbind(logk, logpk))#write.table(rbind(k, kfreq, pk, logk, logpk), paste(FilePath, \\data, sep=), sep = \t, append = TRUE) panel.xyplot(col=blue, logk, logpk, type = c('p', 'r')) lm.r = lm(logpk ~ logk) panel.abline(coef=c(-4.847634, -1.037480)) # --- This gets drawn panel.abline(coef=c(-4.847634, -1.037481)) # --- This doesn't get drawn print(coef(lm.r)) # -4.847634 -1.349699 ; -3.377894 -1.498693 } # end mypanel4loglog myprepanel4loglog- function(x, # x is the variable column in melted data, equals the Network's ages y, # y is the value column in melted data, the degrees ... # Rest of the arguments ) { FreqTable- as.data.frame(table(y)) FreqsVector- sort(FreqTable$Freq) Min- FreqsVector[1] # first element - the lowest value frequency #print(c(Max2, length(y))) list(ylim = c(log(Min / length(y)), 0), xlim = c(0, log(max(y # log(p(k)) is always -ve as p(k) is decimal, so max(log(p(k)) is 0 } # end myprepanel4loglog print(xyplot(value ~ variable | Range, data = md, xlab = log(k); Panel = Range, ylab = log(p(k)), main = log(k) vs. log(p(k)), groups = Range, pch = 20, # dots instead of circles panel = mypanel4loglog, prepanel = myprepanel4loglog, # to set the scale of k and pk scales = list(x = (relation = free), y = (relation = free)), # necessary to make x-axis in each panel adjustable according to k )) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 do Sweave users collaborate with Word users?
See http://biostat.mc.vanderbilt.edu/SweaveConvert for several other approaches. The most solid way I've found is to convert PDF to Word. Frank Alexander Shenkin wrote Thanks for the heads up, Paul. That's good to know. I happen to be in academics, as are my collaborators, but I could imagine running into problems down the road. The good thing seems to be that it's just the users who want to interact with R who need the software. If collaborators are just touching the text, or are inserting graphics of their own, then sounds like things are still good. But, that's just my impression from the brief description I've seen online. Thanks, Allie On 4/9/2012 7:42 AM, Paul Bivand wrote: If you're considering SWord, you should remember that the licence is not the normal R licence and in commercial use will require a commercial licence. While some academic disciplines use Word etc, the issue may be more common outside academia. For those of us where such requirements involve a procurement process, the need to purchase something when other users (and the administration) are happy with their Word/Excel solutions may be an insurmountable barrier. Good luck Paul Bivand Centre for Economic and Social Inclusion (a non-profit organisation) London On 9 April 2012 13:23, Richard M. Heiberger lt;rmh@gt; wrote: You might want to consider SWord, which provides similar facilities for the Word and R user. Word-oriented co-authors can modify the Word part of the document without impacting the R part of the document. SWord is by Thomas Baier thomas@, author of the statconnDCOM interface that is underneath RExcel. See rcom.univie.ac.at for information and download and to sign up on the rcom email list. Rich On Sat, Apr 7, 2012 at 4:54 PM, Alexander Shenkin lt;ashenkin@gt; wrote: Hello All, I'm getting my workflow switched over to Sweave, which is very cool. However, I collaborate with folks (as many of you must as well) who use Word to Track Changes amongst a group while crafting a paper. In the simplest case, there will just be two people (one Sweave user and one Word user) editing a paper. I'm wondering, how do Sweave users go about this? I could convert a sweave file to a .docx easily enough via an intermediary pdf, rtf, html or otherwise. However, once the file has been marked up with changes, the challenge is to migrate those (accepted) changes back to the sweave document. Perhaps the most straightforward way is to manually back-propagate changes, but I imagine that could be a painstaking process. Ideally, I imagine a tool that puts invisible tags in the word document when it is originally produced from Sweave, and is then able to propagate changes back to that sweave file after markup. I'd be pleasantly surprised if such a tool existed. Perhaps there are other ways of making this work. Any thoughts are kindly appreciated! Thanks, Allie __ R-help@ mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.htmllt;http://www.r-project.org/posting-guide.htmlgt; and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@ mailing list https://stat.ethz.ch/mailman/listinfo/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@ mailing list https://stat.ethz.ch/mailman/listinfo/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@ mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. - Frank Harrell Department of Biostatistics, Vanderbilt University -- View this message in context: http://r.789695.n4.nabble.com/How-do-Sweave-users-collaborate-with-Word-users-tp4540061p4544413.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to add 3d-points to bplot {rms} figure?
On Apr 9, 2012, at 6:24 PM, phi771 wrote: Hello! I have created a bplot-figure using this code: *file - 2dcali_red.ttt ux-as.matrix(read.table(file, dec = ,)) mode(ux)-'numeric' vel-ux[,1] ang-ux[,2] x-ux[,3] y-ux[,4] dat- data.frame(ang=ang, x=x,y=y) require(rms) ddist2 - datadist(dat) options(datadist=ddist2) fitn - lrm(ang ~ rcs(x,4) + rcs(y,4), data=dat) predi - Predict(fitn, x, y) bplot(predi,lfun=wireframe, screen = list(z = -40, x = -80), drape=TRUE)* The file 2dcali_red.ttt consists of 4 columns can be found here : http://www.color-space.de/upload/dl/2dcali_red.ttt link The code gives me the bplot-figure, which looks like this: http://r.789695.n4.nabble.com/file/n4544050/Screen_shot_2012-04-10_at_12.04.24_AM.png Now I want to add some 3d data (x,y,z) to this plot (like scatterplot3d), but i don't know how. I also don't know how to extract the x,y,z mesh as coordinates from the bplot. Some help would be very appreciated. Thank you! That is most probably a lattice::wireframe object. You could use 'trellis.focus' or the 'layer' functions in pkg::latticeExtra to attempt panel.3dbars() additions. I haven't actually done this, but I have added panel function results to the lattice::contourplots that are created byrms/Hmisc. The 3d aspect adds a significant level of extra work: http://tolstoy.newcastle.edu.au/R/e2/help/06/10/3274.html Maybe a contourplot would suffice? -- David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] How to get the SS and MS from oneway.test?
Hello everyone: I'm a new member of this group. I have a question about oneway.test. When I use anova(lm()) to analysis the ANOVA, I can get the information about Sum Sq and Mean Sq. (The R code and the results are as follows.) anova(lm(BackCalac~factor(Assay),data=Control)) Analysis of Variance Table Response: BackCalac Df Sum Sq Mean Sq F valuePr(F) factor(Assay) 4 270.846 67.711 56.219 1.345e-10 *** Residuals 20 24.088 1.204 But it default the variances are the same. If the variances aren't equal. I need to use the oneway.test method. Because oneway.test has the option about var.equal=F. Here, I have a question about oneway.test, How can I get SS, and MS information from oneway.test? My R code and the results are as follows. Thank you very much. :) oneway.test(BackCalac~factor(Assay), var.equal=T,data=Control) One-way analysis of means data: BackCalac and factor(Assay) F = 56.2191, num df = 4, denom df = 20, p-value = 1.345e-10 oneway.test(BackCalac~factor(Assay), var.equal=F,data=Control) One-way analysis of means (not assuming equal variances) data: BackCalac and factor(Assay) F = 92.8834, num df = 4.000, denom df = 9.625, p-value = 1.165e-07 -- View this message in context: http://r.789695.n4.nabble.com/How-to-get-the-SS-and-MS-from-oneway-test-tp4544417p4544417.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to add 3d-points to bplot {rms} figure?
On Apr 9, 2012, at 10:39 PM, David Winsemius wrote: On Apr 9, 2012, at 6:24 PM, phi771 wrote: Hello! I have created a bplot-figure using this code: *file - 2dcali_red.ttt ux-as.matrix(read.table(file, dec = ,)) mode(ux)-'numeric' vel-ux[,1] ang-ux[,2] x-ux[,3] y-ux[,4] dat- data.frame(ang=ang, x=x,y=y) require(rms) ddist2 - datadist(dat) options(datadist=ddist2) fitn - lrm(ang ~ rcs(x,4) + rcs(y,4), data=dat) predi - Predict(fitn, x, y) bplot(predi,lfun=wireframe, screen = list(z = -40, x = -80), drape=TRUE)* The file 2dcali_red.ttt consists of 4 columns can be found here : http://www.color-space.de/upload/dl/2dcali_red.ttt link The code gives me the bplot-figure, which looks like this: http://r.789695.n4.nabble.com/file/n4544050/Screen_shot_2012-04-10_at_12.04.24_AM.png Now I want to add some 3d data (x,y,z) to this plot (like scatterplot3d), but i don't know how. I also don't know how to extract the x,y,z mesh as coordinates from the bplot. Some help would be very appreciated. Thank you! That is most probably a lattice::wireframe object. You could use 'trellis.focus' or the 'layer' functions in pkg::latticeExtra to attempt panel.3dbars() additions. I haven't actually done this, but I have added panel function results to the lattice::contourplots that are created byrms/Hmisc. The 3d aspect adds a significant level of extra work: http://tolstoy.newcastle.edu.au/R/e2/help/06/10/3274.html There is also a panel.3d.identify and panel.identify.cloud that might be needed: https://r-forge.r-project.org/scm/viewvc.php/pkg/R/interaction.R?view=markuprevision=675root=lattice My efforts at adding bars to that plot have failed so far. This does not do anything useful: trellis.focus(panel, 1, 1) panel.3dbars(x=-1.5, y=-.1, z=0,distance=20, col=red, xlim=c(-2.4,-1), ylim=c(-0.15, 0.15), zlim=c(-50,100), xlim.scaled=c(-0.5,0.5), ylim.scaled=c(-0.5,0.5), zlim.scaled=c(-0.5,0.5), zero.scaled=25) trellis.unfocus() Maybe a contourplot would suffice? David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to add 3d-points to bplot {rms} figure?
Below. -- Bert On Mon, Apr 9, 2012 at 7:39 PM, David Winsemius dwinsem...@comcast.net wrote: On Apr 9, 2012, at 6:24 PM, phi771 wrote: Hello! I have created a bplot-figure using this code: *file - 2dcali_red.ttt ux-as.matrix(read.table(file, dec = ,)) mode(ux)-'numeric' vel-ux[,1] ang-ux[,2] x-ux[,3] y-ux[,4] dat- data.frame(ang=ang, x=x,y=y) require(rms) ddist2 - datadist(dat) options(datadist=ddist2) fitn - lrm(ang ~ rcs(x,4) + rcs(y,4), data=dat) predi - Predict(fitn, x, y) bplot(predi,lfun=wireframe, screen = list(z = -40, x = -80), drape=TRUE)* The file 2dcali_red.ttt consists of 4 columns can be found here : http://www.color-space.de/upload/dl/2dcali_red.ttt link The code gives me the bplot-figure, which looks like this: http://r.789695.n4.nabble.com/file/n4544050/Screen_shot_2012-04-10_at_12.04.24_AM.png Now I want to add some 3d data (x,y,z) to this plot (like scatterplot3d), but i don't know how. I also don't know how to extract the x,y,z mesh as coordinates from the bplot. Some help would be very appreciated. Thank you! That is most probably a lattice::wireframe object. You could use 'trellis.focus' or the 'layer' functions in pkg::latticeExtra to attempt panel.3dbars() additions. I haven't actually done this, but I have added panel function results to the lattice::contourplots that are created byrms/Hmisc. The 3d aspect adds a significant level of extra work: And is generally a bad idea. Perspective plots of all types are much over-hyped. They are pretty, but the very act of creating artificial perspective for non-realistic shapes can distort perception and lead to misinterpretation of data. Parts of the scene hide other parts; it is difficult to quantitatively judge depths of valleys and heights of hills; etc. And, as David notes, it is very difficult to get right. http://tolstoy.newcastle.edu.au/R/e2/help/06/10/3274.html Maybe a contourplot would suffice? -- and is generally a much better, albeit less sexy, idea. -- David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- 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 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 get the SS and MS from oneway.test?
It is not clear to me how meaningful sums of squares and mean squares are in the context of oneway.test.(). Be that as it may, the short answer is you ***can't*** get at them from the output of oneway.test(). The oneway.test() function simply does not return enough information. If you ***really*** want/need (or think you do!) SS and MS information, you'll simply have to hack the code of oneway.test(). The code looks simple enough. I would advise reading the cited 1951 Biometrics paper by Welch in conjunction with your code-hacking. Good luck. cheers, Rolf Turner On 10/04/12 14:16, lulumecindy wrote: Hello everyone: I'm a new member of this group. I have a question about oneway.test. When I use anova(lm()) to analysis the ANOVA, I can get the information about Sum Sq and Mean Sq. (The R code and the results are as follows.) anova(lm(BackCalac~factor(Assay),data=Control)) Analysis of Variance Table Response: BackCalac Df Sum Sq Mean Sq F valuePr(F) factor(Assay) 4 270.846 67.711 56.219 1.345e-10 *** Residuals 20 24.088 1.204 But it default the variances are the same. If the variances aren't equal. I need to use the oneway.test method. Because oneway.test has the option about var.equal=F. Here, I have a question about oneway.test, How can I get SS, and MS information from oneway.test? My R code and the results are as follows. Thank you very much. :) oneway.test(BackCalac~factor(Assay), var.equal=T,data=Control) One-way analysis of means data: BackCalac and factor(Assay) F = 56.2191, num df = 4, denom df = 20, p-value = 1.345e-10 oneway.test(BackCalac~factor(Assay), var.equal=F,data=Control) One-way analysis of means (not assuming equal variances) data: BackCalac and factor(Assay) F = 92.8834, num df = 4.000, denom df = 9.625, p-value = 1.165e-07 -- View this message in context: http://r.789695.n4.nabble.com/How-to-get-the-SS-and-MS-from-oneway-test-tp4544417p4544417.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to get the SS and MS from oneway.test?
From what I can see the function does not give back this information (either in the print output or the output object structure itself). The thing is that I have never studied the statistic used for this test, so I am not sure how it works and if giving the SS MS is meaningful or not. But if you want to dig in, simply write in R the command: oneway.test The code that seems to be relevant is near the end: } else { m - sum(w.i * m.i)/sum.w.i STATISTIC - sum(w.i * (m.i - m)^2)/((k - 1) * (1 + 2 * (k - 2) * tmp)) PARAMETER - c(k - 1, 1/(3 * tmp)) PVAL - pf(STATISTIC, k - 1, 1/(3 * tmp), lower.tail = FALSE) METHOD - paste(METHOD, (not assuming equal variances)) } So you could simply copy paste the full code of the function and edit it to extract this information... Contact Details:--- Contact me: tal.gal...@gmail.com | 972-52-7275845 Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) | www.r-statistics.com (English) -- On Tue, Apr 10, 2012 at 5:16 AM, lulumecindy cindy852...@gmail.com wrote: Hello everyone: I'm a new member of this group. I have a question about oneway.test. When I use anova(lm()) to analysis the ANOVA, I can get the information about Sum Sq and Mean Sq. (The R code and the results are as follows.) anova(lm(BackCalac~factor(Assay),data=Control)) Analysis of Variance Table Response: BackCalac Df Sum Sq Mean Sq F valuePr(F) factor(Assay) 4 270.846 67.711 56.219 1.345e-10 *** Residuals 20 24.088 1.204 But it default the variances are the same. If the variances aren't equal. I need to use the oneway.test method. Because oneway.test has the option about var.equal=F. Here, I have a question about oneway.test, How can I get SS, and MS information from oneway.test? My R code and the results are as follows. Thank you very much. :) oneway.test(BackCalac~factor(Assay), var.equal=T,data=Control) One-way analysis of means data: BackCalac and factor(Assay) F = 56.2191, num df = 4, denom df = 20, p-value = 1.345e-10 oneway.test(BackCalac~factor(Assay), var.equal=F,data=Control) One-way analysis of means (not assuming equal variances) data: BackCalac and factor(Assay) F = 92.8834, num df = 4.000, denom df = 9.625, p-value = 1.165e-07 -- View this message in context: http://r.789695.n4.nabble.com/How-to-get-the-SS-and-MS-from-oneway-test-tp4544417p4544417.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.