[R] factors and ops
Dear R People I have a variable which is an ID number that is a factor. I would like to look for ID numberbs that are greeater than a particular value. However, factors are a big ugly for these operations. I was messing with the HR data set from the SASmixed package. HR$Patient is a factor like that. Now if I used: subset(HR,as.integer(as.character(Patient)) 214) that will work. If seems to me that there may be a better way. Is that true? R Version 2.3.0 Windows Thanks in advance! sincerely, Erin Hodgess Associate Professor Department of Computer and Mathematical Sciences University of Houston - Downtown mailto: [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Is there a way to draw 3d plot?
Hi all, I have a 2D matrix, which has 100 rows, and 100 columns, I have a 2D matrix, with 100 rows and 100 columns, I want to display it using 3D plot, much like plot3d and mesh/surf functions in matlab. Specifically, in matlab, I just need to do the following: [X, Y]=meshgrid([0:0.01:0.99, 0:0.01:0.99]); % Z is my 2D matrix, surf(X, Y, Z); Note that X and Y are created so that I can associate physical meaning onto the x and y axis of the 3D plot. For example, my 100 rows represent 0, 0.01, 0.02, ... 0.99 here. In Matlab I can also drag in the graphic window and see from different visual angle and perspective of the 3D plot... Are there similar functions in R that (1) show 3D plot; (2) let me manipulate view angles easily? Thanks a lot! [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] factors and ops
Hi On 26 May 2006 at 1:33, Erin Hodgess wrote: Date sent: Fri, 26 May 2006 01:33:34 -0500 From: Erin Hodgess [EMAIL PROTECTED] To: r-help@stat.math.ethz.ch Subject:[R] factors and ops Dear R People I have a variable which is an ID number that is a factor. I would like to look for ID numberbs that are greeater than a particular value. However, factors are a big ugly for these operations. I was messing with the HR data set from the SASmixed package. HR$Patient is a factor like that. Now if I used: subset(HR,as.integer(as.character(Patient)) 214) that will work. AFAIK no other way. If you know you do not want HR$Patient as factor you can use HR$Patient - as.integer(as.character(HR$Patient)) or HR$Patient - as.numeric(as.character(HR$Patient)) and then subset(HR,Patient 214) but beware of attaching data.frame as if it is attached you need to detach it first and than to change factor to numeric and only then to attach it again. HTH Petr If seems to me that there may be a better way. Is that true? R Version 2.3.0 Windows Thanks in advance! sincerely, Erin Hodgess Associate Professor Department of Computer and Mathematical Sciences University of Houston - Downtown mailto: [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Petr Pikal [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Distribution Fitting
Hi, I know this is a bit off-topic, but I am quite puzzled. I am going through several papers about aerosol physics and in this field you often have determine the parameters of a distribution to match your experimental data (one typically uses a Gaussian mixture). However, in many cases people plot a normalized empirical distribution function and then perform some least-square fitting rather than using likelihood functions. As an undergrad, I was told that the former approach is correct only if you have a model for the dynamics (e.g. Ohm law and you perform a least-square fitting), but not if you deal with a distribution and you pick random draws from it (in that case, one should maximize the probability of drawing the data which were actually observed and this leads to likelihood functions). The two approaches do not seem equivalent to me, but I cannot believe that this distinction is ignored in practice... Many thanks Lorenzo __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Is there a way to draw 3d plot?
Michael wrote: Hi all, I have a 2D matrix, which has 100 rows, and 100 columns, I have a 2D matrix, with 100 rows and 100 columns, I want to display it using 3D plot, much like plot3d and mesh/surf functions in matlab. Specifically, in matlab, I just need to do the following: [X, Y]=meshgrid([0:0.01:0.99, 0:0.01:0.99]); % Z is my 2D matrix, surf(X, Y, Z); Note that X and Y are created so that I can associate physical meaning onto the x and y axis of the 3D plot. For example, my 100 rows represent 0, 0.01, 0.02, ... 0.99 here. In Matlab I can also drag in the graphic window and see from different visual angle and perspective of the 3D plot... Are there similar functions in R that (1) show 3D plot; (2) let me manipulate view angles easily? (1) See ?persp (1) *and* (2): See package rgl. Uwe Ligges Thanks a lot! [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Is there a way to draw 3d plot?
On 5/26/06, Uwe Ligges [EMAIL PROTECTED] wrote: Michael wrote: Hi all, I have a 2D matrix, which has 100 rows, and 100 columns, I have a 2D matrix, with 100 rows and 100 columns, I want to display it using 3D plot, much like plot3d and mesh/surf functions in matlab. Specifically, in matlab, I just need to do the following: [X, Y]=meshgrid([0:0.01:0.99, 0:0.01:0.99]); % Z is my 2D matrix, surf(X, Y, Z); Note that X and Y are created so that I can associate physical meaning onto the x and y axis of the 3D plot. For example, my 100 rows represent 0, 0.01, 0.02, ... 0.99 here. In Matlab I can also drag in the graphic window and see from different visual angle and perspective of the 3D plot... Are there similar functions in R that (1) show 3D plot; (2) let me manipulate view angles easily? (1) See ?persp (1) *and* (2): See package rgl. Uwe Ligges Thanks a lot, But a glance at rgl seems requireing shape, etc... and very complicated... Any easier approaches? persp does not allow me to use mouse to rotate [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] missed ylim from plot.default
2006/5/26, Peter Ehlers [EMAIL PROTECTED]: Antonio, Fabio Di Narzo wrote: Hi all. On my R-2.3.0, calling: plot(1:5, rep(1,5), xlim=c(-0.5,7), ylim=c(-0.5,2.8), asp=1) gives to me a plot (tried pdf and X11) whose y axis goes from about -2 to about 4.5. What I've missed? How can I show some pre-specified x and y ranges from a plot (keeping fixed the aspect ratio)? Tnx all, Antonio, Fabio Di Narzo. [[alternative HTML version deleted]] Do you want a wide, but not very tall plot? If so, you could specify the width/height in your X11() call. Peter Ehlers Tnx for your indication, now I see: x and y plotting ranges are influenced by the actual window size (at least when fixing aspect ratio). I think this is obvious... What I missed/forgot above, was that default plot width and height are fixed indipendently from the 'x-ylim' and 'asp' arguments. However now I have a problem: I can write a wrapper to set device width and height depending on the (xlim,ylim,asp) triple. But there's often extra-space in the device due to plot title, axis labels, etc. So, I fear I can't simply set width=height*asp. Some suggestion? (Anyway, by now I only have to produce few plots, so today I can go try by try :-) ). Bests, Antonio. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Indexing vector with repeated values
Hi all, I have a vector which contains many repeated values. Z - c(1,2,3,4,5,1,2,3,4,5,1,2,3,5,5,2,3,4,5) I need to know how many times each number in this vecetor is repeated. I can use a command like this: table (cut(Z,seq(1,5,1))) (1,2] (2,3] (3,4] (4,5] 4 4 3 5 But how can I find a break points for vector with random values and random number sequence intervals? Why result for interval (1,2) is 4, if there is only 3 values of 1? Can I get from vector Z a smal vector Zs 1,2,3,4,5 ? Thanks a lot! Andris Jankevics __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Labels in cluster plots
I'm trying to follow an example from ?clara, which plots two clusters: x - rbind(cbind(rnorm(10,0,8), rnorm(10,0,8)), cbind(rnorm(20,50,8), rnorm(20,50,8))) clarax - clara(x, 2) clusplot(clarax) Suppose I'd like to label these 30 observations on the plot according to a preconceived classification, e.g. c(rep(A,15),rep(B,7),rep(C,8)) Is there a way to plot observations with these pre-specified labels? That is, I'd like the first 15 of 'x' to be labeled as 'A', the next seven as 'B', the next eight as 'C'. (P.S. For principal component plots it was helpfully suggested to use the following as an argument to the plot: pch=c('A','B','C')[factor(c(rep(A,15),rep(B,7),rep(C,8)))]. I can't figure out if there is a similar way here, and what is a general way to put repeated row labels on various graphs...) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Is there a way to draw 3d plot?
2006/5/26, Michael [EMAIL PROTECTED]: On 5/26/06, Uwe Ligges [EMAIL PROTECTED] wrote: Michael wrote: Hi all, I have a 2D matrix, which has 100 rows, and 100 columns, I have a 2D matrix, with 100 rows and 100 columns, I want to display it using 3D plot, much like plot3d and mesh/surf functions in matlab. Specifically, in matlab, I just need to do the following: [X, Y]=meshgrid([0:0.01:0.99, 0:0.01:0.99]); % Z is my 2D matrix, surf(X, Y, Z); Note that X and Y are created so that I can associate physical meaning onto the x and y axis of the 3D plot. For example, my 100 rows represent 0, 0.01, 0.02, ... 0.99 here. In Matlab I can also drag in the graphic window and see from different visual angle and perspective of the 3D plot... Are there similar functions in R that (1) show 3D plot; (2) let me manipulate view angles easily? (1) See ?persp (1) *and* (2): See package rgl. Uwe Ligges Thanks a lot, But a glance at rgl seems requireing shape, etc... and very complicated... Any easier approaches? persp does not allow me to use mouse to rotate [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Have you seen the example in 'rgl' man page? It seems sufficient something like: rgl.surface(x,y,z) In that example, is also showed how to colorize the surface. Antonio. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Indexing vector with repeated values
Andris Jankevics wrote: Hi all, I have a vector which contains many repeated values. Z - c(1,2,3,4,5,1,2,3,4,5,1,2,3,5,5,2,3,4,5) I need to know how many times each number in this vecetor is repeated. I can use a command like this: table (cut(Z,seq(1,5,1))) Why not simply table(Z) ? (1,2] (2,3] (3,4] (4,5] 4 4 3 5 But how can I find a break points for vector with random values and random number sequence intervals? Why result for interval (1,2) is 4, if there is only 3 values of 1? There is no 1 counted, but 4 2s, you have not specified an interval with 1 in it. You were looking for table(cut(Z,seq(0,5,1))) but don't need cut() here at all. Can I get from vector Z a smal vector Zs 1,2,3,4,5 ? seq(min(Z), max(Z)) Uwe Ligges Thanks a lot! Andris Jankevics __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] changing font size in plot(effect())
Dear Peter, Won't that wipe out the other components of axis.text, etc.? Regards, John On Thu, 25 May 2006 06:22:22 -0600 P Ehlers [EMAIL PROTECTED] wrote: John Fox wrote: Dear Emilie, This is, I guess, the effect() function in the effects package. If so, note that the plot.effect() method uses trellis graphics (via the lattice package), not standard R graphics, so you have to control aspects of the plot in a manner consistent with trellis graphics. Unfortunately, you can't just specify the scales argument when you call plot(), since plot.effect() already includes a scales argument when it calls xyplot(). You can, however, set trellis parameters globally, e.g., via something like axis.text - trellis.par.get(axis.text) par.ylab.text - trellis.par.get(par.ylab.text) par.xlab.text - trellis.par.get(par.xlab.text) axis.text$cex - 1.5 par.ylab.text$cex - 1.5 par.xlab.text$cex - 1.5 trellis.par.set(axis.text, axis.text) trellis.par.set(par.ylab.text, par.ylab.text) trellis.par.set(par.xlab.text, par.xlab.text) Perhaps this can be done more efficiently -- I'm far from a trellis whiz -- but the above should work. More generally, in designing the plot method for effect objects, it wasn't my goal to produce publication-quality plots; for that, I suggest that you build a custom graph using the information included in the effect object. I hope this helps, John Possibly a bit more efficient: trellis.par.set(list(axis.text = list(cex = 2), par.ylab.text = list(cex = 1.5), par.xlab.text = list(cex = 0.5))) Also good to know, Emilie: trellis.par.get() to see all the things that can be (re)set. Most are self-explanatory. - Peter Ehlers On Wed, 24 May 2006 13:24:03 -0400 Emilie Berthiaume [EMAIL PROTECTED] wrote: I can't seem to be able to change the font size in an effect display. I've tried the following: par(cex.lab=4) plot(effect (alti,reg8), ylab=detection probability) and plot(effect (alti,reg8), ylab=detection probability, cex=4) but nothing changes. Can anyone help me? thanks. John Fox Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox/ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] factors and ops
Try an ordered factor: subset(HR, ordered(Patient) 214) or HR$Patient - ordered(HR$Patient) subset(HR, Patient 214) You might also check that it orders them in the way you want. ordered(HR$Patient) in the first case or just HR$Patient in the second case will display the data followed by the ordering. On 5/26/06, Erin Hodgess [EMAIL PROTECTED] wrote: Dear R People I have a variable which is an ID number that is a factor. I would like to look for ID numberbs that are greeater than a particular value. However, factors are a big ugly for these operations. I was messing with the HR data set from the SASmixed package. HR$Patient is a factor like that. Now if I used: subset(HR,as.integer(as.character(Patient)) 214) that will work. If seems to me that there may be a better way. Is that true? R Version 2.3.0 Windows Thanks in advance! sincerely, Erin Hodgess Associate Professor Department of Computer and Mathematical Sciences University of Houston - Downtown mailto: [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Vector elements and ratios
Dear useRs, I have two different length vectors: one column (1...m) and one row vector (1...n): 20 40 20 60 5 4 2 Now I have to calculate ratios between column vector elements and each row vector elements: 4 5 10 8 10 20 4 5 20 15 12 30 Thank's in advance for any suggestions, Andrej __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Vector elements and ratios
Andrej Kastrin wrote: Dear useRs, I have two different length vectors: one column (1...m) and one row vector (1...n): 20 40 20 60 5 4 2 Now I have to calculate ratios between column vector elements and each row vector elements: 4 5 10 8 10 20 4 5 20 15 12 30 A - c(20,40,20,60) B - c(5,4,2) A %*% t(1/B) [,1] [,2] [,3] [1,]45 10 [2,]8 10 20 [3,]45 10 [4,] 12 15 30 Thank's in advance for any suggestions, Andrej __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Vector elements and ratios
?outer v1 - c(20,40,20,60) v2 - c(5,4,2) outer(v1,v2,'/') [,1] [,2] [,3] [1,]45 10 [2,]8 10 20 [3,]45 10 [4,] 12 15 30 On 5/26/06, Andrej Kastrin [EMAIL PROTECTED] wrote: Dear useRs, I have two different length vectors: one column (1...m) and one row vector (1...n): 20 40 20 60 5 4 2 Now I have to calculate ratios between column vector elements and each row vector elements: 4 5 10 8 10 20 4 5 20 15 12 30 Thank's in advance for any suggestions, Andrej __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Jim Holtman Cincinnati, OH +1 513 646 9390 (Cell) +1 513 247 0281 (Home) What is the problem you are trying to solve? [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Vector elements and ratios
outer(c(20,40,20,60), c(5,4,2), /) [,1] [,2] [,3] [1,]45 10 [2,]8 10 20 [3,]45 10 [4,] 12 15 30 cheers D On 26/05/06, Andrej Kastrin [EMAIL PROTECTED] wrote: Dear useRs, I have two different length vectors: one column (1...m) and one row vector (1...n): 20 40 20 60 5 4 2 Now I have to calculate ratios between column vector elements and each row vector elements: 4 5 10 8 10 20 4 5 20 15 12 30 Thank's in advance for any suggestions, Andrej __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Simulate Strauss process in 3D
Apple Ho writes: I am searching some information about simulating Strauss process in 3D If you need this immediately, I suggest you follow Brian Ripley's advice. This probably involves downloading a source tar file of the R package and finding the Fortran code in library/spatial/src/. The Fortran should be edited to generate 3D patterns (it's not hard: everywhere you see 'x' and 'y' just include 'z'). Then also edit the R interface to the Fortran. If you can wait about 3 months, we will be releasing an upgrade of the R package `spatstat' that will deal with 3D point patterns, including simulation of the Strauss process in 3D. regards Adrian Baddeley __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Is there a way to draw 3d plot?
On 5/26/2006 3:39 AM, Antonio, Fabio Di Narzo wrote: 2006/5/26, Michael [EMAIL PROTECTED]: On 5/26/06, Uwe Ligges [EMAIL PROTECTED] wrote: Michael wrote: Hi all, I have a 2D matrix, which has 100 rows, and 100 columns, I have a 2D matrix, with 100 rows and 100 columns, I want to display it using 3D plot, much like plot3d and mesh/surf functions in matlab. Specifically, in matlab, I just need to do the following: [X, Y]=meshgrid([0:0.01:0.99, 0:0.01:0.99]); % Z is my 2D matrix, surf(X, Y, Z); Note that X and Y are created so that I can associate physical meaning onto the x and y axis of the 3D plot. For example, my 100 rows represent 0, 0.01, 0.02, ... 0.99 here. In Matlab I can also drag in the graphic window and see from different visual angle and perspective of the 3D plot... Are there similar functions in R that (1) show 3D plot; (2) let me manipulate view angles easily? (1) See ?persp (1) *and* (2): See package rgl. Uwe Ligges Thanks a lot, But a glance at rgl seems requireing shape, etc... and very complicated... Any easier approaches? persp does not allow me to use mouse to rotate [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Have you seen the example in 'rgl' man page? It seems sufficient something like: rgl.surface(x,y,z) In that example, is also showed how to colorize the surface. And an in-development update to rgl adds a number of new routines with interfaces similar to ones from the graphics package. For example, persp3d is probably what Michael wants. If you like playing with experimental versions and don't mind having defaults changing with new releases, etc., you can download a recent build from http://www.stats.uwo.ca/faculty/murdoch/software Bug reports and suggestions about the new stuff would be welcome. Duncan __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Thanks for help
Hello David and Spencer, thank you for your replies! The systemfit-package was exactly what I was looking for. The only little drawback is, that the nonlinear systemfit does not support predifined weights yet (in oder to ensure homegenity of variances when predictios are hetescedastic) Thanks also to Katharina, who mailed me offline and drew my attention to the nls and the nlm function. This answerd my second question. Thomas __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] lme, best model without convergence
Dear R-help list readers, I am fitting mixed models with the lme function of the nlme package. If I get convergence depends on how the method (ML/REM) and which (and how much) parameters will depend randomly on the cluster-variable. How get the bist fit without convergence? I set the parameters msVerbose and returnObject to TRUE: lmeControl(maxIter=5, msMaxIter=200, tolerance=1e-4, niter=50, msTol=1e-5, nlmStepMax=500, ,msVerbose=TRUE ,returnObject=TRUE ) However, the lme-functions does not produce verbose output, nor does it return the best fit if lme is not converging. It returns only an error: Error in lme.formula(y ~ lndbh + I(lndbh^2) + lnh + I(lnh^2), random = ~lndbh + : iteration limit reached without convergence (9) Best regards Thomas __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Computing a reliability index of a statistic with missing data
Thanks Spencer, that is interesting but I must say I'm a bit lost with the terminology. I'll try to catch up but I'm not sure I need a complicated model (MC sounds complicated to me but it may not be...). I plan to use this reliability index just as an indication and I need to compute it in batch for several different charts so I try to keep the statistic as simple as possible but yet efficient. Aziz -Original Message- From: Spencer Graves [mailto:[EMAIL PROTECTED] Sent: May 25, 2006 8:12 PM To: Chaouch, Aziz Cc: R-help@stat.math.ethz.ch Subject: Re: [R] Computing a reliability index of a statistic with missing data Have you considered some kind of binary time series model? 'RSiteSearch(binary time series)' produced 150 hits. One of the first 20 mentioned continuous-time hidden Markov chains (http://finzi.psych.upenn.edu/R/library/repeated/html/chidden.html). I don't know if this will help you or not, but it might be worth examining. hope this helps. Spencer Graves Chaouch, Aziz wrote: Hi All, I'd like to compute a kind of reliability index (RI) that would in a sense stand as a measure of reliability of a statistic (histogram etc) computed on a time serie with missing values. The final goal is that: RI=1 for a perfect reliability RI=0 for a total unreliability (no data at all as an extreme case...) The percentage of missing data is one indication: the more missing data, the less confidence we can have in the statistic. But the distribution of missing data throughout the data serie is important as well: independently of the number of missing data, if available data are regularily spaced in time the RI should be higher than if available data are irregulary spaced. As a measure of sampling regularity, I thought about computing the time to next record and then take its variance over the time interval on which the statistic is computed. The variance of the time to next record would be a measure of sampling regularity so that the final RI could be of the form: RI=1 when n=0 RI~1/n*var(T) with n=% of missing data T=time to next record (in hours) However I need to normalize var(T) to use it to compute the RI. Does someone have an idea on how to do this (or another proposal to compute the RI)? Thanks, Aziz [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Try to debbug R script
Hi, I have a bug and I don't find it. Can u help me please. In attachement : - My script file : CPU_study.txt - My data sample file : a.txt When I use it, the error message is : source(H:/R/_workspace/CPU_study_1.R) Erreur dans [-(`*tmp*`, i, k, value = numeric(0)) : rien à remplacer (= error in ... nothing to replace) I try a lot of thing but it doesn't work. The probleme is line 39 : val - (val + t[(((i*step)-step)+j),k]) Thanks Guillaume LAURENS - trainee student Simulation Industry - EYYSIDV Airbus France This e-mail is intended only for the above addressee. It may contain privileged information. If you are not the addressee you must not copy, distribute, disclose or use any of the information in it. If you have received it in error please delete it and immediately notify the sender. Security Notice: all e-mail, sent to or from this address, may be accessed by someone other than the recipient, for system management and security reasons. This access is controlled under Regulation of Investigatory Powers Act 2000, Lawful Business Practises. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] query: sample size calculation
Dear all, I am doing something wrong. I am trying to apply a formula for sample size calculation as in the book Design and Analysis of Clinical Trials, from Chow et all. There a suggested sample size strategy uses the formula 0.30/0.45=F(0.80,2n,2n)/F(0.025,2n,2n) which gives n=96, as in the book. (here F(alfa,k,n) is the upper (alfa)th quantile of an F distribution with k,n degrees of freedom) I have been trying to get n=96 using the following code in R. val-rep(NA,200) for (n in 1:200) { val[n]- qf(0.8, 2*n, 2*n,lower.tail = TRUE, log.p = FALSE)/qf(0.025, 2*n, 2*n,lower.tail = TRUE, log.p = FALSE); } but val doesnot get any close to 0.30/0.45. Thank you, Manuel __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Try to debbug R script
On 5/26/2006 8:45 AM, LAURENS, Guillaume wrote: Hi, I have a bug and I don't find it. Can u help me please. In attachement : - My script file : CPU_study.txt - My data sample file : a.txt When I use it, the error message is : source(H:/R/_workspace/CPU_study_1.R) Erreur dans [-(`*tmp*`, i, k, value = numeric(0)) : rien à remplacer (= error in ... nothing to replace) I try a lot of thing but it doesn't work. The probleme is line 39 : val - (val + t[(((i*step)-step)+j),k]) The error says that your indexing returned a zero length vector. I'd suggest printing out the values of (((i*step)-step)+j) and k just before this line (or using the debugger to examine them). They probably aren't in the range you were expecting. Duncan Murdoch __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] changing font size in plot(effect())
No, I don't think so, John. I think what I proposed is Deepayan's newer way of setting parameters. I used to do it the other way, too (S-PLUS may still need it), but I find this new way very handy, since it avoids the get/save/alter/set process. Try this: trellis.device() trellis.par.set(list(axis.text = list(col = blue))) # reset colour trellis.par.get(axis.text) # inspect trellis.par.set(list(axis.text = list(cex = 2))) # reset cex trellis.par.get(axis.text) # inspect; colour is still 'blue' Best wishes, Peter John Fox wrote: Dear Peter, Won't that wipe out the other components of axis.text, etc.? Regards, John On Thu, 25 May 2006 06:22:22 -0600 P Ehlers [EMAIL PROTECTED] wrote: John Fox wrote: Dear Emilie, This is, I guess, the effect() function in the effects package. If so, note that the plot.effect() method uses trellis graphics (via the lattice package), not standard R graphics, so you have to control aspects of the plot in a manner consistent with trellis graphics. Unfortunately, you can't just specify the scales argument when you call plot(), since plot.effect() already includes a scales argument when it calls xyplot(). You can, however, set trellis parameters globally, e.g., via something like axis.text - trellis.par.get(axis.text) par.ylab.text - trellis.par.get(par.ylab.text) par.xlab.text - trellis.par.get(par.xlab.text) axis.text$cex - 1.5 par.ylab.text$cex - 1.5 par.xlab.text$cex - 1.5 trellis.par.set(axis.text, axis.text) trellis.par.set(par.ylab.text, par.ylab.text) trellis.par.set(par.xlab.text, par.xlab.text) Perhaps this can be done more efficiently -- I'm far from a trellis whiz -- but the above should work. More generally, in designing the plot method for effect objects, it wasn't my goal to produce publication-quality plots; for that, I suggest that you build a custom graph using the information included in the effect object. I hope this helps, John Possibly a bit more efficient: trellis.par.set(list(axis.text = list(cex = 2), par.ylab.text = list(cex = 1.5), par.xlab.text = list(cex = 0.5))) Also good to know, Emilie: trellis.par.get() to see all the things that can be (re)set. Most are self-explanatory. - Peter Ehlers On Wed, 24 May 2006 13:24:03 -0400 Emilie Berthiaume [EMAIL PROTECTED] wrote: I can't seem to be able to change the font size in an effect display. I've tried the following: par(cex.lab=4) plot(effect (alti,reg8), ylab=detection probability) and plot(effect (alti,reg8), ylab=detection probability, cex=4) but nothing changes. Can anyone help me? thanks. John Fox Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox/ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Try to debbug R script
LAURENS, Guillaume a écrit : Hi, I have a bug and I don't find it. Can u help me please. In attachement : - My script file : CPU_study.txt - My data sample file : a.txt When I use it, the error message is : source(H:/R/_workspace/CPU_study_1.R) Erreur dans [-(`*tmp*`, i, k, value = numeric(0)) : rien à remplacer (= error in ... nothing to replace) I try a lot of thing but it doesn't work. The probleme is line 39 : val - (val + t[(((i*step)-step)+j),k]) Have you try the debug library ? With the function mtrace() you can follow the work step by step and in the same time with the function bp() you can add break points... Etienne Thanks Guillaume LAURENS - trainee student Simulation Industry - EYYSIDV Airbus France This e-mail is intended only for the above addressee. It may contain privileged information. If you are not the addressee you must not copy, distribute, disclose or use any of the information in it. If you have received it in error please delete it and immediately notify the sender. Security Notice: all e-mail, sent to or from this address, may be accessed by someone other than the recipient, for system management and security reasons. This access is controlled under Regulation of Investigatory Powers Act 2000, Lawful Business Practises. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] query: sample size calculation
[EMAIL PROTECTED] wrote: Dear all, I am doing something wrong. I am trying to apply a formula for sample size calculation as in the book Design and Analysis of Clinical Trials, from Chow et all. There a suggested sample size strategy uses the formula 0.30/0.45=F(0.80,2n,2n)/F(0.025,2n,2n) which gives n=96, as in the book. (here F(alfa,k,n) is the upper (alfa)th quantile of an F distribution with k,n degrees of freedom) I have been trying to get n=96 using the following code in R. val-rep(NA,200) for (n in 1:200) { val[n]- qf(0.8, 2*n, 2*n,lower.tail = TRUE, log.p = FALSE)/qf(0.025, 2*n, 2*n,lower.tail = TRUE, log.p = FALSE); } Don't you want 'lower.tail = FALSE'? Peter Ehlers but val doesnot get any close to 0.30/0.45. Thank you, Manuel __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] changing font size in plot(effect())
Dear Peter, OK -- I see that this works as advertized. Thanks, John On Fri, 26 May 2006 07:12:52 -0600 P Ehlers [EMAIL PROTECTED] wrote: No, I don't think so, John. I think what I proposed is Deepayan's newer way of setting parameters. I used to do it the other way, too (S-PLUS may still need it), but I find this new way very handy, since it avoids the get/save/alter/set process. Try this: trellis.device() trellis.par.set(list(axis.text = list(col = blue))) # reset colour trellis.par.get(axis.text) # inspect trellis.par.set(list(axis.text = list(cex = 2))) # reset cex trellis.par.get(axis.text) # inspect; colour is still 'blue' Best wishes, Peter John Fox wrote: Dear Peter, Won't that wipe out the other components of axis.text, etc.? Regards, John On Thu, 25 May 2006 06:22:22 -0600 P Ehlers [EMAIL PROTECTED] wrote: John Fox wrote: Dear Emilie, This is, I guess, the effect() function in the effects package. If so, note that the plot.effect() method uses trellis graphics (via the lattice package), not standard R graphics, so you have to control aspects of the plot in a manner consistent with trellis graphics. Unfortunately, you can't just specify the scales argument when you call plot(), since plot.effect() already includes a scales argument when it calls xyplot(). You can, however, set trellis parameters globally, e.g., via something like axis.text - trellis.par.get(axis.text) par.ylab.text - trellis.par.get(par.ylab.text) par.xlab.text - trellis.par.get(par.xlab.text) axis.text$cex - 1.5 par.ylab.text$cex - 1.5 par.xlab.text$cex - 1.5 trellis.par.set(axis.text, axis.text) trellis.par.set(par.ylab.text, par.ylab.text) trellis.par.set(par.xlab.text, par.xlab.text) Perhaps this can be done more efficiently -- I'm far from a trellis whiz -- but the above should work. More generally, in designing the plot method for effect objects, it wasn't my goal to produce publication-quality plots; for that, I suggest that you build a custom graph using the information included in the effect object. I hope this helps, John Possibly a bit more efficient: trellis.par.set(list(axis.text = list(cex = 2), par.ylab.text = list(cex = 1.5), par.xlab.text = list(cex = 0.5))) Also good to know, Emilie: trellis.par.get() to see all the things that can be (re)set. Most are self-explanatory. - Peter Ehlers On Wed, 24 May 2006 13:24:03 -0400 Emilie Berthiaume [EMAIL PROTECTED] wrote: I can't seem to be able to change the font size in an effect display. I've tried the following: par(cex.lab=4) plot(effect (alti,reg8), ylab=detection probability) and plot(effect (alti,reg8), ylab=detection probability, cex=4) but nothing changes. Can anyone help me? thanks. John Fox Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox/ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html John Fox Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox/ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] MART(tm) vs. gbm
Hello, I tried two different implementations of the stochastic gradient boosting (Friedman 2002) : the MART(tm) with R tool (http://www-stat.stanford.edu/~jhf/R-MART.html) and the gbm R package. To me, both seemed fairly comparable, except maybe regarding the different loss criterion proposed and the fact that the gbm tool is sightly more convenient to use. However, it seemed that the MART with R tool systematically outperforms the gbm tool in terms of goodness of fit (whatever the way of choosing the best iteration for the gbm package). I tried to find out if there were specific options that could have explained it but nothing came out. See below for an example of how I compare both implementations. Did any one had the same experience, and can anyone give me hints about such performance differences or tell me if I am missing something obvious? Thank you in advance,Manuel Here are the arguments and options I used for comparison purposes, working on a 1600 records * 15 variables dataset : # the MART with R tool lx - mart( as.matrix(x), y, c(1,1,1,1,1,1,1,1,1,1,1,1,1,2,2) niter=1000, tree.size=6, learn.rate=0.01, loss.cri=2 # gaussian ) # for gbm gbm1 - gbm(y ~ v1 + v2 + v3 + v4 + v5+ v6+ v7+ v8 + v9 + v10 + v11 + v12 + v13 + v14 + v15, data=data, var.monotone=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0), distribution=gaussian, n.trees=1000, shrinkage=0.01, interaction.depth=6, bag.fraction = 0.5, train.fraction = 0.5, n.minobsinnode = 10, cv.folds = 1, keep.data=TRUE) # I then do predictions on the same dataset, and further perform goodness of fit comparisons #... __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] I cannot load the package Rcmdr
I can load the package in other R file, but I can not do it in the file. How to solve the question? local({pkg - select.list(sort(.packages(all.available = TRUE))) + if(nchar(pkg)) library(pkg, character.only=TRUE)}) Loading required package: tcltk Loading required package: car Error in parse(file, n, text, prompt) : syntax error in * Error: .onAttach failed in 'attachNamespace' Error: package/namespace load failed for 'Rcmdr' [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] MART(tm) vs. gbm
Hello, I have been using two different implementations of the stochastic gradient boosting (Friedman 2002) : MART(tm) with R and the gbm package. Both are fairly comparable except that the MART with R systematically strongly (depending on the dataset though) outperforms the gbm tool in terms of goodness of fit. For instance, a # gbm package gbm1 - gbm(Y~X2+X3+X4+X5+X6, data=data, var.monotone=c(0,0,0,0,0), # 0: no monotone restrictions distribution=gaussian, # bernoulli, adaboost, gaussian, # poisson, and coxph available n.trees=3000,# number of trees shrinkage=0.005, # shrinkage or learning rate, # 0.001 to 0.1 usually work interaction.depth=6, # 1: additive model, 2: two-way interactions, etc. bag.fraction = 0.5, # subsampling fraction, 0.5 is probably best train.fraction = 0.5,# fraction of data for training, # first train.fraction*N used for training n.minobsinnode = 10, # minimum total weight needed in each node cv.folds = 5,# do 5-fold cross-validation keep.data=TRUE, # keep a copy of the dataset with the object verbose=TRUE)# print out progress # MART with R X - as.matrix(cbind(data$X2,as.numeric(data$X3), as.numeric(data$X4),as.numeric(data$X5),data$X6)) Y - data$Y mart(X, Y, c(1,2,2,2,1) , niter=3000, tree.size=6,learn.rate=0.005, loss.cri=2 #gaussian too ) leads to very different goodnesses of fit (I can provide the dataset if needed). Did anyone already encountered this, is there an explanation, am I missing something obvious in the argument settings? Thank you in advance, Manuel __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Function as.Date leading to error implying that strptime requires 3 arguments
Thank you to Prof Ripley, Gabor Grothendieck, Dirk Eddelbuettel and others who emailed with suggestions and comments. Dr Ripley's suggestion that I reinstall appears to have fixed the problem, though I confess I have upgraded only to 2.3.0 rather than the beta 2.3.1. I use Windows XP Pro, Service Pack 2.0 -- my apologies for not mentioning this. While this instance was not serious, in itself, is it a symptom of something more troubling? I worry that some feature of my usage has led to this problem. However, if it were something I did regularly, then it should recur. At that time, I will certainly have to try to correct whatever it happens to be. In the meantime, however, if anyone has any suggestions about what I might have done to create this problem I'd be pleased to hear them. I do (or rather, did) have multiple versions of R installed, and I do make extensive use of several libraries, including the Bioconductor set. My suspicion is that this would not occur if I were to remove old versions when I upgrade. Cheers, Rob -Original Message- From: Prof Brian Ripley [mailto:[EMAIL PROTECTED] Sent: Friday, May 19, 2006 11:39 PM To: Rob Balshaw Cc: r-help@stat.math.ethz.ch Subject: Re: [R] Function as.Date leading to error implying that strptime requires 3 arguments Your system (unstated) is corrupted. R runs its examples as part of it test suite, so this is not an error in R 2.2.1 (which is not current: see the posting guide which asked you to update *before* posting). One possibility is that you have strptime from a much earlier version of R in use somehow. I suggest you remove all versions of R from your system and reinstall one latest version (2.3.1 beta). On Fri, 19 May 2006, Rob Balshaw wrote: I'm using R V 2.2.1. When I try an example from the as.Date help page, I get an error. x - c(1jan1960, 2jan1960, 31mar1960, 30jul1960) z - as.Date(x, %d%b%Y) Error in strptime(x, format) : 2 arguments passed to 'strptime' which requires 3 Any suggestions would be appreciated. Thanks, Rob __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Is there a way to draw 3d plot?
There is also the rotate.persp function in the TeachingDemos package. It creates a set of slider bars that you can move with the mouse to change the different options to persp including the angles. It is not quite as convenient as clicking in the plot and dragging like rgl allows, but it does have the advantage that when you have the plot you like, you can see what the current settings are and recreate the same plot using persp and those settings rather than having to remember how you rotated the graph interactively (you can query the current rotation in rgl as well and set it by command in future plots). -- Gregory (Greg) L. Snow Ph.D. Statistical Data Center Intermountain Healthcare [EMAIL PROTECTED] (801) 408-8111 -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Michael Sent: Friday, May 26, 2006 1:25 AM To: Uwe Ligges Cc: R-help@stat.math.ethz.ch Subject: Re: [R] Is there a way to draw 3d plot? On 5/26/06, Uwe Ligges [EMAIL PROTECTED] wrote: Michael wrote: Hi all, I have a 2D matrix, which has 100 rows, and 100 columns, I have a 2D matrix, with 100 rows and 100 columns, I want to display it using 3D plot, much like plot3d and mesh/surf functions in matlab. Specifically, in matlab, I just need to do the following: [X, Y]=meshgrid([0:0.01:0.99, 0:0.01:0.99]); % Z is my 2D matrix, surf(X, Y, Z); Note that X and Y are created so that I can associate physical meaning onto the x and y axis of the 3D plot. For example, my 100 rows represent 0, 0.01, 0.02, ... 0.99 here. In Matlab I can also drag in the graphic window and see from different visual angle and perspective of the 3D plot... Are there similar functions in R that (1) show 3D plot; (2) let me manipulate view angles easily? (1) See ?persp (1) *and* (2): See package rgl. Uwe Ligges Thanks a lot, But a glance at rgl seems requireing shape, etc... and very complicated... Any easier approaches? persp does not allow me to use mouse to rotate [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to multiply a constant to a matrix?
I still can't see why this is a problem. If a 1x1 matrix should be treated as a scalar, then it can just be wrapped in drop(), and the arithmetic will be computed correctly by R. Are there any cases where this cannot be done? More specifically, are there any matrix algebra expressions where, depending on the particular dimensions of the variables used, drop() must be used in some cases, and not in other cases? A related but different behavior is the default dropping dimensions with extent equal to one by indexing operations. This can be problematic because if one is not careful, incorrect results can be obtained for particular values used in the expression. For example, consider the following, in which we are trying to compute the cross product of some columns of x with some rows of y. If x has n rows and y has n columns, then the result should always be an nxn matrix. However, if we are not careful with using drop=F in the indexing expressions, we can inadvertently end up with a 1x1 inner product matrix result for the case where we just use one column of x and one row of y. The solution to this is to always use drop=F in indexing in situations where this can occur. x - matrix(1:9, ncol=3) y - matrix(-(1:9), ncol=3) i - 1:2 x[,i] %*% y[i,] [,1] [,2] [,3] [1,] -9 -24 -39 [2,] -12 -33 -54 [3,] -15 -42 -69 i - 1:3 x[,i] %*% y[i,] [,1] [,2] [,3] [1,] -30 -66 -102 [2,] -36 -81 -126 [3,] -42 -96 -150 # i has just one element -- the expression without drop=F # no longer computes an outer product i - 2 x[,i] %*% y[i,] [,1] [1,] -81 x[,i,drop=F] %*% y[i,,drop=F] [,1] [,2] [,3] [1,] -8 -20 -32 [2,] -10 -25 -40 [3,] -12 -30 -48 Cannot all cases in the situations you mention be handled in an analogous manner, by always wrapping appropriate quadratic expressions in drop(), or are there some cases where the result of the quadratic expression must be treated as a matrix, and other cases where the result of the quadratic expression must be treated as a scalar? -- Tony Plate Michael wrote: imagine when you have complicated matrix algebra computation using R, you cannot prevent some middle-terms become quadratic and absorbs into one scalar, right? if R cannot intelligently determine this, and you have to manually add drop everywhere, do you think it is reasonable? On 5/23/06, Patrick Burns [EMAIL PROTECTED] wrote: I think drop(B/D) * solve(A) would be a more transparent approach. It isn't that R can not do what you want, it is that it is saving you from shooting yourself in the foot in your attempt. What you are doing is not really a matrix computation. Patrick Burns [EMAIL PROTECTED] +44 (0)20 8525 0696 http://www.burns-stat.com (home of S Poetry and A Guide for the Unwilling S User) Michael wrote: This is very strange: I want compute the following in R: g = B/D * solve(A) where B and D are quadratics so they are just a scalar number, e.g. B=t(a) %*% F %*% a; I want to multiply B/D to A^(-1), but R just does not allow me to do that and it keeps complaining that nonconformable array, etc. I tried the following two tricks and they worked: as.numeric(B/D) * solve(A) diag(as.numeric(B/D), 5, 5) %*% solve (A) But if R cannot intelligently do scalar and matrix multiplication, it is really problemetic. It basically cannot be used to do computations, since in complicated matrix algebras, you have to distinguish where is scalar, and scalars obtained from quadratics cannot be directly used to multiply another matrix, etc. It is going to a huge mess... Any thoughts? [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] upgrade or recompile
Dear list, My current R 2.2.0 in my linux desktop was installed by compiling sources, now I want to upgrade it to 2.3.0, is it possible to upgrade it without re-compiling sources? Thanks, mike [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] multiple comparisons of time series data
I am interested in a statistical comparison of multiple (5) time series' generated from modeling software (Hydrologic Simulation Program Fortran). The model output simulates daily bacteria concentration in a stream. The multiple time series' are a result of varying our representation of the stream within the model. Our main question is: Do the different methods used to represent a stream produce different results at a statistically significant level? We want to compare each otput time series to determine if there is a difference before looking into the cause within the model. In a previous study, the Kolmogorov-Smirnov k-sample test was used to compare multiple time series'. I am unsure about the strength of the Kolmogorov-Smirnov test and I have set out to determine if there are any other tests to compare multiple time series'. I know htat R has the ks.test but I am unsure how this test handles multiple comparisons. Is there something similar to a pairwise.t.test with a bonferroni corection, only with time series data? Does R currently (v 2.3.0) have a comparison test that takes into account the strong serial correlation of time series data? Kyle Hall Graduate Research Assistant Biological Systems Engineering Virginia Tech __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution
Hi, I'm still strugling with this copula model but this seems to be the way to go. I'm now trying to model the marginal distributions and and for wind direction I use a mixture of 2 von mises. I'd like to estimate all the parameters (m1,m1,kappa1,kappa2,p) by maximizing the likelihood but I don't know how to define the likelihood (or log-likelihood) of a mixture of 2 Von Mises to use it with the function fitDistr in the MASS package. Can you help me define this likelihood function and use it through the fitDistr function? Thanks, Aziz -Original Message- From: Dimitrios Rizopoulos [mailto:[EMAIL PROTECTED] Sent: May 12, 2006 4:35 PM To: Chaouch, Aziz Subject: RE: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution look at the following code: library(copula) par(mfrow = c(2, 2)) x - mvdc(normalCopula(sin(0.5 * pi /2)), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) x - mvdc(frankCopula(5.736276), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) x - mvdc(gumbelCopula(2), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) x - mvdc(claytonCopula(2), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) the values of the association parameter I've chosen in each copula correspond to Kendall's tau 0.5; assuming also standard normal marginal distributions look at the different shapes you get! If possible try something similar for you case (i.e., using von Mises and Weibull marginals) and check if the association shape for a specific copula is more appropriate for your application. If this is not possible fit models assumig different copulas and check which one provides a better fit to your data. I hope it helps. Best, Dimitris Quoting Chaouch, Aziz [EMAIL PROTECTED]: Hi Dimitris, I'm not sure to understand your suggestion. How would you build that contour plot for a particular copula and on what is computed the Kendall's tau? Thanks, Aziz -Original Message- From: Dimitris Rizopoulos [mailto:[EMAIL PROTECTED] Sent: May 12, 2006 9:57 AM To: Chaouch, Aziz; [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Subject: Re: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution the choice of the copula is, in fact, a model selection problem. First, you could have a look at the contour plots of different copulas (preferably for the same value of Kendall's tau), and decide if some of them assume a more appropriate association structure for your application, compared to the others. Afterwards, you may fit various copula functions, check the fit on the data, compute AIC (since these are typically not nested models), etc. regarding the Von Mises distribution and if could be used in mvdc(), that I don't know. It'd be better to contact the copula package maintainer and ask. I hope it helps. Best, Dimitirs Dimitris Rizopoulos Ph.D. Student Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven, Belgium Tel: +32/(0)16/336899 Fax: +32/(0)16/337015 Web: http://www.med.kuleuven.be/biostat/ http://www.student.kuleuven.be/~m0390867/dimitris.htm - Original Message - From: Chaouch, Aziz [EMAIL PROTECTED] To: Dimitris Rizopoulos [EMAIL PROTECTED]; [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Sent: Friday, May 12, 2006 3:13 PM Subject: RE: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution Thanks a lot! I wasn't aware of that copula package and it could well be appropriate to use it for my application. However if I read the copula help correctly, I still need to know what kind of copula to use to link the distribution of wind speeds and directions. Is there a way to determine this in R? Moreover can I use the Von Mises distribution from the circular or CircStats package into the mvdc function of the copula package or does the mvdc function only recognize distributions available natively within R? Thanks again to all, your help is highly appreciated for a newbie like me! Regards, Aziz -Original Message- From: Dimitris Rizopoulos [mailto:[EMAIL PROTECTED] Sent: May 12, 2006 3:01 AM To: Philip He; Chaouch, Aziz Cc: r-help@stat.math.ethz.ch Subject: Re: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution - Original Message - From: Philip He [EMAIL PROTECTED] To: Chaouch, Aziz [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Sent: Thursday, May 11, 2006 11:21 PM Subject: Re: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution On 5/11/06, Chaouch, Aziz [EMAIL PROTECTED]
Re: [R] upgrade or recompile
Mike Wolfgang wrote: Dear list, My current R 2.2.0 in my linux desktop was installed by compiling sources, now I want to upgrade it to 2.3.0, is it possible to upgrade it without re-compiling sources? Thanks, No, except if installing a binary distribution. BTW: You might want to try an R-2.3.1 beta release instead anyway. Uwe Ligges mike [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Maximum likelihood estimate of bivariatevonmises-weibulldistribution
Hi Aziz, I am attaching a file that contains the functions that I wrote to compute the MLE for a binary vonMises mixture, using the EM algorithm. It also contains a simulation example. Let me know if you have any trouble. Hope this is helpful, Ravi. -- Ravi Varadhan, Ph.D. Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University Ph: (410) 502-2619 Fax: (410) 614-9625 Email: [EMAIL PROTECTED] Webpage:http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html -- -Original Message- From: [EMAIL PROTECTED] [mailto:r-help- [EMAIL PROTECTED] On Behalf Of Chaouch, Aziz Sent: Friday, May 26, 2006 12:02 PM To: Dimitrios Rizopoulos Cc: r-help@stat.math.ethz.ch Subject: Re: [R] Maximum likelihood estimate of bivariatevonmises- weibulldistribution Hi, I'm still strugling with this copula model but this seems to be the way to go. I'm now trying to model the marginal distributions and and for wind direction I use a mixture of 2 von mises. I'd like to estimate all the parameters (m1,m1,kappa1,kappa2,p) by maximizing the likelihood but I don't know how to define the likelihood (or log-likelihood) of a mixture of 2 Von Mises to use it with the function fitDistr in the MASS package. Can you help me define this likelihood function and use it through the fitDistr function? Thanks, Aziz -Original Message- From: Dimitrios Rizopoulos [mailto:[EMAIL PROTECTED] Sent: May 12, 2006 4:35 PM To: Chaouch, Aziz Subject: RE: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution look at the following code: library(copula) par(mfrow = c(2, 2)) x - mvdc(normalCopula(sin(0.5 * pi /2)), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) x - mvdc(frankCopula(5.736276), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) x - mvdc(gumbelCopula(2), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) x - mvdc(claytonCopula(2), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) the values of the association parameter I've chosen in each copula correspond to Kendall's tau 0.5; assuming also standard normal marginal distributions look at the different shapes you get! If possible try something similar for you case (i.e., using von Mises and Weibull marginals) and check if the association shape for a specific copula is more appropriate for your application. If this is not possible fit models assumig different copulas and check which one provides a better fit to your data. I hope it helps. Best, Dimitris Quoting Chaouch, Aziz [EMAIL PROTECTED]: Hi Dimitris, I'm not sure to understand your suggestion. How would you build that contour plot for a particular copula and on what is computed the Kendall's tau? Thanks, Aziz -Original Message- From: Dimitris Rizopoulos [mailto:[EMAIL PROTECTED] Sent: May 12, 2006 9:57 AM To: Chaouch, Aziz; [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Subject: Re: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution the choice of the copula is, in fact, a model selection problem. First, you could have a look at the contour plots of different copulas (preferably for the same value of Kendall's tau), and decide if some of them assume a more appropriate association structure for your application, compared to the others. Afterwards, you may fit various copula functions, check the fit on the data, compute AIC (since these are typically not nested models), etc. regarding the Von Mises distribution and if could be used in mvdc(), that I don't know. It'd be better to contact the copula package maintainer and ask. I hope it helps. Best, Dimitirs Dimitris Rizopoulos Ph.D. Student Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven, Belgium Tel: +32/(0)16/336899 Fax: +32/(0)16/337015 Web: http://www.med.kuleuven.be/biostat/ http://www.student.kuleuven.be/~m0390867/dimitris.htm - Original Message - From: Chaouch, Aziz [EMAIL PROTECTED] To: Dimitris Rizopoulos [EMAIL PROTECTED]; [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Sent: Friday, May 12, 2006 3:13 PM Subject: RE: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution Thanks a lot! I wasn't aware of that copula package and it could well be appropriate to use it for my application.
[R] how to pick a value from AR equation
i need to compute several (hundreds) of regression AR and PP.test for untary roots. is there an easy way to pick the values of interest from the output of these operations? i would like to fill a matrix (in a for cycle) with all the values of coefficients, standard deviations, PP statistics and the relevant critical value for the series i have. thanks in advance L -- credo nella ragione umana, e nella liberta' e nella giustizia che dalla ragione scaturiscono. (sciascia) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] I cannot load the package Rcmdr
Dear Zhang Jian, Is it possible that you're using a version of the Rcmdr package that's incompatible with your version of R? (I won't have access to email this weekend, so will be unable to reply until Monday to your response.) I hope this helps, John On Fri, 26 May 2006 09:52:31 -0400 zhang jian [EMAIL PROTECTED] wrote: I can load the package in other R file, but I can not do it in the file. How to solve the question? local({pkg - select.list(sort(.packages(all.available = TRUE))) + if(nchar(pkg)) library(pkg, character.only=TRUE)}) Loading required package: tcltk Loading required package: car Error in parse(file, n, text, prompt) : syntax error in * Error: .onAttach failed in 'attachNamespace' Error: package/namespace load failed for 'Rcmdr' [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html John Fox Department of Sociology McMaster University Hamilton, Ontario, Canada http://socserv.mcmaster.ca/jfox/ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] combinatorial programming problem
Hola! I am programming a class (S3) symarray for storing the results of functions symmetric in its k arguments. Intended use is for association indices for more than two variables, for instance coresistivity against antibiotics. There is one programming problem I haven't solved, making an inverse of the index function indx() --- se code below. It could for instance return the original k indexes in strictly increasing order, to make indx() formally invertible. Any ideas? Kjetil Code: # Implementing an S3 class for symarrays with array rank r for dimension # [k, k, ..., k] with k=r repeated r times. We do not store the diagonal. # Storage requirement is given by {r, k}= choose(k, r) # where r=array rank, k=maximum index symarray - function(data=NA, dims=c(1,1)){ r - dims[1] k - dims[2] if(r k) stop(symarray needs dimension larger than array rank) len - choose(k, r) out - data[1:len] attr(out, dims) - dims class(out) - symarray out } # Index calculation: indx - function(inds, k){ r - length(inds) if(r==1) return(inds) else { if(inds[1]==1) { return( indx(inds[-1]-1, k-1 ) ) } else { return( indx(c(inds[1]-1, seq(from=k-r+2, by=1, to=k)), k) + indx( inds[-1]-inds[1], k-inds[1] )) } } } # end indx # Methods for assignment and indexing: [.symarray - function(x, inds, drop=FALSE){ dims - attr(x, dims) k - dims[2] inds - indx(inds, k) res - NextMethod([, x) res } [-.symarray - function(x, inds, value){ dims - attr(x, dims) k - dims[2] inds - indx(inds, k) res - NextMethod([-, x) res } __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Returned mail: Data format error
Hello, Your mail to [EMAIL PROTECTED] was caught by the SpamAssassin filter running on the bcs.org.uk mail system. To confirm that your mail is genuine, please click this link, or paste it into your browser: https://bcsnet.bcs.org.uk/approve.php?c=66e3c8769c1611445312f4 You will not have to do this again for any mail sent to this recipient ([EMAIL PROTECTED]). Thank you. -- British Computer Society - www.bcs.org.uk Email Services from gradwell dot com - www.gradwell.com __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to pick a value from AR equation
Lorenzo, in order to help you more, you should consider sending to the list some relevant code with the functions and packages you are using. Rogerio. - Original Message - From: Lorenzo Bencivelli [EMAIL PROTECTED] To: r-help@stat.math.ethz.ch Sent: Friday, May 26, 2006 1:21 PM Subject: [R] how to pick a value from AR equation i need to compute several (hundreds) of regression AR and PP.test for untary roots. is there an easy way to pick the values of interest from the output of these operations? i would like to fill a matrix (in a for cycle) with all the values of coefficients, standard deviations, PP statistics and the relevant critical value for the series i have. thanks in advance L -- credo nella ragione umana, e nella liberta' e nella giustizia che dalla ragione scaturiscono. (sciascia) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to pick a value from AR equation
ok, i'll try to be more precise. both functions i need (ar() and PP.test()) are in the package stats. what i would like to do: omega-array(0, dim=c(N+1,7)) omega[1,]-(series,mean,std,max,min,ar(1)coeff,stdar(1)) for (i in (2:N+1)){ omega[i,1]-i-1 omega[i,2]-mean(y[(i-1),]) (...) omega[i,6]- (???) omega[i,7]- (???)} for the PP.test, i would like to do something similar. N is the number of the series for which i would like to compute the AR(1) model and the unitary root test. thanks in advance. L -- credo nella ragione umana, e nella liberta' e nella giustizia che dalla ragione scaturiscono. (sciascia) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] I cannot load the package Rcmdr
On Fri, 26 May 2006, John Fox wrote: Dear Zhang Jian, Is it possible that you're using a version of the Rcmdr package that's incompatible with your version of R? That ought to have generated an error much earlier, but of course the dependencies in an old version of Rcmdr might not have been prescient. I would suggest running update.packages() to be sure. (I won't have access to email this weekend, so will be unable to reply until Monday to your response.) I hope this helps, John On Fri, 26 May 2006 09:52:31 -0400 zhang jian [EMAIL PROTECTED] wrote: I can load the package in other R file, but I can not do it in the file. How to solve the question? local({pkg - select.list(sort(.packages(all.available = TRUE))) + if(nchar(pkg)) library(pkg, character.only=TRUE)}) Loading required package: tcltk Loading required package: car Error in parse(file, n, text, prompt) : syntax error in * Error: .onAttach failed in 'attachNamespace' Error: package/namespace load failed for 'Rcmdr' Please look at the section on debugging in the `Writing R Extensions' manual. (If you do not see it, you need to update your R to = 2.3.0.) At a minimum, running traceback() at that point would have helped both you and us. (I am guessing this occurs in the call to Commander(), but would like not to be guessing.) Also, please check out the posting guide and tell us the information requested, e.g. sessionInfo(), and Sys.getlocale() in case it is a charset issue. -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to get a series of results by loop
I donot know how to do it remainly. Perhaps my question is not clear. For example, test=data.frame(sp=paste(c(sp),1:100,sep=),x=rnorm(100),y=rnorm(100)) # I want to get the data in the circle of radius = 1 in every point, and I want to add a sign (e.g. sp1,sp2) in order that I can distinguish different groups # if I want to choose a point, I can do it like this r = 1 # radius attach(test) point100=test[(x-x[100])^2+(y-y[100])^2=r^2,] result100=data.frame (note=paste(c(p),rep(100,length(point100$sp)),sep=),point100) detach(test) But if I want to do all points in test, I donot know how to do it. I try to do it by using for, but it doesnot work. Can you give me some advices? thanks! jian zhang [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to pick a value from AR equation
Ciao, you didn't provide lots of details so I'll try to answer your questions with examples you can find typing ?PP.test and ?ar. In general, when you run a command in R you get an output. You just need to understand what kind of object you get. For the PP.test x - rnorm(1000) pp - PP.test(x) The command 'str' is very useful str(pp) List of 5 $ statistic: Named num -33.2 ..- attr(*, names)= chr Dickey-Fuller $ parameter: Named num 7 ..- attr(*, names)= chr Truncation lag parameter $ p.value : num 0.01 $ method : chr Phillips-Perron Unit Root Test $ data.name: chr x - attr(*, class)= chr htest It's a list. So you can select the element you need, e.g. pp$statistic Dickey-Fuller -33.22576 or as.numeric(pp$statistic) [1] -33.22576 so you get rid of the attributes. Now, same thing for 'ar'. fit.ar - ar(lh) str(fit.ar) List of 14 $ order : int 3 $ ar : num [1:3] 0.6534 -0.0636 -0.2269 $ var.pred: num 0.196 $ x.mean : num 2.4 $ aic : Named num [1:17] 18.307 0.996 0.538 0.000 1.490 ... ..- attr(*, names)= chr [1:17] 0 1 2 3 ... $ n.used : int 48 $ order.max : num 16 $ partialacf : num [1:16, 1, 1] 0.576 -0.223 -0.227 0.103 -0.076 ... $ resid : Time-Series [1:48] from 1 to 48: NA NA NA -0.200 -0.169 ... $ method : chr Yule-Walker $ series : chr lh $ frequency : num 1 $ call: language ar(x = lh) $ asy.var.coef: num [1:3, 1:3] 0.02156 -0.01518 0.00482 -0.01518 0.03117 ... - attr(*, class)= chr ar This is for what concerns your questions. I'd like to give you a hint for a simpler programming code. You use the first row of your matrix as header omega-array(0, dim=c(N+1,7)) omega[1,]-(series,mean,std,max,min,ar(1)coeff,stdar(1)) you can simply type omega-array(0, dim=c(N,7)) colnames(omega) - c(series,mean,std,max,min,ar(1)coeff,stdar(1)) so your matrix has Nx7 empty cells and avoid an headache trying to adjust the index 'i' in the loop 'for (i in (2:N+1))'. Also, if you really need to specify the 'id' of the record omega[i,1]-i-1 you just type before the loop omega[,1] - 1:N or rownames(omega) - 1:N In summary, omega-array(0, dim=c(N,7)) colnames(omega) - c(series,mean,std,max,min,ar(1)coeff,stdar(1)) omega[,1] - 1:N for (i in 1:N){ omega[i,2]-mean(y[i,]) (...) omega[i,6]- (???) # see before omega[i,7]- (???)} # see before Finally, if 'y' is a prespecified matrix I suggest you to take a look at ?rowMeans, ?rowSums, and ?apply in order to execute 'mean' or 'sd' or any other function applied to 'y[i,]' at once instead of N times Hope this helps, Marco Geraci --- Lorenzo Bencivelli [EMAIL PROTECTED] wrote: ok, i'll try to be more precise. both functions i need (ar() and PP.test()) are in the package stats. what i would like to do: omega-array(0, dim=c(N+1,7)) omega[1,]-(series,mean,std,max,min,ar(1)coeff,stdar(1)) for (i in (2:N+1)){ omega[i,1]-i-1 omega[i,2]-mean(y[(i-1),]) (...) omega[i,6]- (???) omega[i,7]- (???)} for the PP.test, i would like to do something similar. N is the number of the series for which i would like to compute the AR(1) model and the unitary root test. thanks in advance. L -- credo nella ragione umana, e nella liberta' e nella giustizia che dalla ragione scaturiscono. (sciascia) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] find position to last number in a vector
Hello Is there a function to extract the position (i.e. row number or more desirable the value from another column with the same rownumber ) of last number in a vector? This is done in a loop with 1000s of columns. e.g. vector/column.n: na,na,10,na,2,na,na,1,na,na fixed column:10,20,30,40,50,60,70,80,85,100 result: 8 and/or 80 Thanks... Best Regards Anders __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] curve peeling
R2.3.0, WinXP. I am trying to fit the model Y(t) = C1*exp(-exp(a)*t) + C2*exp(-exp(b)*t) - (C1+C2)*exp(-exp(c)*t) to some pharmacokinetic data. The data derives from a single oral dose of a drug. Does anyone know how to use peeling and/or other methods to obtain good starting estimates of the parameters in the triexponential? The book Mixed Effects Models in S and SPLUS describes self-starting models for the biexponential and first-order one-compartment PK models, but not for the triexponential. Much thanks in advance, Greg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] find position to last number in a vector
Is this what you want? x [1] NA NA 10 NA 2 NA NA 1 NA NA y [1] 10 20 30 40 50 60 70 80 85 100 y[max(which(!is.na(x)))] [1] 80 Andy From: Anders Bjørgesæter Hello Is there a function to extract the position (i.e. row number or more desirable the value from another column with the same rownumber ) of last number in a vector? This is done in a loop with 1000s of columns. e.g. vector/column.n: na,na,10,na,2,na,na,1,na,na fixed column:10,20,30,40,50,60,70,80,85,100 result: 8 and/or 80 Thanks... Best Regards Anders __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] pdf output incompatible with latest ghostscript
I just upgraded to the latest ghostscript, AFPL 8.54. When I run the R output through gs, I get AFPL Ghostscript 8.54: Set UseCUEColor for UseDeviceIndependentColor to work properly. followed by a segfault. the resulting ghostscript file is incomplete (no %%EOF). probably a ghostscript bug, but caused by R output. regards, /iaw __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Maximum likelihood estimate of bivariatevonmises-weibulldistribution
Thanks Ravi! That's probably more than what I need. I suppose I'll have to get back to you for some more explanations but for now I'm trying to apply the log-likelihood (likelihood) that you defined to the fitdistr function or use the log-likelihood with optim(). Thanks again, that's really helpful! Regards, Aziz -Original Message- From: Ravi Varadhan [mailto:[EMAIL PROTECTED] Sent: May 26, 2006 12:18 PM To: Chaouch, Aziz Cc: r-help@stat.math.ethz.ch Subject: RE: [R] Maximum likelihood estimate of bivariatevonmises-weibulldistribution Hi Aziz, I am attaching a file that contains the functions that I wrote to compute the MLE for a binary vonMises mixture, using the EM algorithm. It also contains a simulation example. Let me know if you have any trouble. Hope this is helpful, Ravi. -- Ravi Varadhan, Ph.D. Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University Ph: (410) 502-2619 Fax: (410) 614-9625 Email: [EMAIL PROTECTED] Webpage:http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html -- -Original Message- From: [EMAIL PROTECTED] [mailto:r-help- [EMAIL PROTECTED] On Behalf Of Chaouch, Aziz Sent: Friday, May 26, 2006 12:02 PM To: Dimitrios Rizopoulos Cc: r-help@stat.math.ethz.ch Subject: Re: [R] Maximum likelihood estimate of bivariatevonmises- weibulldistribution Hi, I'm still strugling with this copula model but this seems to be the way to go. I'm now trying to model the marginal distributions and and for wind direction I use a mixture of 2 von mises. I'd like to estimate all the parameters (m1,m1,kappa1,kappa2,p) by maximizing the likelihood but I don't know how to define the likelihood (or log-likelihood) of a mixture of 2 Von Mises to use it with the function fitDistr in the MASS package. Can you help me define this likelihood function and use it through the fitDistr function? Thanks, Aziz -Original Message- From: Dimitrios Rizopoulos [mailto:[EMAIL PROTECTED] Sent: May 12, 2006 4:35 PM To: Chaouch, Aziz Subject: RE: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution look at the following code: library(copula) par(mfrow = c(2, 2)) x - mvdc(normalCopula(sin(0.5 * pi /2)), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) x - mvdc(frankCopula(5.736276), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) x - mvdc(gumbelCopula(2), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) x - mvdc(claytonCopula(2), c(norm, norm), list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7)) the values of the association parameter I've chosen in each copula correspond to Kendall's tau 0.5; assuming also standard normal marginal distributions look at the different shapes you get! If possible try something similar for you case (i.e., using von Mises and Weibull marginals) and check if the association shape for a specific copula is more appropriate for your application. If this is not possible fit models assumig different copulas and check which one provides a better fit to your data. I hope it helps. Best, Dimitris Quoting Chaouch, Aziz [EMAIL PROTECTED]: Hi Dimitris, I'm not sure to understand your suggestion. How would you build that contour plot for a particular copula and on what is computed the Kendall's tau? Thanks, Aziz -Original Message- From: Dimitris Rizopoulos [mailto:[EMAIL PROTECTED] Sent: May 12, 2006 9:57 AM To: Chaouch, Aziz; [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Subject: Re: [R] Maximum likelihood estimate of bivariate vonmises-weibulldistribution the choice of the copula is, in fact, a model selection problem. First, you could have a look at the contour plots of different copulas (preferably for the same value of Kendall's tau), and decide if some of them assume a more appropriate association structure for your application, compared to the others. Afterwards, you may fit various copula functions, check the fit on the data, compute AIC (since these are typically not nested models), etc. regarding the Von Mises distribution and if could be used in mvdc(), that I don't know. It'd be better to contact the copula package maintainer and ask. I hope it helps. Best, Dimitirs Dimitris Rizopoulos Ph.D. Student Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven,
Re: [R] find position to last number in a vector
I am sure there are several ways of doing this but how about something like this: y - c(NA,NA,10,NA,2,NA,NA,1,NA,NA) x - c(10,20,30,40,50,60,70,80,85,100) tail(x[!is.na(y)], 1) or rev(x[!is.na(y)])[1] Rudimentary checks indicate that the second expression is much faster than the first. Also, if y is contains NAs only, the first expression returns numeric(0) and the second one NA. Hope this helps, Andy __ Andy Jaworski 518-1-01 Process Laboratory 3M Corporate Research Laboratory - E-mail: [EMAIL PROTECTED] Tel: (651) 733-6092 Fax: (651) 736-3122 Anders Bjørgesæter anders.bjorgesat To [EMAIL PROTECTED]r-help@stat.math.ethz.ch Sent by: cc [EMAIL PROTECTED] at.math.ethz.ch Subject [R] find position to last number in a vector 05/26/2006 01:30 PM Hello Is there a function to extract the position (i.e. row number or more desirable the value from another column with the same rownumber ) of last number in a vector? This is done in a loop with 1000s of columns. e.g. vector/column.n: na,na,10,na,2,na,na,1,na,na fixed column:10,20,30,40,50,60,70,80,85,100 result: 8 and/or 80 Thanks... Best Regards Anders __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to get a series of results by loop
Try this: set.seed(1) r - .1 n - 5 test - data.frame(sp = paste(sp, 1:n, sep = ), x= rnorm(n), y = rnorm(n)) r - 1 # radius f - function(i) cbind(orig = test[i,1], subset(test, (x - x[i])^2 + (y - y[i])^2 r^2)) do.call(rbind, lapply(1:n, f)) On 5/26/06, zhang jian [EMAIL PROTECTED] wrote: I donot know how to do it remainly. Perhaps my question is not clear. For example, test=data.frame(sp=paste(c(sp),1:100,sep=),x=rnorm(100),y=rnorm(100)) # I want to get the data in the circle of radius = 1 in every point, and I want to add a sign (e.g. sp1,sp2) in order that I can distinguish different groups # if I want to choose a point, I can do it like this r = 1 # radius attach(test) point100=test[(x-x[100])^2+(y-y[100])^2=r^2,] result100=data.frame (note=paste(c(p),rep(100,length(point100$sp)),sep=),point100) detach(test) But if I want to do all points in test, I donot know how to do it. I try to do it by using for, but it doesnot work. Can you give me some advices? thanks! jian zhang [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] PC rotation question
I wrote: On p. 48 of Statistics Complements to the 3rd MASS edition, http://www.stats.ox.ac.uk/pub/MASS3/VR3stat.pdf I read that the orthogonal rotations of Z Lambda^-1 remain uncorrelated, where Z is the PC and Lambda is the diag matrix of singular values. However, the example below that text is A - loadings(ir.pca) %*% diag(ir.pca$sdev) If ir.pca$sdev are the singular values, should that be diag(1 / ir.pca$sdev), or is it some discrepancy between S+ and R that I'm missing? To dwell some more on this and in hopes to get replies I did a small experimantation (below) that makes me to suspect that the correct syntax is A - loadings(ir.pca) %*% diag(1/ir.pca$sdev). Any comments? library(MASS) a - 1; b - 0.2; c - -.3; d - .7 x - scale(mvrnorm(n=1000,c(0,0,0),matrix(c(a,b,c, b,a,d, c,d,a),3,3))) e - eigen(cov(x)) cat(expect identity corr for orthogonal rotations:\n) zs - e$vectors %*% diag(1/sqrt(e$values)) pdx - x %*% varimax(zs, normalize = FALSE)$loadings print.table(cor(pdx), digits=2) cat(expect zero: , prcomp(x)$sdev - sqrt(e$values), \n) cat(Now zero corr between projections is not preserved:\n) zs - e$vectors %*% diag(sqrt(e$values)) pdx - x %*% varimax(zs, normalize = FALSE)$loadings print.table(cor(pdx), digits=2) Output: expect identity corr for orthogonal rotations: [,1] [,2] [,3] [1,] 1.0e+00 -3.0e-16 -4.5e-16 [2,] -3.0e-16 1.0e+00 2.3e-16 [3,] -4.5e-16 2.3e-16 1.0e+00 expect zero: 4.440892e-16 2.220446e-16 2.775558e-16 Now zero corr between projections is not preserved: [,1] [,2] [,3] [1,] 1.00 -0.35 -0.88 [2,] -0.35 1.00 -0.11 [3,] -0.88 -0.11 1.00 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] missed ylim from plot.default
Try the following function to see if it does what you want. The basic syntax for your example would be: tmp - squishplot( xlim=c(-0.5,7), ylim=c(-0.5,2.8), asp=1 ) plot(1:5, rep(1,5), xlim=c(-0.5,7), ylim=c(-0.5,2.8)) par(tmp) # reset plotting region for future plots The function definition is: squishplot - function(xlim,ylim,asp=1){ if(length(xlim) 2) stop('xlim must be a vector of length 2') if(length(ylim) 2) stop('ylim must be a vector of length 2') tmp - par(c('plt','pin','xaxs','yaxs')) if( tmp$xaxs == 'i' ){ # not extended axis range xlim - range(xlim) } else { # extended range tmp.r - diff(range(xlim)) xlim - range(xlim) + c(-1,1)*0.04*tmp.r } if( tmp$yaxs == 'i' ){ # not extended axis range ylim - range(ylim) } else { # extended range tmp.r - diff(range(ylim)) ylim - range(ylim) + c(-1,1)*0.04*tmp.r } tmp2 - (ylim[2]-ylim[1])/(xlim[2]-xlim[1]) tmp.y - tmp$pin[1] * tmp2 * asp if(tmp.y tmp$pin[2]){ # squish vertically par(pin=c(tmp$pin[1], tmp.y)) par(plt=c(tmp$plt[1:2], par('plt')[3:4])) } else { # squish horizontally tmp.x - tmp$pin[2]/tmp2/asp par(pin=c(tmp.x, tmp$pin[2])) par(plt=c(par('plt')[1:2], tmp$plt[3:4])) } return(invisible(tmp['plt'])) } Let me know of any improvments you can think of (including a better name) and how it works. -- Gregory (Greg) L. Snow Ph.D. Statistical Data Center Intermountain Healthcare [EMAIL PROTECTED] (801) 408-8111 -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Antonio, Fabio Di Narzo Sent: Friday, May 26, 2006 1:29 AM To: Peter Ehlers Cc: R-help@stat.math.ethz.ch Subject: Re: [R] missed ylim from plot.default 2006/5/26, Peter Ehlers [EMAIL PROTECTED]: Antonio, Fabio Di Narzo wrote: Hi all. On my R-2.3.0, calling: plot(1:5, rep(1,5), xlim=c(-0.5,7), ylim=c(-0.5,2.8), asp=1) gives to me a plot (tried pdf and X11) whose y axis goes from about -2 to about 4.5. What I've missed? How can I show some pre-specified x and y ranges from a plot (keeping fixed the aspect ratio)? Tnx all, Antonio, Fabio Di Narzo. [[alternative HTML version deleted]] Do you want a wide, but not very tall plot? If so, you could specify the width/height in your X11() call. Peter Ehlers Tnx for your indication, now I see: x and y plotting ranges are influenced by the actual window size (at least when fixing aspect ratio). I think this is obvious... What I missed/forgot above, was that default plot width and height are fixed indipendently from the 'x-ylim' and 'asp' arguments. However now I have a problem: I can write a wrapper to set device width and height depending on the (xlim,ylim,asp) triple. But there's often extra-space in the device due to plot title, axis labels, etc. So, I fear I can't simply set width=height*asp. Some suggestion? (Anyway, by now I only have to produce few plots, so today I can go try by try :-) ). Bests, Antonio. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] R.oo question
This is a simple R.oo question but I, thankfully, hope that someone would explain it to me so I would better understand this work frame. I create this class: setConstructorS3(MyExample, function(param=0) { print(paste(called with param=, param)) extend(Object(), MyExample, .param = param ); }) From what is printed out, who made the second call to the class with the default param? MyExample(1) [1] called with param= 1 [1] called with param= 0# - this is line in question [1] MyExample: 0x02831708 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] R newbie attempting to plot data
Hi, I just started using R and am having trouble with the below error: I type: df - read.csv(/home/rex/Desktop/mytable.csv) which gives me what I want: ... 639 2006-05-26 16:46:54 4 16 640 2006-05-26 17:05:36 5 17 641 2006-05-26 17:30:48 6 17 But now I try: plot(df[4],df[3]) Error in pmatch(x, table, duplicates.ok) : argument is not of mode character Any ideas? I tried read.table and importing the data from MySQL too; same error. Thanks, Rex [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R newbie attempting to plot data
Rex Eastbourne wrote: Hi, I just started using R and am having trouble with the below error: I type: df - read.csv(/home/rex/Desktop/mytable.csv) which gives me what I want: ... 639 2006-05-26 16:46:54 4 16 640 2006-05-26 17:05:36 5 17 641 2006-05-26 17:30:48 6 17 But now I try: plot(df[4],df[3]) Error in pmatch(x, table, duplicates.ok) : argument is not of mode character You probably intended: plot(df[,4],df[,3]) Any ideas? I tried read.table and importing the data from MySQL too; same error. Thanks, Rex [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R.oo question
On 5/26/06, Omar Lakkis [EMAIL PROTECTED] wrote: This is a simple R.oo question but I, thankfully, hope that someone would explain it to me so I would better understand this work frame. I create this class: setConstructorS3(MyExample, function(param=0) { print(paste(called with param=, param)) extend(Object(), MyExample, .param = param ); }) From what is printed out, who made the second call to the class with the default param? MyExample(1) [1] called with param= 1 [1] called with param= 0# - this is line in question [1] MyExample: 0x02831708 That is because of a new US government rule requiring that one instance of every R.oo class is forwarded to the NSA. Seriously, you will see that this happens once and only once for each class defined this way and it normally happens when you create your first instance (otherwise before that). The first call to the constructor creates a so called static instance of the class. The static instance is for instance used when you type 'MyExample' to get the API coupled to class MyExample. Another example: Object Object public $(name) public $-(name, value) public [[(name) ... public objectSize(...) public print(...) public save(file=NULL, path=NULL, compress=TRUE, ...) } Technical details: The static instance is not always needed, but quite often. If needed, it has to be create at some stage and I found that having extend() to do this is quite convenient. The alternative would be to let you do it explicitly, e.g. getStaticInstance(Object). Note that there is no central registry/database keeping track of the classes defined by R.oo/Object, but all is just plain S3 classes making use of the S3/UseMethod dispatching mechanisms - that's all. There is a low-level way to figure out if the call to the constructor is for the static instance or not, but normally it is easier not to output anything in the constructor. If you want to know the low-level way, I'll have have to get back to, because I don't know it off the top of my head. I might provide a method for this if there is a need for it, e.g. hasStaticInstance(MyExample). Cheers Henrik (the author) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Too many open files
This may be more of an OS question ... I have this call r = get.hist.quote(symbol, start= format(start, %Y-%m-%d), end= format(end, %Y-%m-%d)) which does a url request in a loop and my program runs out of file handlers after few hundred rotations. The error message is: 'Too many open files'. Other than increasing the file handlers assigned to my process, is there a way to cleanly release and reuse these connections? I run on Mac OS X __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Too many open files
Omar Lakkis [EMAIL PROTECTED] writes: This may be more of an OS question ... I have this call r = get.hist.quote(symbol, start= format(start, %Y-%m-%d), end= format(end, %Y-%m-%d)) which does a url request in a loop and my program runs out of file handlers after few hundred rotations. The error message is: 'Too many open files'. Other than increasing the file handlers assigned to my process, is there a way to cleanly release and reuse these connections? Inside your loop you need to close the connection object created by url(). for (i in 1:500) { con - url(urls[i]) ## ... stuff here ... close(con) } R only allows you to have a fixed number of open connections at one time, and they do not get closed automatically when they go out of scope. These commands may help make clear how things work... showConnections() description class mode text isopen can read can write f = url(http://www.r-project.org;, open=r) showConnections() descriptionclass mode text isopen can read can write 3 http://www.r-project.org; url r text opened yesno rm(f) showConnections() descriptionclass mode text isopen can read can write 3 http://www.r-project.org; url r text opened yesno f - getConnection(3) close(f) f Error in summary.connection(x) : invalid connection showConnections() description class mode text isopen can read can write + seth __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Too many open files
Try closeAllConnections() On 5/26/06, Omar Lakkis [EMAIL PROTECTED] wrote: This may be more of an OS question ... I have this call r = get.hist.quote(symbol, start= format(start, %Y-%m-%d), end= format(end, %Y-%m-%d)) which does a url request in a loop and my program runs out of file handlers after few hundred rotations. The error message is: 'Too many open files'. Other than increasing the file handlers assigned to my process, is there a way to cleanly release and reuse these connections? I run on Mac OS X __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Building V2.3.0 on Tru64 V5.1B
Hi, I am going to be attempting a build of R for the alphas (our main need is to collate and present graphically, performance metrics we gather on our applications). We have: CC=cc, CXX=cxx, Fortran V5.51 (f77, f95), GNU make 3.80, Perl 5.8.4 I am hesitant to rebuild GCC on this platform, and would like to know if others have success in building R with native C, C++ and Fortran. I plan on using the following configure string: CPPFLAGS=-I/opt/gnu/include -I/opt/libjpeg/include -I/opt/libpng/include -I/opt/ncurses/include CFLAGS=-std1 LDFLAGS=-L/opt/gnu/lib -L/opt/libjpeg/lib -L/opt/libpng/lib -L/opt/ncurses/lib -lncurses -ljpeg -lpng -rpath /opt/libpng/lib:/opt/libjpeg/lib:/opt/ncurses/lib:/opt/gnu/lib:/usr/shl ib R_PAPERSIZE=letter CC=cc CXX=cxx ./configure prefix=/opt/R --without-readline --without-blas I expect to update the config.site file with the above values. If someone has built R-2.3.0 on Tru64 V5.1B, I would appreciate hearing how. Thank you, Narendra Ravi __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Trouble passing list or non-list to function using ...
Hello, Simply put, I'm trying to call a function testme with value age=NA. I wish to use dotlist-list(...) inside the function and have dotlist become: $age [1] NA I'm modifying existing code and need to minimize changing that code so it's easiest to conform how I call the existing function. My sample code fragment, results, and R.version information are listed below. I've been searching for prior questions regarding this phenomena, but haven't quite figured out how to do this. Thank you in advance, -jason ##BEGIN TEST CODE ### testme-function(...){ dotlist-list(...) dotlist } nm-age testme(age=NA) testme(nm=NA) testme(age=NA,age2=NA) testlist-list(age,age2) for (i in 1:length(testlist)){ print(testme(testlist[i])) } ##END TEST CODE # ##BEGIN RESULTS # testme-function(...){ + + dotlist-list(...) + + dotlist + + } nm-age testme(age=NA) *This is what I'm really after $age [1] NA testme(nm=NA) $nm [1] NA #trying w/ 2 vars testme(age=NA,age2=NA) $age [1] NA $age2 [1] NA testlist-list(age,age2) for (i in 1:length(testlist)){ + + print(testme(testlist[i])) + + } [[1]] [[1]][[1]] [1] age [[1]] [[1]][[1]] [1] age2 ##END RESULTS ### platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major 2 minor 3.0 year 2006 month 04 day24 svn rev37909 language R version.string Version 2.3.0 (2006-04-24) [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Trouble passing list or non-list to function using ...
Its not clear from your post which of the tests does not give you the desired output and what the desired output is in that case. Just guessing but maybe you meant to use testlist[[i]] instead of testlist[i] in the loop? On 5/26/06, Jason Barnhart [EMAIL PROTECTED] wrote: Hello, Simply put, I'm trying to call a function testme with value age=NA. I wish to use dotlist-list(...) inside the function and have dotlist become: $age [1] NA I'm modifying existing code and need to minimize changing that code so it's easiest to conform how I call the existing function. My sample code fragment, results, and R.version information are listed below. I've been searching for prior questions regarding this phenomena, but haven't quite figured out how to do this. Thank you in advance, -jason ##BEGIN TEST CODE ### testme-function(...){ dotlist-list(...) dotlist } nm-age testme(age=NA) testme(nm=NA) testme(age=NA,age2=NA) testlist-list(age,age2) for (i in 1:length(testlist)){ print(testme(testlist[i])) } ##END TEST CODE # ##BEGIN RESULTS # testme-function(...){ + + dotlist-list(...) + + dotlist + + } nm-age testme(age=NA) *This is what I'm really after $age [1] NA testme(nm=NA) $nm [1] NA #trying w/ 2 vars testme(age=NA,age2=NA) $age [1] NA $age2 [1] NA testlist-list(age,age2) for (i in 1:length(testlist)){ + + print(testme(testlist[i])) + + } [[1]] [[1]][[1]] [1] age [[1]] [[1]][[1]] [1] age2 ##END RESULTS ### platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major 2 minor 3.0 year 2006 month 04 day24 svn rev37909 language R version.string Version 2.3.0 (2006-04-24) [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html