[R] predicting values into the future
Hi: I have usually used the GROWTH() excel function to do this but now want to see if I can do this with R. I want to predict values into the future, possibly with the predict.arima Function. I have the following weekly fish weight averages: weight - c(2.1,2.4,2.8,3.6,4.1,5.2,6.3) week - c(1,2,3,4,5,6,7) I would like to predict what the weight will be by week 10 based on my weight values and make a line plot of all the weights(including the predicted values). I have two questions: 1- Should the predicted values be linear or exponential? 2- Is the predict.arima function appropriate to do this? Thanks in advance. Felipe D. Carrillo Supervisory Fishery Biologist Department of the Interior US Fish Wildlife Service California, USA __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Question about arima.sim()
Hi, I tried to simulate an ARIMA model by using arima.sim(), say arima.sim(n=100,list(order=c(1,0,1),ar=0.6,ma=0.9,sd=1), but the acf and pacf of simulated data using acf() and pacf() are so much different from the theoritcal acf and pacf. For instance, in my case, ar=0.6 and ma=0.9, so the acf for all lags should be greater than 0 based on the theoritical calculation, but the acf of simulated data always has some negative values which bothers me. Can any one tell me why? and how can I get the simulated data close to what I expect? Thanks very much JC _ [[elided Hotmail spam]] plorer 8. [[elided Hotmail spam]] [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] predicting values into the future
Here is a bit of an exploration of your data but first a couple of notes. * the information about Excel is probably a bit superfluous here. Some of us have no idea about Excel, and rather hope it can stay that way. * With such a short series, you don't stand much chance of fitting a time series model such as with arima. It's clearly not stationary, too. If you had multiple growth curves you may stand some chance of fitting a correlated model, but with just one, I don't think so. For now, I think you just may have to make the hopeful assumption of independence. You might like to look at this. weightData - data.frame(weight = c(2.1,2.4,2.8,3.6,4.1,5.2,6.3), week = 1:7) plot(weight ~ week, weightData) plot(log(weight) ~ week, weightData) ### clearly the log plot seems to linearise things. ### Try an non-linear regression: wModel - nls(weight ~ alpha + beta*exp(gamma*week), weightData, start = c(alpha = 0.0, beta = 1, gamma = 0.2), trace = TRUE) you should look at the residuals from this to see if the assumptions look reasonable. With only 7, you can't see much, though. now suppose you want to predict for another 3 weeks: newData - data.frame(week = 1:10) newData$pweight - predict(wModel, newData) plot(pweight ~ week, newData, pch = 4, col = red, ylab = Weight, xlab = Week) with(weightData, points(week, weight)) looks OK to me (thought fish cannot keep on growing exponentailly forever - this is clearly a model with limitations and you have to be careful when pushing it too far). finally predict on a more continuous scale and add in the result as a blue line. lData - data.frame(week = seq(1, 10, len = 1000)) with(lData, lines(week, predict(wModel, lData), col = blue)) Now that we have over-analysed this miniscule data set to blazes, perhaps it's time for a beer! __ Bill Venables http://www.cmis.csiro.au/bill.venables/ -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Felipe Carrillo Sent: Sunday, 5 April 2009 4:13 PM To: r-h...@stat.math.ethz.ch Subject: [R] predicting values into the future Hi: I have usually used the GROWTH() excel function to do this but now want to see if I can do this with R. I want to predict values into the future, possibly with the predict.arima Function. I have the following weekly fish weight averages: weight - c(2.1,2.4,2.8,3.6,4.1,5.2,6.3) week - c(1,2,3,4,5,6,7) I would like to predict what the weight will be by week 10 based on my weight values and make a line plot of all the weights(including the predicted values). I have two questions: 1- Should the predicted values be linear or exponential? 2- Is the predict.arima function appropriate to do this? Thanks in advance. Felipe D. Carrillo Supervisory Fishery Biologist Department of the Interior US Fish Wildlife Service California, USA __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Filling in Gapped Interval of a Data Frame As Series
Dear all, I have a data frame that looks like this: dat V1 V2 1 -8 100 2 -6 20.2 3 -1 1.5 4 2 3.00 5 3 78.8 6 5 33.2 7 6 44.5 8 7 5.00 Now I want to fill the V1 column in series (note that it is gapped), with the corresponding value in V2 as 0. So in the end we would like to have: newdat V1 V2 1 -8 100 2 -6 20.2 3 -7 0 4 -6 0 5 -5 0 6-4 0 7-3 0 8-2 0 9-1 1.5 100 0 111 0 12 2 3.00 13 3 78.8 14 4 0 15 5 33.2 16 6 44.5 17 7 5.00 Is there a way to do it in R? - Gundala Viswanath Jakarta - Indonesia __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Filling in Gapped Interval of a Data Frame As Series
one way is the following: dat. - data.frame(V1 = min(dat$V1):max(dat$V1), V2 = 0) newdat - merge(dat, dat., by = V1, all.y = TRUE, suffixes = c(, .)) newdat$V2[na.ind] - newdat$V2.[na.ind - is.na(newdat$V2)] newdat[-3] I hope it helps. Best, Dimitris Gundala Viswanath wrote: Dear all, I have a data frame that looks like this: dat V1 V2 1 -8 100 2 -6 20.2 3 -1 1.5 4 2 3.00 5 3 78.8 6 5 33.2 7 6 44.5 8 7 5.00 Now I want to fill the V1 column in series (note that it is gapped), with the corresponding value in V2 as 0. So in the end we would like to have: newdat V1 V2 1 -8 100 2 -6 20.2 3 -7 0 4 -6 0 5 -5 0 6-4 0 7-3 0 8-2 0 9-1 1.5 100 0 111 0 12 2 3.00 13 3 78.8 14 4 0 15 5 33.2 16 6 44.5 17 7 5.00 Is there a way to do it in R? - Gundala Viswanath Jakarta - Indonesia __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Dimitris Rizopoulos Assistant Professor Department of Biostatistics Erasmus University Medical Center Address: PO Box 2040, 3000 CA Rotterdam, the Netherlands Tel: +31/(0)10/7043478 Fax: +31/(0)10/7043014 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] axis colours
Thanks a lot. That's what I was looking for. Must have missed the fg bit in the ? par. Cheers, Umesh On Sun, Apr 5, 2009 at 3:19 AM, Duncan Murdoch murd...@stats.uwo.ca wrote: On 04/04/2009 5:13 PM, Umesh Srinivasan wrote: Hi, Is there a way to use par to change the colours of the axes lines (not the labels)? I've done this: par(bg = black) plot(x, y, col = yellow, pch = 16) but I have to use axes = F within the plot command, and then use axis (1, col = yellow) axis (2, col = yellow) and so on for axes 3 4. This does not help since the new axes are not of the same length as the default 'box' that you get around a plot. You can use box() to draw a box. So you can get yellow axes and box using plot(1,1, axes=F) axis(1, col=yellow) axis(2, col=yellow) box(col=yellow) par() doesn't have an option for this. plot() will use the default par(col) setting for the box, but that also affects the points being drawn. Duncan Murdoch [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Which model to keep (negative BIC)
Quoting cladoo.26 cladoo...@laposte.net: Hi, My questions concern the function 'mclustBIC' which compute BIC for a range of clusters of several models on the given data and the other function 'mclustModel' which choose the best model and the best number of cluster accordind to the results of the previous cited function. 1) When trying the following example (see ?mclustModel), I get negative BIC computed by 'mclustBIC', and the best model according to the results of 'mclustModel' is the one with the highest BIC (i.e. the closer to zero). irisBIC - mclustBIC(iris[,-5]) plot(irisBIC) mclustModel(iris[,-5], irisBIC) Because I don't find anything about this point, could someone confirm that when the BIC are positive, we try to the minimize the criterion (the model with the smallest BIC is the best one) but when the BIC are negative we look for the higher BIC (the model with a the BIC closest to zero is the best one) ? The mclust package seems to be using a definition of BIC that is the negative of the usual one, i.e. the bic() function in the mclust package returns 2 * loglik - nparams * log(n) where loglik is the log likelihood, n is the number of observations and nparams is the number of parameters. BIC is normally defined as -2 * loglik + nparams * log(n) and the optimal model is the one with the minimum BIC. However in this case, you want to maximize it. 2) Does the $G argument from the output of 'mclustModel' represent the best number of clusters for the chosen model ? According to the documentation it does, and you can verify from your plot that the VEV model with 2 components has maximum BIC Many thanks, this is my first post on R help, but I often consult the forum for 4 years. Cladoo --- This message and its attachments are strictly confidenti...{{dropped:8}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Eclipse and StatET Howto (also added Subversion, Rtools)
G'day Peter, On Sun, 05 Apr 2009 11:26:40 +0200 Peter Dalgaard p.dalga...@biostat.ku.dk wrote: Berwin A Turlach wrote: G'day Dirk, On Sat, 4 Apr 2009 20:27:22 -0500 Dirk Eddelbuettel e...@debian.org wrote: On 4 April 2009 at 14:37, Ken-JP wrote: Yes, I have x11-common installed, and dpkg -S /etc/X11/rgb.txt shows not found for me. This is on Ubuntu 8.10 amd64. Same 8.10 for both amd64 and i386 where I checked -- both have the file. It could be a leftover from an earlier install of mine, or an accidental deletion at your end. My guess that it is the former. :) Some time ago I wiped Debian of my laptop and installed Kubuntu 8.10 fresh (i386 flavour); so no left overs from previous versions and on my laptop the results are: ber...@berwin5:~$ apt-file find /etc/X11/rgb.txt ber...@berwin5:~$ dpkg -S /etc/X11/rgb.txt dpkg: /etc/X11/rgb.txt not found. Hum.. Fedora 9 doesn't have it either. It does have /usr/share/X11/rgb.txt though, so please check whether it has moved. Apparently it has: ber...@berwin5:~$ apt-file find rgb.txt [...] x11-common: /usr/share/X11/rgb.txt [...] However, ber...@berwin5:~$ ls -l /usr/share/X11/rgb.txt lrwxrwxrwx 1 root root 16 Jan 14 03:28 /usr/share/X11/rgb.txt - /etc/X11/rgb.txt ber...@berwin5:~$ more /usr/share/X11/rgb.txt /usr/share/X11/rgb.txt: No such file or directory Cheers, Berwin __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Eclipse and StatET Howto (also added Subversion, Rtools)
Berwin A Turlach wrote: G'day Dirk, On Sat, 4 Apr 2009 20:27:22 -0500 Dirk Eddelbuettel e...@debian.org wrote: On 4 April 2009 at 14:37, Ken-JP wrote: Yes, I have x11-common installed, and dpkg -S /etc/X11/rgb.txt shows not found for me. This is on Ubuntu 8.10 amd64. Same 8.10 for both amd64 and i386 where I checked -- both have the file. It could be a leftover from an earlier install of mine, or an accidental deletion at your end. My guess that it is the former. :) Some time ago I wiped Debian of my laptop and installed Kubuntu 8.10 fresh (i386 flavour); so no left overs from previous versions and on my laptop the results are: ber...@berwin5:~$ apt-file find /etc/X11/rgb.txt ber...@berwin5:~$ dpkg -S /etc/X11/rgb.txt dpkg: /etc/X11/rgb.txt not found. Hum.. Fedora 9 doesn't have it either. It does have /usr/share/X11/rgb.txt though, so please check whether it has moved. I'm curious as to why (only) R/tcltk would be confused by this sort of thing. -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - (p.dalga...@biostat.ku.dk) FAX: (+45) 35327907 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] threshold distribution
On Sun, 2009-04-05 at 01:13 -0400, jim holtman wrote: Here is what I get from using 'fitdistr' in R to fit to a lognormal. The resulting density plot from the distribution seems to be a reason match to the data. x - scan() 1: 0.80010 0.72299 0.69893 0.99597 0.89200 0.69312 0.73613 1.13559 9: 0.85009 0.85804 0.73324 1.04826 0.84002 14: 1.76330 0.71980 0.89416 0.89450 0.98670 0.83571 0.73833 0.66549 22: 0.93641 0.80418 0.95285 0.76876 0.82588 27: 1.09394 1.00195 1.14976 0.80008 1.11947 1.09484 0.81494 0.68696 35: 0.82364 0.84390 0.71402 0.80293 1.02873 40: Read 39 items plot(density(x)) library(MASS) fitdistr(x, 'lognormal') meanlogsdlog -0.134806360.19118861 ( 0.03061468) ( 0.02164785) lines(dlnorm(x, -.1348, .1911), col='red') Hi Jim I agree with your solution but my plot result not fine. I obtain same result. fitdistr(x, 'lognormal') meanlogsdlog -0.134806360.19118861 ( 0.03061468) ( 0.02164785) In plot when I use points (blue) and curve (green) the fit o lognormal and density(data) is fine but when I use lines (red) the fit is bad (in attach I send a PDF of my output) Do you know why this happen? Thanks in advance -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil P.S. my script is x - scan() 0.80010 0.72299 0.69893 0.99597 0.89200 0.69312 0.73613 1.13559 0.85009 0.85804 0.73324 1.04826 0.84002 1.76330 0.71980 0.89416 0.89450 0.98670 0.83571 0.73833 0.66549 0.93641 0.80418 0.95285 0.76876 0.82588 1.09394 1.00195 1.14976 0.80008 1.11947 1.09484 0.81494 0.68696 0.82364 0.84390 0.71402 0.80293 1.02873 require(MASS) fitdistr(x, 'lognormal') pdf(adj.pdf) plot(density(x)) lines(dlnorm(x, -.1348, .1911), col='red') points(x,dlnorm(x, -.1348, .1911), col='blue') curve(dlnorm(x, -.1348, .1911), col='green',add=T) dev.off() adj.pdf Description: Adobe PDF document __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] showing values in levelplot or image
I am trying to visualize a 2D matrix, with some auxiliary labels associated with each rectangle in the chart. The image command and levelplot in the lattice package map data to colors, but I couldn't find any option to specify values I want to show. Is there an easy way to do this? Thanks, Ken [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] predicting values into the future
Try this: library(forecast) weight - as.numeric(weight) fc - forecast(weight) plot(fc) With this data it chooses a model with multiplicative error and trend and no seasonality. On Sun, Apr 5, 2009 at 2:13 AM, Felipe Carrillo mazatlanmex...@yahoo.com wrote: Hi: I have usually used the GROWTH() excel function to do this but now want to see if I can do this with R. I want to predict values into the future, possibly with the predict.arima Function. I have the following weekly fish weight averages: weight - c(2.1,2.4,2.8,3.6,4.1,5.2,6.3) week - c(1,2,3,4,5,6,7) I would like to predict what the weight will be by week 10 based on my weight values and make a line plot of all the weights(including the predicted values). I have two questions: 1- Should the predicted values be linear or exponential? 2- Is the predict.arima function appropriate to do this? Thanks in advance. Felipe D. Carrillo Supervisory Fishery Biologist Department of the Interior US Fish Wildlife Service California, USA __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Filling in Gapped Interval of a Data Frame As Series
Try this: library(zoo) Lines - -8 100 + -6 20.2 + -1 1.5 +2 3.00 +3 78.8 +5 33.2 +6 44.5 +7 5.00 DF - read.table(textConnection(Lines)) z - zoo(DF$V2, DF$V1) zz - merge(z, zoo(, seq(min(time(z)), max(time(z, fill = 0); zz -8-7-6-5-4-3-2-1 0 1 2 3 4 100.0 0.0 20.2 0.0 0.0 0.0 0.0 1.5 0.0 0.0 3.0 78.8 0.0 5 6 7 33.2 44.5 5.0 If you want to create a data frame from that try this: newDF - data.frame(Time = time(zz), Value = coredata(zz)) although you might be best off to just leave it as a zoo object so that you can make use of other zoo methods too. plot(zz, type = o) See the three zoo vignettes (PDF documents) for more info. On Sun, Apr 5, 2009 at 3:22 AM, Gundala Viswanath gunda...@gmail.com wrote: Dear all, I have a data frame that looks like this: dat V1 V2 1 -8 100 2 -6 20.2 3 -1 1.5 4 2 3.00 5 3 78.8 6 5 33.2 7 6 44.5 8 7 5.00 Now I want to fill the V1 column in series (note that it is gapped), with the corresponding value in V2 as 0. So in the end we would like to have: newdat V1 V2 1 -8 100 2 -6 20.2 3 -7 0 4 -6 0 5 -5 0 6 -4 0 7 -3 0 8 -2 0 9 -1 1.5 10 0 0 11 1 0 12 2 3.00 13 3 78.8 14 4 0 15 5 33.2 16 6 44.5 17 7 5.00 Is there a way to do it in R? - Gundala Viswanath Jakarta - Indonesia __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] unicode only works with a second one
Hi Greg and Paul, I tried several things, but I did not succeed: * I could not find the library(EBImage) on CRAN in Austria to open an png image in R. * I could not import the image via pixmap (read.pnm) as described on http://wiki.r-project.org/rwiki/doku.php?id=tips:graphics-misc:translucency because my GIMP cannot write pnm format. * I could not manage use the 'grImport' package to trace the svg: readPicture(Aries.svg) Fehler in readPicture(Aries.svg) : Version mismatch: RGML file needs to be recreated with PostScriptTrace() * I gave up modifiying the svg code from wikipedia to make it an R array (structure) as you greg described it above. If you have any hint for me please let me know. I am willing to contribute something to TeachingDemos (although I am not sure if this is not a license problem as I trace the (public domain) images from wikimedia. Otherwise I am happy with the Hershey fonts so far. Thomas __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Eclipse and StatET Howto (also added Subversion, Rtools)
Peter Dalgaard wrote: Hum.. Fedora 9 doesn't have it either. It does have /usr/share/X11/rgb.txt though, so please check whether it has moved. I'm curious as to why (only) R/tcltk would be confused by this sort of thing. If you check here: https://bugs.launchpad.net/ubuntu/+source/xorg/+bug/300935 You will find many other unhappy apps: emacs xterm vnc freeNX netpbm x3270 stage and any others that uses names in place of rgb codes for X11. Like one of the posters said, I don't care for arguing whether or not the file should be made obsolete or not. All I know is, removing it breaks some very basic programs. - Ken -- View this message in context: http://www.nabble.com/Eclipse-and-StatET-Howto-%28also-added-Subversion%2C-Rtools%29-tp22764049p22892507.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] help with formula and data= argument
Hi, It seems to me that you could write a generic function myplot() and have different methods for each class of object (like plot does). Either S3 or S4 classes would do I think. Then it is only a matter of making each method work separately. In particular, the method for a formula could contain the with(data, ) construct to make all the variables accessible to the text() function. HTH, baptiste On 4 Apr 2009, at 23:27, Derek Ogle wrote: Sorry for posting this twice, but I still have not solved this problem and am hoping for some assistance. I am attempting to write a function that is flexible enough to respond to the user providing a formula (with a data= argument) or not (similar to plot(x,y) versus plot(y~x,data=data)). I have found a method to work with this in a simple case but am having trouble determining how to find a variable from within the data= argument that is not part of the formula. The following code illustrates (I know that plotrix::thigmophobe.labels() does what this function does) my problem ... myplot - function(x,y=NULL,data=NULL,label=NULL) { if (class(x)==formula) { mf - model.frame(x,data=data) x - mf[,2] y - mf[,1] } if (is.null(y)) stop(Y-axis variable is missing) plot(x,y) if (!is.null(label)) text(x,y,label) } # dummy data df - data .frame(x=runif(10),y=runif(10),grp=factor(rep(c(Yes,No),each=5)) ) # both calls work as expected with(df,myplot(x,y)) myplot(y~x,data=df) # only first works as I would hope with(df,myplot(x,y,label=grp)) myplot(y~x,data=df,label=grp) # this works but is clumsy myplot(y~x,data=df,label=df$grp) Any help with how to make this function recognize the grp variable in df when using the formula without having to type df$grp when supplying it to the label= argument would be greatly appreciated. Thank you in advance. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. _ Baptiste Auguié School of Physics University of Exeter Stocker Road, Exeter, Devon, EX4 4QL, UK Phone: +44 1392 264187 http://newton.ex.ac.uk/research/emag __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Time series forecasting
Dear all: I'm a newbie and an amateur seeking help with forecasting the next in a non-stationary time series, with constraints of 1 (low) and 27 (high) applicable to all. What I need help with is the solution concept. The series has 439 observations as of last week. I'd like to analyze obs 1 - 30 (which are historical and therefore invariate), to solve for 31. The history: Obs 12 21 31 416 59 66 77 811 911 101 1112 1214 1313 142 154 165 1714 186 194 207 215 228 237 2415 2511 263 274 286 298 304 31?? (a known) For backtesting of forecasting accuracy, I can use either a sliding window ( 1 - 30 to solve for 31, 2 - 31 to solve for 32, 3 - 32 to solve for 33, etc.) OR a cumulative window (1 - 30 to solve for 31, 1 - 31 to solve for 32, 1 - 32 to solve for 33, etc.), whichever works better. I can also supply different windows if deemed appropriate, e.g., 50 or 75 or 100 obs, whatever, in either configuration. The 30 obs window is selected for this list query only so as not to take up too much message space. Query: How would you solve for ob 31 in the above series, with the constraints stated? (If you need a longer history, say, 50 obs or more, I can supply it off-list.) I've tried all the relevant Excel functions with no success, and suspect that the solution lies in some combination of them. Here I defer to the collective wisdom of you all. Once the correct concept is established, I can proceed to set it up in R for this and other similar series. Many TIA and regards, Perry E. Gary Tokyo [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] help with formula and data= argument
Try adding this line: label - eval(substitute(label), df, environment(formula)) On Sat, Apr 4, 2009 at 6:27 PM, Derek Ogle do...@northland.edu wrote: Sorry for posting this twice, but I still have not solved this problem and am hoping for some assistance. I am attempting to write a function that is flexible enough to respond to the user providing a formula (with a data= argument) or not (similar to plot(x,y) versus plot(y~x,data=data)). I have found a method to work with this in a simple case but am having trouble determining how to find a variable from within the data= argument that is not part of the formula. The following code illustrates (I know that plotrix::thigmophobe.labels() does what this function does) my problem ... myplot - function(x,y=NULL,data=NULL,label=NULL) { if (class(x)==formula) { mf - model.frame(x,data=data) x - mf[,2] y - mf[,1] } if (is.null(y)) stop(Y-axis variable is missing) plot(x,y) if (!is.null(label)) text(x,y,label) } # dummy data df - data.frame(x=runif(10),y=runif(10),grp=factor(rep(c(Yes,No),each=5)) ) # both calls work as expected with(df,myplot(x,y)) myplot(y~x,data=df) # only first works as I would hope with(df,myplot(x,y,label=grp)) myplot(y~x,data=df,label=grp) # this works but is clumsy myplot(y~x,data=df,label=df$grp) Any help with how to make this function recognize the grp variable in df when using the formula without having to type df$grp when supplying it to the label= argument would be greatly appreciated. Thank you in advance. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] comparing columns in a dataframe
Hi, Have you looked at the compare package? It might do what you want (I just remember seeing its description on R News recently but I've never used it), d - data.frame(x=1:10,y=sin(1:10),z=factor(letters[1:10])) d1 - d d1$x[2:3] - jitter(d$x[2:3] ) d2 - subset(d1, !(z %in% c(a,g))) compare(d,d2,allowAll=T)-test str(test) Just an idea, baptiste On 4 Apr 2009, at 23:45, markle...@verizon.net wrote: Hi: you've got to create a setdiff in both directions in order to get the lone ones in each column because setdiff is not commutative meaning that setdiff(a,b) does not equal setdiff(b,a). once you do that, then ( setdiff1 + setdiff2 - intersect ) should equal the union. if it doesn't, that would be weird and more investigation would need to be done. On Apr 4, 2009, Bob Green bgr...@dyson.brisnet.org.au wrote: hello, I am hoping for some advice regarding comparing variables from 3 versions of a spreadsheet which have been combined into a single dataframe. The aim is to identify which rows have been changed. The dataframe contains 177 rows of data (each cell contains text). 'intersect' produced a file with 35 rows, 'union' a file with 303 rows and 'setdiff' a file with 130 rows Below is the code that I have started with. Ideally I would like to identify the actual row numbers where there is difference in the variables (either pairwise or between 3 variables). x - read.csv(c://rec_compare.csv,header=T, as.is=TRUE) u - union(x$rm1, x$redc1) write.csv(u,c:/union_test.csv) i - intersect(x$rm1, x$redc1) write.csv(i,c:/intersect_test.csv) sd - setdiff(x$rm1, x$redc1) write.csv(sd,c:/setdiff_test.csv) Any suggestions are appreciated. regards Bob __ [1]r-h...@r-project.org mailing list [2]https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide [3]http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. References 1. mailto:R-help@r-project.org 2. https://stat.ethz.ch/mailman/listinfo/r-help 3. http://www.R-project.org/posting-guide.html __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. _ Baptiste Auguié School of Physics University of Exeter Stocker Road, Exeter, Devon, EX4 4QL, UK Phone: +44 1392 264187 http://newton.ex.ac.uk/research/emag __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] extract the p value of F statistics from the lm class
Dear R users I have run an regression and want to extract the p value of the F statistics, but I can find a way to do that. x-summary(lm(log(RV2)~log(IV.m),data=b)) Call: lm(formula = log(RV2) ~ log(IV.m), data = b[[11]]) Residuals: Min 1Q Median 3Q Max -0.26511 -0.09718 -0.01326 0.11095 0.29777 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) -0.3059 0.1917 -1.5950.121 log(IV.m) 0.9038 0.1065 8.488 1.38e-09 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 0.1435 on 31 degrees of freedom Multiple R-squared: 0.6991, Adjusted R-squared: 0.6894 F-statistic: 72.04 on 1 and 31 DF, p-value: 1.379e-09 names(x) [1] call terms residuals [4] coefficients aliased sigma [7] dfr.squared adj.r.squared [10] fstatisticcov.unscaled x$fstatistic valuenumdfdendf 72.04064 1.0 31.0 But can not find the p value of F statistics. Thanks Ted -- View this message in context: http://www.nabble.com/extract-the-p-value-of-F-statistics-from-the-lm-class-tp22891475p22891475.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] extract the p value of F statistics from the lm class
On 05/04/2009 4:18 AM, tedzzx wrote: Dear R users I have run an regression and want to extract the p value of the F statistics, but I can find a way to do that. x-summary(lm(log(RV2)~log(IV.m),data=b)) Call: lm(formula = log(RV2) ~ log(IV.m), data = b[[11]]) Residuals: Min 1Q Median 3Q Max -0.26511 -0.09718 -0.01326 0.11095 0.29777 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) -0.3059 0.1917 -1.5950.121 log(IV.m) 0.9038 0.1065 8.488 1.38e-09 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 0.1435 on 31 degrees of freedom Multiple R-squared: 0.6991, Adjusted R-squared: 0.6894 F-statistic: 72.04 on 1 and 31 DF, p-value: 1.379e-09 names(x) [1] call terms residuals [4] coefficients aliased sigma [7] dfr.squared adj.r.squared [10] fstatisticcov.unscaled x$fstatistic valuenumdfdendf 72.04064 1.0 31.0 But can not find the p value of F statistics. If you're looking for something like that, the two places to look are: - the man page ?summary.lm (which gives the answer) - unclass(x) will display the object without the fancy printing, so you can see that the man page is accurate. (Sometimes man pages are incomplete, and this way is needed, but not in this case.) Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] extract the p value of F statistics from the lm class
Hi, what about the following: ## some test data x - 1:10 y - x + rnorm(x) ## model and summary m - lm(y~x) sm - summary(m) sm # str(sm) # sm$fstatistic ## and now: the manual case 1 - pf(sm$fstatistic[1], sm$fstatistic[2], sm$fstatistic[3]) Hope it helps, ThPe tedzzx schrieb: Dear R users I have run an regression and want to extract the p value of the F statistics, but I can find a way to do that. x-summary(lm(log(RV2)~log(IV.m),data=b)) Call: lm(formula = log(RV2) ~ log(IV.m), data = b[[11]]) Residuals: Min 1Q Median 3Q Max -0.26511 -0.09718 -0.01326 0.11095 0.29777 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) -0.3059 0.1917 -1.5950.121 log(IV.m) 0.9038 0.1065 8.488 1.38e-09 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 0.1435 on 31 degrees of freedom Multiple R-squared: 0.6991, Adjusted R-squared: 0.6894 F-statistic: 72.04 on 1 and 31 DF, p-value: 1.379e-09 names(x) [1] call terms residuals [4] coefficients aliased sigma [7] dfr.squared adj.r.squared [10] fstatisticcov.unscaled x$fstatistic valuenumdfdendf 72.04064 1.0 31.0 But can not find the p value of F statistics. Thanks Ted __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] extract the p value of F statistics from the lm class
On 05/04/2009 4:18 AM, tedzzx wrote: Dear R users I have run an regression and want to extract the p value of the F statistics, but I can find a way to do that. x-summary(lm(log(RV2)~log(IV.m),data=b)) Call: lm(formula = log(RV2) ~ log(IV.m), data = b[[11]]) Residuals: Min 1Q Median 3Q Max -0.26511 -0.09718 -0.01326 0.11095 0.29777 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) -0.3059 0.1917 -1.5950.121 log(IV.m) 0.9038 0.1065 8.488 1.38e-09 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 0.1435 on 31 degrees of freedom Multiple R-squared: 0.6991, Adjusted R-squared: 0.6894 F-statistic: 72.04 on 1 and 31 DF, p-value: 1.379e-09 names(x) [1] call terms residuals [4] coefficients aliased sigma [7] dfr.squared adj.r.squared [10] fstatisticcov.unscaled x$fstatistic valuenumdfdendf 72.04064 1.0 31.0 But can not find the p value of F statistics. Sorry, I misread your question: the p-value of that statistic isn't returned. You just need to calculate it yourself, as f - x$fstatistic pf(f[1], f[2], f[3], lower=FALSE) Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] extract the p value of F statistics from the lm class
On 05-Apr-09 08:18:27, tedzzx wrote: Dear R users I have run an regression and want to extract the p value of the F statistics, but I can find a way to do that. x-summary(lm(log(RV2)~log(IV.m),data=b)) Call: lm(formula = log(RV2) ~ log(IV.m), data = b[[11]]) Residuals: Min 1Q Median 3Q Max -0.26511 -0.09718 -0.01326 0.11095 0.29777 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) -0.3059 0.1917 -1.5950.121 log(IV.m) 0.9038 0.1065 8.488 1.38e-09 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 0.1435 on 31 degrees of freedom Multiple R-squared: 0.6991, Adjusted R-squared: 0.6894 F-statistic: 72.04 on 1 and 31 DF, p-value: 1.379e-09 names(x) [1] call terms residuals [4] coefficients aliased sigma [7] dfr.squared adj.r.squared [10] fstatisticcov.unscaled x$fstatistic valuenumdfdendf 72.04064 1.0 31.0 But can not find the p value of F statistics. Thanks Ted Maybe you were looking in the wrong place. A few lines above the output from x$fstatistic x$fstatistic valuenumdfdendf 72.04064 1.0 31.0 you will find F-statistic: 72.04 on 1 and 31 DF, p-value: 1.379e-09 and therefore will find the P-value. However, maybe that is not the question you really wanted to ask. If that is what I think it may be, you could 1: Observe that x$fstatistic is a vector with 3 values which are: value of F; numerator df; demoninator df 2: Note (from ?pf) pf(q, df1, df2, ncp, lower.tail = TRUE, log.p = FALSE) 3: Therefore do pf(x$fstatistic[1],x$fstatistic[2],x$fstatistic[3],lower.tail=FALSE) # [1] 1.378626e-09 Note that the P-value is not in the list of values returned by lm() although $fstatistic is one of the values. The computation of the P-value in the displayed output from summary.lm() is done by the 'print' method for summary.lm() (just as in [3] above). Hoping this helps, Ted. E-Mail: (Ted Harding) ted.hard...@manchester.ac.uk Fax-to-email: +44 (0)870 094 0861 Date: 05-Apr-09 Time: 13:12:45 -- XFMail -- __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] threshold distribution
You didn't properly specify the x-axis coordinates in your call to lines(): plot(density(x)) lines(x,dlnorm(x, -.1348, .1911), col='red') points(x,dlnorm(x, -.1348, .1911), col='blue') curve(dlnorm(x, -.1348, .1911), col='green',add=T) On Sun, Apr 5, 2009 at 6:40 AM, Bernardo Rangel Tura t...@centroin.com.br wrote: On Sun, 2009-04-05 at 01:13 -0400, jim holtman wrote: Here is what I get from using 'fitdistr' in R to fit to a lognormal. The resulting density plot from the distribution seems to be a reason match to the data. x - scan() 1: 0.80010 0.72299 0.69893 0.99597 0.89200 0.69312 0.73613 1.13559 9: 0.85009 0.85804 0.73324 1.04826 0.84002 14: 1.76330 0.71980 0.89416 0.89450 0.98670 0.83571 0.73833 0.66549 22: 0.93641 0.80418 0.95285 0.76876 0.82588 27: 1.09394 1.00195 1.14976 0.80008 1.11947 1.09484 0.81494 0.68696 35: 0.82364 0.84390 0.71402 0.80293 1.02873 40: Read 39 items plot(density(x)) library(MASS) fitdistr(x, 'lognormal') meanlog sdlog -0.13480636 0.19118861 ( 0.03061468) ( 0.02164785) lines(dlnorm(x, -.1348, .1911), col='red') Hi Jim I agree with your solution but my plot result not fine. I obtain same result. fitdistr(x, 'lognormal') meanlog sdlog -0.13480636 0.19118861 ( 0.03061468) ( 0.02164785) In plot when I use points (blue) and curve (green) the fit o lognormal and density(data) is fine but when I use lines (red) the fit is bad (in attach I send a PDF of my output) Do you know why this happen? Thanks in advance -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil P.S. my script is x - scan() 0.80010 0.72299 0.69893 0.99597 0.89200 0.69312 0.73613 1.13559 0.85009 0.85804 0.73324 1.04826 0.84002 1.76330 0.71980 0.89416 0.89450 0.98670 0.83571 0.73833 0.66549 0.93641 0.80418 0.95285 0.76876 0.82588 1.09394 1.00195 1.14976 0.80008 1.11947 1.09484 0.81494 0.68696 0.82364 0.84390 0.71402 0.80293 1.02873 require(MASS) fitdistr(x, 'lognormal') pdf(adj.pdf) plot(density(x)) lines(dlnorm(x, -.1348, .1911), col='red') points(x,dlnorm(x, -.1348, .1911), col='blue') curve(dlnorm(x, -.1348, .1911), col='green',add=T) dev.off() __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Mike Lawrence Graduate Student Department of Psychology Dalhousie University Looking to arrange a meeting? Check my public calendar: http://tinyurl.com/mikes-public-calendar ~ Certainty is folly... I think. ~ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Quick Question - MLE for Geometric Brownian Motion Process with Jumps
Hi All, I am relatively new to R and having a great experience with it but have come up with a little road block. I am tying to run a maxlik regression and keep getting the error, NA in the initial gradient My Code is below: gbmploglik-function(param){ mu-param[1] sigma-param[2] lamda-param[3] nu-param[4] gama-param[5] logLikVal- - n*lamda - .5*n*log(2*pi) + sum(log(sum(for(j in 1:10)(cat((lamda^j/factorial(j))*(1/((sigma^2+j*gama^2)^.5)*exp( - (combinedlrph1-mu-j*nu)^2/2*(sigma^2+j*gama^2 logLikVal } rescbj- maxLik(gbmploglik, grad = NULL, hess = NULL, start=c(1,1,1,1,1), method = Newton-Raphson) I was wondering if there is something that I have to do with the grad= and maybe put something other then NULL. The other issue is that there might be something wrong with the loglikelihood function, because of the loop that I put in it. If you have any suggestion or see something incorrect that I am doing please let me know. Thanks JP __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] problem with lattice tiff or bitmap: character size and color
Hi all, I am trying to make tiff files of lattice plots at a resolution greater than 300 dpi required by a journal (PLoS ONE). I have tried both the tiff and bitmap functions. tiff keeps panel colors but reduces axes and tick labels so they are nearly invisible. bitmap maintains correct label size but only produces greyscale. Regular plots work fine with tiff; the problem is only with lattice plots. How can I keep the right sizes using tiff or, alternatively, how can I get bitmap to produce colors? I have made this test function using an xyplot example to demonstrate the problem. `test` - function(tiff=F, bitmap=F) { ### requires states ### states - data.frame(state.x77, state.name = dimnames(state.x77)[[1]], ### state.region = state.region) if(bitmap) bitmap(bitmaptest.tiff, width=17.15, height=17.15, units=cm, res=1200, pointsize=10, type=tifflzw,bg=white) if(tiff) tiff(file=tifftest.tiff,width=17.15,height=17.15,units=cm, res=1200, pointsize=10, compression = lzw) ### tiff(file=tifftest.tiff,width=17.15,height=17.15,units=cm, ###res=72, pointsize=10, compression = lzw) plot(xyplot(Murder ~ Population | state.region, data = states, groups = state.name, panel = function(x, y, subscripts, groups) ltext(x = x, y = y, label = groups[subscripts], cex=0.5, fontfamily = HersheySans))) if(bitmap) dev.off() if(tiff) dev.off() } R.Version() $platform [1] x86_64-redhat-linux-gnu $arch [1] x86_64 $os [1] linux-gnu $system [1] x86_64, linux-gnu $status [1] $major [1] 2 $minor [1] 8.1 $year [1] 2008 $month [1] 12 $day [1] 22 $`svn rev` [1] 47281 $language [1] R $version.string [1] R version 2.8.1 (2008-12-22) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] unicode only works with a second one
It's in the Bioconductor repository. http://bioconductor.org/packages/2.3/bioc/html/EBImage.html On Apr 5, 2009, at 6:14 AM, Thomas Steiner wrote: EBImage David Winsemius, MD Heritage Laboratories West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Code of sample() in C
Dear All, I would like to use the function sample() in a program written in C. Is there somewhere the code of sample() written in C? Thanks in advance, Paul __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] threshold distribution
The data suggest a lognormal threshold to me (and perhaps to the originator, based on his title). ## x - c(0.80010, 0.72299, 0.69893, 0.99597, 0.89200, 0.69312, 0.73613, 1.13559, 0.85009, 0.85804, 0.73324, 1.04826, 0.84002, 1.76330, 0.71980, 0.89416, 0.89450, 0.98670, 0.83571, 0.73833, 0.66549, 0.93641, 0.80418, 0.95285, 0.76876, 0.82588, 1.09394, 1.00195, 1.14976, 0.80008, 1.11947, 1.09484, 0.81494, 0.68696, 0.82364, 0.84390, 0.71402, 0.80293, 1.02873) # plot the original density x.threshold - 0 qqnorm(log(x - x.threshold), datax = TRUE) qqline(log(x - x.threshold), datax = TRUE) # plot the lognormal density with a threshold x.threshold - 0.63 qqnorm(log(x - x.threshold), datax = TRUE) qqline(log(x - x.threshold), datax = TRUE) # compute and plot the associated log likelihoods x.threshold - 0 loglikelihood - 1/(x-x.threshold) * dnorm(log(x - x.threshold), mean(log(x - x.threshold)), sd(log(x - x.threshold))) x.threshold - 0.63 loglikelihood.63 - 1/(x-x.threshold) * dnorm(log(x - x.threshold), mean(log(x - x.threshold)), sd(log(x - x.threshold))) x.minus.x.threshold - x - x.threshold plot(x.minus.x.threshold, loglikelihood.63, xlim = c(0, 2), ylim =c(0, 3.5)) points(x, loglikelihood, col=red) # compare loglikelihood ratio with chi-square sum.loglikelihood - sum(loglikelihood) print(sum.loglikelihood) sum.loglikelihood.63 - sum(loglikelihood.63) print(sum.loglikelihood.63) log.likelihiood.ratio - sum.loglikelihood.63 - sum.loglikelihood print(log.likelihiood.ratio) significant.difference - qchisq(p = 0.95, df = 1)/2 print(significant.difference) # Charles Annis, P.E. charles.an...@statisticalengineering.com phone: 561-352-9699 eFax: 614-455-3265 http://www.StatisticalEngineering.com -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Mike Lawrence Sent: Sunday, April 05, 2009 8:13 AM To: t...@centroin.com.br Cc: r-help Subject: Re: [R] threshold distribution You didn't properly specify the x-axis coordinates in your call to lines(): plot(density(x)) lines(x,dlnorm(x, -.1348, .1911), col='red') points(x,dlnorm(x, -.1348, .1911), col='blue') curve(dlnorm(x, -.1348, .1911), col='green',add=T) On Sun, Apr 5, 2009 at 6:40 AM, Bernardo Rangel Tura t...@centroin.com.br wrote: On Sun, 2009-04-05 at 01:13 -0400, jim holtman wrote: Here is what I get from using 'fitdistr' in R to fit to a lognormal. The resulting density plot from the distribution seems to be a reason match to the data. x - scan() 1: 0.80010 0.72299 0.69893 0.99597 0.89200 0.69312 0.73613 1.13559 9: 0.85009 0.85804 0.73324 1.04826 0.84002 14: 1.76330 0.71980 0.89416 0.89450 0.98670 0.83571 0.73833 0.66549 22: 0.93641 0.80418 0.95285 0.76876 0.82588 27: 1.09394 1.00195 1.14976 0.80008 1.11947 1.09484 0.81494 0.68696 35: 0.82364 0.84390 0.71402 0.80293 1.02873 40: Read 39 items plot(density(x)) library(MASS) fitdistr(x, 'lognormal') meanlog sdlog -0.13480636 0.19118861 ( 0.03061468) ( 0.02164785) lines(dlnorm(x, -.1348, .1911), col='red') Hi Jim I agree with your solution but my plot result not fine. I obtain same result. fitdistr(x, 'lognormal') meanlog sdlog -0.13480636 0.19118861 ( 0.03061468) ( 0.02164785) In plot when I use points (blue) and curve (green) the fit o lognormal and density(data) is fine but when I use lines (red) the fit is bad (in attach I send a PDF of my output) Do you know why this happen? Thanks in advance -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil P.S. my script is x - scan() 0.80010 0.72299 0.69893 0.99597 0.89200 0.69312 0.73613 1.13559 0.85009 0.85804 0.73324 1.04826 0.84002 1.76330 0.71980 0.89416 0.89450 0.98670 0.83571 0.73833 0.66549 0.93641 0.80418 0.95285 0.76876 0.82588 1.09394 1.00195 1.14976 0.80008 1.11947 1.09484 0.81494 0.68696 0.82364 0.84390 0.71402 0.80293 1.02873 require(MASS) fitdistr(x, 'lognormal') pdf(adj.pdf) plot(density(x)) lines(dlnorm(x, -.1348, .1911), col='red') points(x,dlnorm(x, -.1348, .1911), col='blue') curve(dlnorm(x, -.1348, .1911), col='green',add=T) dev.off() __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Mike Lawrence Graduate Student Department of Psychology Dalhousie University Looking to arrange a meeting? Check my public calendar: http://tinyurl.com/mikes-public-calendar ~ Certainty is folly... I think. ~ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] loop problem for extract coefficients
Dear R users, I have problem with extracting coefficients from a object. Here, X (predictor)and Y (response) are two matrix , I am regressing X ( dimensions 10 x 20) on each of columns of Y[,1] (10 x 1) and want to store the coefficient values. I have performed a Elastic Net regression and I want to store the coeffcients in each iteration. I got an error message . I do not know where is the problem Please help me. Thanks *Code:* --- library(elasticnet) X-matrix(rnorm(200),ncol=20) Y-matrix(rnorm(200),ncol=20) loop - 20 size - 20 enres-matrix(nrow = size, ncol = loop) fit-matrix(nrow = size, ncol = loop) store-matrix(nrow = size, ncol = loop) for(j in 1: 10) print (paste(j,/200,sep=)) { enres-enet(x=X,y=Y[,j],lambda=1,normalize=TRUE,intercept=TRUE) fit-predict.enet(enres, X, type=coefficients) store[,j]-fit$coefficients } library(elasticnet) Loading required package: lars X-matrix(rnorm(200),ncol=20) Y-matrix(rnorm(200),ncol=20) loop - 20 size - 20 enres-matrix(nrow = size, ncol = loop) fit-matrix(nrow = size, ncol = loop) store-matrix(nrow = size, ncol = loop) for(j in 1: 10) + print (paste(j,/200,sep=)) [1] 1/200 [1] 2/200 [1] 3/200 [1] 4/200 [1] 5/200 [1] 6/200 [1] 7/200 [1] 8/200 [1] 9/200 [1] 10/200 { + enres-enet(x=X,y=Y[,j],lambda=1,normalize=TRUE,intercept=TRUE) + fit-predict.enet(enres, X, type=coefficients) + store[,j]-fit$coefficients + } *Error in store[, j] - fit$coefficients : number of items to replace is not a multiple of replacement length * [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Quick Question - MLE for Geometric Brownian Motion Process with Jumps
John-Paul Taylor johnpaul.taylor at ryerson.ca writes: I am tying to run a maxlik regression and keep getting the error, NA in the initial gradient My Code is below: gbmploglik-function(param){ mu-param[1] sigma-param[2] lamda-param[3] nu-param[4] gama-param[5] logLikVal- - n*lamda - .5*n*log(2*pi) + sum(log(sum(for(j in 1:10)(cat((lamda^j/factorial(j))*(1/((sigma^2+j*gama^2)^.5)*exp( - (combinedlrph1-mu-j*nu)^2/2*(sigma^2+j*gama^2 logLikVal } rescbj- maxLik(gbmploglik, grad = NULL, hess = NULL, start=c(1,1,1,1,1), method = Newton-Raphson) I was wondering if there is something that I have to do with the grad= and maybe put something other then NULL. The other issue is that there might be something wrong with the loglikelihood function, because of the loop that I put in it. I have not used that function, but the error message is rather clear in telling you that your start values are not good. Try to plot the function and the gradient at that point. Sometime moving by a very small amount helps. Dieter __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Code of sample() in C
Hi Paul, It is in the main/src/random.c file of the source code. HTH! Best wishes, Ranjan On Sun, 5 Apr 2009 15:35:25 +0100 Paul Smith phh...@gmail.com wrote: Dear All, I would like to use the function sample() in a program written in C. Is there somewhere the code of sample() written in C? Thanks in advance, Paul __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Eclipse and StatET Howto (also added Subversion, Rtools)
Ken-JP wrote: Peter Dalgaard wrote: Hum.. Fedora 9 doesn't have it either. It does have /usr/share/X11/rgb.txt though, so please check whether it has moved. I'm curious as to why (only) R/tcltk would be confused by this sort of thing. If you check here: https://bugs.launchpad.net/ubuntu/+source/xorg/+bug/300935 You will find many other unhappy apps: emacs xterm vnc freeNX netpbm x3270 stage and any others that uses names in place of rgb codes for X11. Like one of the posters said, I don't care for arguing whether or not the file should be made obsolete or not. All I know is, removing it breaks some very basic programs. - Ken OK, now I get it. Apparently, some people believe that because yellow does not internationalize, everyone should use #00 which is equally cryptic in all languages Still not sure what the exact Tk angle is. Looks like tk8.5 has internalized the color table. Other applications have just copied rgb.txt to some place where they can find them. -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - (p.dalga...@biostat.ku.dk) FAX: (+45) 35327907 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] number of zeros in a matrix -row by row
Hi. I have an n x m matrix M some of who's entries are zeros. I want to know how many zeros there are in each row-perhaps stored in a 1 x n vector which lists the number of zeros for each row of M. Before I had a vector V and I was able to get the number of zeros in V by doing length(V[ V==0]) but when I try something similar for M, like M[ M==0] it creates a vector not a matrix and so this does not work. Does anyone have a solution to this? Thank you. -- View this message in context: http://www.nabble.com/number-of-zeros-in-a-matrix--row-by-row-tp22893147p22893147.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] mtext in barplot
Hii, Can anybody help me to put a text under the barplots. I will describe the percental between six grouped barplots. I tried to do it with mtext but without success. Here is my code: test -read.table(file=D:/mobile.txt) pdf(file = D:/mobil126.pdf, width = 6.67, height = 5, onefile = TRUE, family = Helvetica, title = R Graphics Output, fonts = NULL, version = 1.1, bg=white, pointsize=10) barplot(as.matrix(test), main=OLSR,xlab=Hops, col=c(skyblue1,salmon),width- c(1,1),names=c(2-1,4-2,3-1,2-1,2-1,1-2), legend = rownames(x), beside=TRUE) legend(topright, c(OLSR),cex=0.8,ncol =1.5,col = c(red),bg=c(lightskyblue1)) dev.off() I tired to do it with mtext(side=1,at=x, text =c(Mean, rere), col = red), line = 1, cex = 0.75) but without success, I get this example from a R tutorial... Can anybody help me please ? greetings, johnh -- View this message in context: http://www.nabble.com/mtext-in-barplot-tp22893563p22893563.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] number of zeros in a matrix -row by row
if 'M' is your matrix, then try this: rowSums(M == 0) I hope it helps. Best, Dimitris onyourmark wrote: Hi. I have an n x m matrix M some of who's entries are zeros. I want to know how many zeros there are in each row-perhaps stored in a 1 x n vector which lists the number of zeros for each row of M. Before I had a vector V and I was able to get the number of zeros in V by doing length(V[ V==0]) but when I try something similar for M, like M[ M==0] it creates a vector not a matrix and so this does not work. Does anyone have a solution to this? Thank you. -- Dimitris Rizopoulos Assistant Professor Department of Biostatistics Erasmus University Medical Center Address: PO Box 2040, 3000 CA Rotterdam, the Netherlands Tel: +31/(0)10/7043478 Fax: +31/(0)10/7043014 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] mtext in barplot
johnhj wrote: Hii, Can anybody help me to put a text under the barplots. I will describe the percental between six grouped barplots. I tried to do it with mtext but without success. Here is my code: test -read.table(file=D:/mobile.txt) pdf(file = D:/mobil126.pdf, width = 6.67, height = 5, onefile = TRUE, family = Helvetica, title = R Graphics Output, fonts = NULL, version = 1.1, bg=white, pointsize=10) barplot(as.matrix(test), main=OLSR,xlab=Hops, col=c(skyblue1,salmon),width- c(1,1),names=c(2-1,4-2,3-1,2-1,2-1,1-2), legend = rownames(x), beside=TRUE) legend(topright, c(OLSR),cex=0.8,ncol =1.5,col = c(red),bg=c(lightskyblue1)) dev.off() I tired to do it with mtext(side=1,at=x, text =c(Mean, rere), col = red), line = 1, cex = 0.75) but without success, I get this example from a R tutorial... What is x? The bars in a barplot are at locations 1:n. Uwe Ligges Can anybody help me please ? greetings, johnh __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Eclipse and StatET Howto (also added Subversion, Rtools)
On 5 April 2009 at 03:33, Ken-JP wrote: | If you check here: | | https://bugs.launchpad.net/ubuntu/+source/xorg/+bug/300935 | | You will find many other unhappy apps: | emacs | xterm | vnc | freeNX | netpbm | x3270 | stage | | and any others that uses names in place of rgb codes for X11. | | Like one of the posters said, I don't care for arguing whether or not the | file should be made obsolete or not. All I know is, removing it breaks some | very basic programs. As I read, this was a bug in new Ubuntu installations, and it seems to have been corrected. Dirk -- Three out of two people have difficulties with fractions. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Time series forecasting
Hi Perry, my impression after a very cursory glance: this looks like noise. Perhaps you should think a little more about your series - what kind of seasonality there could be (is this weekly data? or monthly?), whether the peaks and troughs could be due to some kind of external driver, whether you really have count data, that kind of thing. Until then, there is little else to do than to use a very simple method, e.g. forecast the last observation (random walk) or the mean of the observations (white noise), or the median. All of these benchmarks can be surprisingly hard to beat. If you have seasonality but no external influence, you could look at smoothing methods, they are nice to interpret and usually perform very well. I'd recommend Hyndman et al., Forecasting with Exponential Smoothing: The State Space Approach and the accompanying forecast R package, mainly with the ets() function. You could also look at arima(). I fitted an ARIMA model to your data, and as expected, it returned a simple mean (not that I would recommend blindly fitting ARIMA to just any data): Call: arima(x = foo[, 2]) Coefficients: intercept 7.2333 s.e. 0.8009 sigma^2 estimated as 19.25: log likelihood = -86.93, aic = 177.85 And for count data, you could use some variants of ARIMA, e.g., INAR. HTH, Stephan pg...@gol.com schrieb: Dear all: I'm a newbie and an amateur seeking help with forecasting the next in a non-stationary time series, with constraints of 1 (low) and 27 (high) applicable to all. What I need help with is the solution concept. The series has 439 observations as of last week. I'd like to analyze obs 1 - 30 (which are historical and therefore invariate), to solve for 31. The history: Obs 12 21 31 416 59 66 77 811 911 101 1112 1214 1313 142 154 165 1714 186 194 207 215 228 237 2415 2511 263 274 286 298 304 31?? (a known) For backtesting of forecasting accuracy, I can use either a sliding window ( 1 - 30 to solve for 31, 2 - 31 to solve for 32, 3 - 32 to solve for 33, etc.) OR a cumulative window (1 - 30 to solve for 31, 1 - 31 to solve for 32, 1 - 32 to solve for 33, etc.), whichever works better. I can also supply different windows if deemed appropriate, e.g., 50 or 75 or 100 obs, whatever, in either configuration. The 30 obs window is selected for this list query only so as not to take up too much message space. Query: How would you solve for ob 31 in the above series, with the constraints stated? (If you need a longer history, say, 50 obs or more, I can supply it off-list.) I've tried all the relevant Excel functions with no success, and suspect that the solution lies in some combination of them. Here I defer to the collective wisdom of you all. Once the correct concept is established, I can proceed to set it up in R for this and other similar series. Many TIA and regards, Perry E. Gary Tokyo [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] threshold distribution
In my haste I forgot to take the logs of the likelihoods. How embarrassing. The conclusion is unchanged - the data support a non-zero threshold. Here's the corrected code: x.threshold - 0 loglikelihood - log(1/(x-x.threshold) * dnorm(log(x - x.threshold), mean(log(x - x.threshold)), sd(log(x - x.threshold sum.loglikelihood - sum(loglikelihood) print(sum.loglikelihood) x.threshold - 0.63 loglikelihood.63 - log(1/(x-x.threshold) * dnorm(log(x - x.threshold), mean(log(x - x.threshold)), sd(log(x - x.threshold sum.loglikelihood.63 - sum(loglikelihood.63) print(sum.loglikelihood.63) x.minus.x.threshold - x - x.threshold plot(loglikelihood.63 ~ x.minus.x.threshold, xlim = c(0, 2), ylim =c(-10, 2)) points(x, loglikelihood, col=red) log.likelihiood.ratio - sum.loglikelihood.63 - sum.loglikelihood print(log.likelihiood.ratio) significant.difference - qchisq(p = 0.95, df = 1)/2 print(significant.difference) Charles Annis, P.E. charles.an...@statisticalengineering.com phone: 561-352-9699 eFax: 614-455-3265 http://www.StatisticalEngineering.com -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Charles Annis, P.E. Sent: Sunday, April 05, 2009 10:39 AM To: 'Mike Lawrence'; t...@centroin.com.br Cc: 'r-help' Subject: Re: [R] threshold distribution The data suggest a lognormal threshold to me (and perhaps to the originator, based on his title). ## x - c(0.80010, 0.72299, 0.69893, 0.99597, 0.89200, 0.69312, 0.73613, 1.13559, 0.85009, 0.85804, 0.73324, 1.04826, 0.84002, 1.76330, 0.71980, 0.89416, 0.89450, 0.98670, 0.83571, 0.73833, 0.66549, 0.93641, 0.80418, 0.95285, 0.76876, 0.82588, 1.09394, 1.00195, 1.14976, 0.80008, 1.11947, 1.09484, 0.81494, 0.68696, 0.82364, 0.84390, 0.71402, 0.80293, 1.02873) # plot the original density x.threshold - 0 qqnorm(log(x - x.threshold), datax = TRUE) qqline(log(x - x.threshold), datax = TRUE) # plot the lognormal density with a threshold x.threshold - 0.63 qqnorm(log(x - x.threshold), datax = TRUE) qqline(log(x - x.threshold), datax = TRUE) # compute and plot the associated log likelihoods x.threshold - 0 loglikelihood - 1/(x-x.threshold) * dnorm(log(x - x.threshold), mean(log(x - x.threshold)), sd(log(x - x.threshold))) x.threshold - 0.63 loglikelihood.63 - 1/(x-x.threshold) * dnorm(log(x - x.threshold), mean(log(x - x.threshold)), sd(log(x - x.threshold))) x.minus.x.threshold - x - x.threshold plot(x.minus.x.threshold, loglikelihood.63, xlim = c(0, 2), ylim =c(0, 3.5)) points(x, loglikelihood, col=red) # compare loglikelihood ratio with chi-square sum.loglikelihood - sum(loglikelihood) print(sum.loglikelihood) sum.loglikelihood.63 - sum(loglikelihood.63) print(sum.loglikelihood.63) log.likelihiood.ratio - sum.loglikelihood.63 - sum.loglikelihood print(log.likelihiood.ratio) significant.difference - qchisq(p = 0.95, df = 1)/2 print(significant.difference) # Charles Annis, P.E. charles.an...@statisticalengineering.com phone: 561-352-9699 eFax: 614-455-3265 http://www.StatisticalEngineering.com -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Mike Lawrence Sent: Sunday, April 05, 2009 8:13 AM To: t...@centroin.com.br Cc: r-help Subject: Re: [R] threshold distribution You didn't properly specify the x-axis coordinates in your call to lines(): plot(density(x)) lines(x,dlnorm(x, -.1348, .1911), col='red') points(x,dlnorm(x, -.1348, .1911), col='blue') curve(dlnorm(x, -.1348, .1911), col='green',add=T) On Sun, Apr 5, 2009 at 6:40 AM, Bernardo Rangel Tura t...@centroin.com.br wrote: On Sun, 2009-04-05 at 01:13 -0400, jim holtman wrote: Here is what I get from using 'fitdistr' in R to fit to a lognormal. The resulting density plot from the distribution seems to be a reason match to the data. x - scan() 1: 0.80010 0.72299 0.69893 0.99597 0.89200 0.69312 0.73613 1.13559 9: 0.85009 0.85804 0.73324 1.04826 0.84002 14: 1.76330 0.71980 0.89416 0.89450 0.98670 0.83571 0.73833 0.66549 22: 0.93641 0.80418 0.95285 0.76876 0.82588 27: 1.09394 1.00195 1.14976 0.80008 1.11947 1.09484 0.81494 0.68696 35: 0.82364 0.84390 0.71402 0.80293 1.02873 40: Read 39 items plot(density(x)) library(MASS) fitdistr(x, 'lognormal') meanlog sdlog -0.13480636 0.19118861 ( 0.03061468) ( 0.02164785) lines(dlnorm(x, -.1348, .1911), col='red') Hi Jim I agree with your solution but my plot result not fine. I obtain same result. fitdistr(x, 'lognormal') meanlog sdlog -0.13480636 0.19118861 ( 0.03061468) ( 0.02164785) In plot when I use points (blue) and curve (green) the fit o lognormal and density(data) is fine but when I use lines (red) the fit is bad (in attach I send a PDF of my output) Do you know why this happen? Thanks in advance -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil P.S. my script is x - scan()
[R] Problem with Dynamo-Package
Good day, I am facing a problem when I am installing the dynamo-package and loading it. After I installed the package, I received the following warning message: In file.create(f.tg) : cannot create file 'C:\PROGRA~2\R\R-28~1.1/doc/html/packages.html', reason 'Permission denied' and when I load the package, an error message pops up saying that the application failed to start because libgsl.dll was not found. Re-installing the application may fix the problem Can anyone help to solve this problem. Mohammad [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] showing values in levelplot or image
I have used ltext in a panel function like in the code segment below: print(levelplot(reslot ~ aisle + quad, reslot.df, col.regions=topo.colors(max(quad.POG$values)), panel=function(x, y, z, ...){ panel.levelplot(x, y, z, ...) panel.abline(v=seq(1.5, by=1, length=length(aisleOrder)), h=seq(1.5, by=1, length=26), col='gray') ltext(x, y, z, col='red', font=1) }, On Sun, Apr 5, 2009 at 5:42 AM, Ken Wilson kwil...@gmail.com wrote: I am trying to visualize a 2D matrix, with some auxiliary labels associated with each rectangle in the chart. The image command and levelplot in the lattice package map data to colors, but I couldn't find any option to specify values I want to show. Is there an easy way to do this? Thanks, Ken [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Jim Holtman Cincinnati, OH +1 513 646 9390 What is the problem that you are trying to solve? __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Question: R software
I need the commond in R for double integral over the range - infinityxy infinity . Could you please help me. Thanks and regards Salya [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Question: R software
check package adapt, http://cran.r-project.org/web/packages/adapt/index.html Best, Dimitris salia a wrote: I need the commond in R for double integral over the range - infinityxy infinity . Could you please help me. Thanks and regards Salya [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Dimitris Rizopoulos Assistant Professor Department of Biostatistics Erasmus University Medical Center Address: PO Box 2040, 3000 CA Rotterdam, the Netherlands Tel: +31/(0)10/7043478 Fax: +31/(0)10/7043014 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Code of sample() in C
Thanks, Ranjan. Got it! I am now wondering whether there is some simpler way of implementing the following in C: sample(1:50,5) Paul On Sun, Apr 5, 2009 at 4:10 PM, Ranjan Maitra mai...@iastate.edu wrote: Hi Paul, It is in the main/src/random.c file of the source code. HTH! Best wishes, Ranjan On Sun, 5 Apr 2009 15:35:25 +0100 Paul Smith phh...@gmail.com wrote: Dear All, I would like to use the function sample() in a program written in C. Is there somewhere the code of sample() written in C? Thanks in advance, Paul __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] lattice different colors in different areas
In making grid graphs, one can make the background semi-transparent with a line like grid.rect(gp=gpar(lty=0, fill=rgb(.5, .5, 0,. 25))) and then make the area where points and lines are plotted white with lines like pushViewport(plotViewport(c(5,4,3,1))) pushViewport(dataViewport(year, m, name=plotRegion)) grid.rect(gp=gpar(fill=white)) (The area where the labels, title, legend, etc. are remains semi-transparent.) I need to use the features of lattice for some graphs but want to keep the same color theme. In lattice, the background can be changed with the lines below, but it makes the entire graph this color. bsettings=trellis.par.get(background) bsettings$col=rgb(.5, .5, 0, .25) trellis.par.set(background, bsettings) How do I now change just the plot region back to white? Thanks. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] inverting a table
Is there an easy way to invert a table? (not to solve for the inverted matrix, just swap rows for columns vice versa). I've gone through my data manipulation bible (Phil Spector's book), but to no avail. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] [OT ?] rant (was : Re: Conversions From standard to metric units)
On Fri, 03-Apr-2009 at 07:48PM -0400, Duncan Murdoch wrote: On 03/04/2009 5:37 PM, Emmanuel Charpentier wrote: Le vendredi 03 avril 2009 à 14:17 -0400, stephen sefick a écrit : I am starting to use R for almost any sort of calculation that I need. I am a biologist that works in the states, and there is often a need to convert from standard units to metric units. rant US/Imperial units are *not* standard units. But they are fun: you should see the arguments you can have about whether imperial fluid ounces are the same volume as US fluid ounces. (They're not: US ounces are bigger. But not big enough so that their gallons catch up!) Hence the word standard. :-) Duncan Murdoch The former metric system is now called Système International (International System) for a reason, which is *not* gallocentrism of a few 6e7 frogs, but rather laziness of about 5.6e9 losers who refuse to load their memories with meaningless conversion factors... /rant Emmanuel Charpentier who has served his time with pounds per cubic feet, furlongs per fortnight, BTU and other figments of British/American sadistic imagination, thank you very much... /rant # Again, didn't work the first time... Is there a package in R for this already? If not I believe that I am going to write some of the most often used in function form. My question is should I include this in my StreamMetabolism package. It is not along the same theme lines, but could loosely fit. The reason that I ask is that I don't want to clutter CRAN with a small package containing some conversion functions because I am to lazy to source them into R every time that I use them, but I also don't want the StreamMetabolism package to turn into StephenMisc Fuctions. Thoughts, comments, or suggestions would be appreciated. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- ~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~. ___Patrick Connolly {~._.~} Great minds discuss ideas _( Y )_ Average minds discuss events (:_~*~_:) Small minds discuss people (_)-(_) . Eleanor Roosevelt ~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] inverting a table
?t # I knew it worked with matrices so seemed reasonable to see if it worked on tables. t(with(warpbreaks, table(wool, tension)) ) wool tension A B L 9 9 M 9 9 H 9 9 with(warpbreaks, table(wool, tension)) tension wool L M H A 9 9 9 B 9 9 9 On Apr 5, 2009, at 3:56 PM, Donald Braman wrote: Is there an easy way to invert a table? (not to solve for the inverted matrix, just swap rows for columns vice versa). I've gone through my data manipulation bible (Phil Spector's book), but to no avail. [[alternative HTML version deleted]] David Winsemius, MD Heritage Laboratories West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] edit for Time series data in data frames
Dear R People: I can edit a data frame with time series as x.df - edit(as.matrix(x.df)) Is there a better way, please? Thanks, erin -- Erin Hodgess Associate Professor Department of Computer and Mathematical Sciences University of Houston - Downtown mailto: erinm.hodg...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Code of sample() in C
I presume you mean sampling without replacement? The following (adapted from the file you asked for) will do it, but you will need to incorporate Rmath as standalone in order to get runif to work. This function will give you k indices from n, sampled WOR. Thus, n = 50 for you and k = 5. You will get a vector y of length 5 (a pointer of int actually) which will contain these indices. Thus you will define a vector z (of length 50) which is 1:50, and then using the function, your SRWOR sample will be z[y[i]] where i = 0, 1,...4. I haven't tried this function out much myself, so YMMV. HTH! Best wishes, Ranjan #include stdlib.h #ifndef USING_RLIB #define MATHLIB_STANDALONE 1 /*It is essential to have this before the call to the Rmath's header file because this decides the definitions to be set. */ #endif #include Rmath.h /* Note that use of this function involves a prior call to the Rmath library to get the seeds in place. It is assumed that this is done in the calling function. */ /* Equal probability sampling; without-replacement case */ /* Adapted from the R function called SampleNoReplace */ /** * Stores k out of n indices sampled at random without replacement * in y. */ int srswor(int n, int k, int *y) { if (k n) { return 1; } else { const double len = (double) n; int i; int* x = malloc(n * sizeof(int)); if (!x) { return 2; } for (i = 0; i n; ++i) x[i] = i; for (i = 0; i k; ++i) { const int j = (int)(len * runif(0.0, 1.0)); y[i] = x[j]; x[j] = x[--n]; } free(x); } return 0; } On Sun, 5 Apr 2009 20:11:04 +0100 Paul Smith phh...@gmail.com wrote: Thanks, Ranjan. Got it! I am now wondering whether there is some simpler way of implementing the following in C: sample(1:50,5) Paul On Sun, Apr 5, 2009 at 4:10 PM, Ranjan Maitra mai...@iastate.edu wrote: Hi Paul, It is in the main/src/random.c file of the source code. HTH! Best wishes, Ranjan On Sun, 5 Apr 2009 15:35:25 +0100 Paul Smith phh...@gmail.com wrote: Dear All, I would like to use the function sample() in a program written in C. Is there somewhere the code of sample() written in C? Thanks in advance, Paul __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Souce macros help
I was trying to understand some of the source in optimi.c and in the SANN source I see: SETCADR(OS-R_gcall, x); PROTECT_WITH_INDEX(s = eval(OS-R_gcall, OS-R_env), ipx); REPROTECT(s = coerceVector(s, REALSXP), ipx); if(LENGTH(s) != n) error(_(candidate point in optim evaluated to length %d not %d), LENGTH(s), n); I think I need a little help it it is not too much to ask. Admitedly I could search for the definition of each of these macros and after some time decipher their meaning. But in an effort to save some time I am appealing to this group. First, SETCADR(OS-R_gcall, x); I was unable to find thie function call in the sources it seems to be called everywhere so searches turn up many hits. Second, PROTECT_WITH_INDEX(s = eval(OS-R_gcall, OS-R_env), ipx); I am assuming that this acutall calls the function pointed to by R_gcall. But, I am not sure where ipx fits in and what PROTECT_WITH_INDEX does. I read the wirting exstensions documentation and this specific macro is dealt with in section 5.9.1 but I still am having a hard time understanding what is happening. Perhaps the added call to eval is confiusing me. I am assuing that this makes the function call. Third, REPROTECT(s = coerceVector(s, REALSXP), ipx); Again the vector ipx shows up. Would someone please help me to understand what this statement is doing? Thank you for your time and patience. Kevin __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] number of zeros in a matrix -row by row
M - matrix(0:15, 8, 4, byrow=TRUE) rowSums(M == 0) [1] 1 0 0 0 1 0 0 0 Bill Venables http://www.cmis.csiro.au/bill.venables/ -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of onyourmark Sent: Sunday, 5 April 2009 10:00 PM To: r-help@r-project.org Subject: [R] number of zeros in a matrix -row by row Hi. I have an n x m matrix M some of who's entries are zeros. I want to know how many zeros there are in each row-perhaps stored in a 1 x n vector which lists the number of zeros for each row of M. Before I had a vector V and I was able to get the number of zeros in V by doing length(V[ V==0]) but when I try something similar for M, like M[ M==0] it creates a vector not a matrix and so this does not work. Does anyone have a solution to this? Thank you. -- View this message in context: http://www.nabble.com/number-of-zeros-in-a-matrix--row-by-row-tp22893147p22893147.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] loop problem for extract coefficients
Perhaps your loop should be more than just a print statement. That works fine! You need to place the print statement after the '{', not before it. fit$coefficients is a 21 x 20 array (the rows are lablelled 0 to 20) and you are trying to put it in the jth *column* of a 20 x 20 matrix. Not surprisingly it does not fit. Perhaps you need to know more about just what it is this software is doing. Bill Venables http://www.cmis.csiro.au/bill.venables/ -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Alex Roy Sent: Monday, 6 April 2009 12:49 AM To: r-help@r-project.org Subject: [R] loop problem for extract coefficients Dear R users, I have problem with extracting coefficients from a object. Here, X (predictor)and Y (response) are two matrix , I am regressing X ( dimensions 10 x 20) on each of columns of Y[,1] (10 x 1) and want to store the coefficient values. I have performed a Elastic Net regression and I want to store the coeffcients in each iteration. I got an error message . I do not know where is the problem Please help me. Thanks *Code:* --- library(elasticnet) X-matrix(rnorm(200),ncol=20) Y-matrix(rnorm(200),ncol=20) loop - 20 size - 20 enres-matrix(nrow = size, ncol = loop) fit-matrix(nrow = size, ncol = loop) store-matrix(nrow = size, ncol = loop) for(j in 1: 10) print (paste(j,/200,sep=)) { enres-enet(x=X,y=Y[,j],lambda=1,normalize=TRUE,intercept=TRUE) fit-predict.enet(enres, X, type=coefficients) store[,j]-fit$coefficients } library(elasticnet) Loading required package: lars X-matrix(rnorm(200),ncol=20) Y-matrix(rnorm(200),ncol=20) loop - 20 size - 20 enres-matrix(nrow = size, ncol = loop) fit-matrix(nrow = size, ncol = loop) store-matrix(nrow = size, ncol = loop) for(j in 1: 10) + print (paste(j,/200,sep=)) [1] 1/200 [1] 2/200 [1] 3/200 [1] 4/200 [1] 5/200 [1] 6/200 [1] 7/200 [1] 8/200 [1] 9/200 [1] 10/200 { + enres-enet(x=X,y=Y[,j],lambda=1,normalize=TRUE,intercept=TRUE) + fit-predict.enet(enres, X, type=coefficients) + store[,j]-fit$coefficients + } *Error in store[, j] - fit$coefficients : number of items to replace is not a multiple of replacement length * [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] lattice different colors in different areas
Hi William Deese wrote: In making grid graphs, one can make the background semi-transparent with a line like grid.rect(gp=gpar(lty=0, fill=rgb(.5, .5, 0,. 25))) and then make the area where points and lines are plotted white with lines like pushViewport(plotViewport(c(5,4,3,1))) pushViewport(dataViewport(year, m, name=plotRegion)) grid.rect(gp=gpar(fill=white)) (The area where the labels, title, legend, etc. are remains semi-transparent.) I need to use the features of lattice for some graphs but want to keep the same color theme. In lattice, the background can be changed with the lines below, but it makes the entire graph this color. bsettings=trellis.par.get(background) bsettings$col=rgb(.5, .5, 0, .25) trellis.par.set(background, bsettings) How do I now change just the plot region back to white? Like this ... ? xyplot(1:10 ~ 1:10, panel=function(...) { grid.rect(gp=gpar(fill=white)) panel.xyplot(...) }, par.settings=list(background=list(col=rgb(.5, .5, 0, .25 Paul Thanks. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Dr Paul Murrell Department of Statistics The University of Auckland Private Bag 92019 Auckland New Zealand 64 9 3737599 x85392 p...@stat.auckland.ac.nz http://www.stat.auckland.ac.nz/~paul/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] unicode only works with a second one
Hi Thomas Steiner wrote: Hi Greg and Paul, I tried several things, but I did not succeed: * I could not find the library(EBImage) on CRAN in Austria to open an png image in R. * I could not import the image via pixmap (read.pnm) as described on http://wiki.r-project.org/rwiki/doku.php?id=tips:graphics-misc:translucency because my GIMP cannot write pnm format. * I could not manage use the 'grImport' package to trace the svg: readPicture(Aries.svg) Fehler in readPicture(Aries.svg) : Version mismatch: RGML file needs to be recreated with PostScriptTrace() Yep, you need to convert to PostScript (ImageMagick or InkScape ought to do it) before you can import it. See http://www.stat.auckland.ac.nz/~paul/R/grImport/import.pdf for a thorough description of how the package works. Paul * I gave up modifiying the svg code from wikipedia to make it an R array (structure) as you greg described it above. If you have any hint for me please let me know. I am willing to contribute something to TeachingDemos (although I am not sure if this is not a license problem as I trace the (public domain) images from wikimedia. Otherwise I am happy with the Hershey fonts so far. Thomas -- Dr Paul Murrell Department of Statistics The University of Auckland Private Bag 92019 Auckland New Zealand 64 9 3737599 x85392 p...@stat.auckland.ac.nz http://www.stat.auckland.ac.nz/~paul/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Code of sample() in C
G'day Ranjan, On Sun, 5 Apr 2009 17:00:56 -0500 Ranjan Maitra mai...@iastate.edu wrote: [...[ I haven't tried this function out much myself, so YMMV. There is an obvious problem in this code since len is not decreased but n is when an index is sampled. The last line in the last for loop should be x[j] = x[--len]; Otherwise this algorithm will produce samples in which an index can be repeated. And the probabilities with which an index is sampled would be non-trivial to determine; they would definitely not be uniform. HTH. Cheers, Berwin #include stdlib.h #ifndef USING_RLIB #define MATHLIB_STANDALONE 1 /*It is essential to have this before the call to the Rmath's header file because this decides the definitions to be set. */ #endif #include Rmath.h /* Note that use of this function involves a prior call to the Rmath library to get the seeds in place. It is assumed that this is done in the calling function. */ /* Equal probability sampling; without-replacement case */ /* Adapted from the R function called SampleNoReplace */ /** * Stores k out of n indices sampled at random without replacement * in y. */ int srswor(int n, int k, int *y) { if (k n) { return 1; } else { const double len = (double) n; int i; int* x = malloc(n * sizeof(int)); if (!x) { return 2; } for (i = 0; i n; ++i) x[i] = i; for (i = 0; i k; ++i) { const int j = (int)(len * runif(0.0, 1.0)); y[i] = x[j]; x[j] = x[--n]; } free(x); } return 0; } === Full address = Berwin A TurlachTel.: +65 6516 4416 (secr) Dept of Statistics and Applied Probability+65 6516 6650 (self) Faculty of Science FAX : +65 6872 3919 National University of Singapore 6 Science Drive 2, Blk S16, Level 7 e-mail: sta...@nus.edu.sg Singapore 117546http://www.stat.nus.edu.sg/~statba __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] sasweave, and R 2.5.1 vs 2.8.1
Looking for advice on a problem with Sweave, sasweave, and version 2.5.1 versus 2.8.1 of R. I'm running all this on WinXP. I have both 2.5.1 and 2.8.1 installed. I have sasweave set up to run R 2.5.1, as this line in the runR.bat file shows: set RTERM=C:\Program Files\R\R-2.5.1\bin\Rterm.exe And things work fine. For example, this minimal file compiles OK: \documentclass{article} \usepackage{Sweave} \usepackage{amsmath} \begin{document} options, echo=FALSE, hide=TRUE= setwd(C:/data/profdev/msbiostatistics/tamu/stat659/homework/homework8) options(SweaveSyntax = SweaveSyntaxNoweb) @ \Sexpr{3+2} \end{document} But if I change to R version 2.8.1, by making the change in runR.bat in my sasweave directory: set RTERM=C:\Program Files\R\R-2.8.1\bin\Rterm.exe then none of my \Sexpr{} expressions work. Near as I can tell, everything else works fine. Anyone else experienced anything similar, and can advise? Thanks. --Chris Ryan __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to sort and plot data?
hi Erin, Thanks for your reply to my problems. I tried and it works, but it sorted all the column. In my case, I want it to sort d user_id column and another information in other column will follow it. Now, after I sorted, the information which website viewed by user was wrong. I want to as how we can sort or filter in excel. Thank you. regards, hema On Fri, Apr 3, 2009 at 4:44 PM, Umesh Srinivasan umesh.sriniva...@gmail.com wrote: Hi, There is definitely a more elegant way of doing this which I don't know about (without a for loop), but try this: mat - matrix(NA, nrow = max(user_id), ncol = 2) mat[,1] - 1:max(user_id) # 1st column of matrix is the user ID for (i in 1:max(user_id)){ temp1 - subset(data, user_id = i) temp2 - unique(temp1$website) mat[2,i] - length(temp2) } The matrix will give you user id and number of sites visited, provided user id ranges from 1 to the number of users. There must be a way to do this using table, but I cant figure it out. Cheers, Umesh On Fri, Apr 3, 2009 at 1:42 PM, Dieter Menne dieter.me...@menne-biomed.dewrote: Hem wrote: user_id website time 20google0930 21yahoo0935 20facebook1000 25facebook1015 61google0940 ... My problem is how to sort the data? So that, I can get information about one user_id viewed how many website perday? Maybe you were looking at the wrong item, because what you want is not sorting, but a table. Check the documentation of table or ftable. Dieter -- View this message in context: http://www.nabble.com/how-to-sort-and-plot-data--tp22861661p22863918.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] predicting values into the future
Thank you Bill and Gabor: I do have a few years of fish sizes data (from larvae to juvenile). If I melt them together how can I create the best model to predict future weights. My weeks go up to 16 right now, but I want to predict weights for week 17 to 20. Based on both of you examples it was easy to fit a model with one year but I can't figure out how to create a model with all the years together. Perhaps my data isn't organized the correct way. See the data below: mydf - read.table(textConnection(first second third fourth fifth sixth seventh 0.08 0.003 0.1 0.003 0.003 0.003 0.002 0.3 0.11 0.6 0.06 0.13 0.02 0.12 1.78 0.33 2.1 0.26 0.69 0.23 0.34 2.91 0.63 4 0.76 1.51 0.87 1.08 4.4 1.51 6.2 1.63 2.57 1.27 1.91 5.5 2.34 8.17 2.49 3.2 2.65 2.8 6.69 3.3 9.64 3.6 5 3.91 4.17 7.8 5 12.1 5.5 6.2 4.9 5.2 9.1 6.2 15 7 7.7 6.3 6.7 10.6 7.7 16.5 8.5 9.2 7.8 8.2 12.1 9.2 18 10 11.6 9.4 9.7 13.6 12.2 19.5 12.9 13.3 10.9 11.2 16.8 13.7 22.2 14.4 15 12.9 13.2 18 15.7 23.7 15.9 16 13.9 14.2 20 16.7 26 16.9 16 17.2 17.9 21 18 27.2 18.1 17 18.3 19.1),header=T) mydf mydf$week - 1:16 mydf - melt(mydf) colnames(mydf) - c('week','year','weights') attach(mydf) --- On Sun, 4/5/09, bill.venab...@csiro.au bill.venab...@csiro.au wrote: From: bill.venab...@csiro.au bill.venab...@csiro.au Subject: RE: [R] predicting values into the future To: mazatlanmex...@yahoo.com, r-h...@stat.math.ethz.ch Date: Sunday, April 5, 2009, 12:09 AM Here is a bit of an exploration of your data but first a couple of notes. * the information about Excel is probably a bit superfluous here. Some of us have no idea about Excel, and rather hope it can stay that way. * With such a short series, you don't stand much chance of fitting a time series model such as with arima. It's clearly not stationary, too. If you had multiple growth curves you may stand some chance of fitting a correlated model, but with just one, I don't think so. For now, I think you just may have to make the hopeful assumption of independence. You might like to look at this. weightData - data.frame(weight = c(2.1,2.4,2.8,3.6,4.1,5.2,6.3), week = 1:7) plot(weight ~ week, weightData) plot(log(weight) ~ week, weightData) ### clearly the log plot seems to linearise things. ### Try an non-linear regression: wModel - nls(weight ~ alpha + beta*exp(gamma*week), weightData, start = c(alpha = 0.0, beta = 1, gamma = 0.2), trace = TRUE) you should look at the residuals from this to see if the assumptions look reasonable. With only 7, you can't see much, though. now suppose you want to predict for another 3 weeks: newData - data.frame(week = 1:10) newData$pweight - predict(wModel, newData) plot(pweight ~ week, newData, pch = 4, col = red, ylab = Weight, xlab = Week) with(weightData, points(week, weight)) looks OK to me (thought fish cannot keep on growing exponentailly forever - this is clearly a model with limitations and you have to be careful when pushing it too far). finally predict on a more continuous scale and add in the result as a blue line. lData - data.frame(week = seq(1, 10, len = 1000)) with(lData, lines(week, predict(wModel, lData), col = blue)) Now that we have over-analysed this miniscule data set to blazes, perhaps it's time for a beer! __ Bill Venables http://www.cmis.csiro.au/bill.venables/ -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Felipe Carrillo Sent: Sunday, 5 April 2009 4:13 PM To: r-h...@stat.math.ethz.ch Subject: [R] predicting values into the future Hi: I have usually used the GROWTH() excel function to do this but now want to see if I can do this with R. I want to predict values into the future, possibly with the predict.arima Function. I have the following weekly fish weight averages: weight - c(2.1,2.4,2.8,3.6,4.1,5.2,6.3) week - c(1,2,3,4,5,6,7) I would like to predict what the weight will be by week 10 based on my weight values and make a line plot of all the weights(including the predicted values). I have two questions: 1- Should the predicted values be linear or exponential? 2- Is the predict.arima function appropriate to do this? Thanks in advance. Felipe D. Carrillo Supervisory Fishery Biologist Department of the Interior US Fish Wildlife Service California, USA __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read
Re: [R] predicting values into the future
Thank you Bill and Gabor: I do have a few years of fish sizes data (from larvae to juvenile). If I melt them together how can I create the best model to predict future weights. My weeks go up to 16 right now, but I want to predict weights for week 17 to 20. Based on both of you examples it was easy to fit a model with one year but I can't figure out how to create a model with all the years together.See the data below: mydf - read.table(textConnection(first second third fourth fifth sixth seventh 0.08 0.003 0.1 0.003 0.003 0.003 0.002 0.3 0.11 0.6 0.06 0.13 0.02 0.12 1.78 0.33 2.1 0.26 0.69 0.23 0.34 2.91 0.63 4 0.76 1.51 0.87 1.08 4.4 1.51 6.2 1.63 2.57 1.27 1.91 5.5 2.34 8.17 2.49 3.2 2.65 2.8 6.69 3.3 9.64 3.6 5 3.91 4.17 7.8 5 12.1 5.5 6.2 4.9 5.2 9.1 6.2 15 7 7.7 6.3 6.7 10.6 7.7 16.5 8.5 9.2 7.8 8.2 12.1 9.2 18 10 11.6 9.4 9.7 13.6 12.2 19.5 12.9 13.3 10.9 11.2 16.8 13.7 22.2 14.4 15 12.9 13.2 18 15.7 23.7 15.9 16 13.9 14.2 20 16.7 26 16.9 16 17.2 17.9 21 18 27.2 18.1 17 18.3 19.1),header=T) mydf mydf$week - 1:16 mydf - melt(mydf) colnames(mydf) - c('week','year','weights') attach(mydf) * With such a short series, you don't stand much chance of fitting a time series model such as with arima. It's clearly not stationary, too. If you had multiple growth curves you may stand some chance of fitting a correlated model, but with just one, I don't think so. For now, I think you just may have to make the hopeful assumption of independence. You might like to look at this. weightData - data.frame(weight = c(2.1,2.4,2.8,3.6,4.1,5.2,6.3), week = 1:7) plot(weight ~ week, weightData) plot(log(weight) ~ week, weightData) ### clearly the log plot seems to linearise things. ### Try an non-linear regression: wModel - nls(weight ~ alpha + beta*exp(gamma*week), weightData, start = c(alpha = 0.0, beta = 1, gamma = 0.2), trace = TRUE) you should look at the residuals from this to see if the assumptions look reasonable. With only 7, you can't see much, though. now suppose you want to predict for another 3 weeks: newData - data.frame(week = 1:10) newData$pweight - predict(wModel, newData) plot(pweight ~ week, newData, pch = 4, col = red, ylab = Weight, xlab = Week) with(weightData, points(week, weight)) looks OK to me (thought fish cannot keep on growing exponentailly forever - this is clearly a model with limitations and you have to be careful when pushing it too far). finally predict on a more continuous scale and add in the result as a blue line. lData - data.frame(week = seq(1, 10, len = 1000)) with(lData, lines(week, predict(wModel, lData), col = blue)) Now that we have over-analysed this miniscule data set to blazes, perhaps it's time for a beer! __ Bill Venables http://www.cmis.csiro.au/bill.venables/ Hi: I have usually used the GROWTH() excel function to do this but now want to see if I can do this with R. I want to predict values into the future, possibly with the predict.arima Function. I have the following weekly fish weight averages: weight - c(2.1,2.4,2.8,3.6,4.1,5.2,6.3) week - c(1,2,3,4,5,6,7) I would like to predict what the weight will be by week 10 based on my weight values and make a line plot of all the weights(including the predicted values). I have two questions: 1- Should the predicted values be linear or exponential? 2- Is the predict.arima function appropriate to do this? Thanks in advance. Felipe D. Carrillo Supervisory Fishery Biologist Department of the Interior US Fish Wildlife Service California, USA __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Need help in calculating studentized residuals/leverage values of non-linear model [nls()]
Hi there, I hope I can get advice regarding the calculation of leverage values or studentized residual values of a non-linear regression model. It seems like rstudent() does not work on a nls object. Many thanks in advance! Best regards, Xingli __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] applying function to non zero elements of a sparse Matrix
The class Matrix (esp. dgCMatrix) is extraordinary. However, I still have growin pains with it. One common need (I believe) is to apply the same function to all non-zero elements of the sparse matrix. Is there a simple, elegant way? -g __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] number of zeros in a matrix -row by row
Dear Dimitris, Thanks very much. I tried dim(M) and get: [1] 5030 2142 which I thought meant that my matrix M is actually a matrix, but I tried str(M) and get: int [1:5030, 1:2142] 2 1 2 0 1 0 2 1 2 -1 ... When I tried rowSums(M == 0) I got: Error in rowSums[M == 0] : object of type 'closure' is not subsettable What am I doing wrong? Thanks so much. Dimitris Rizopoulos-4 wrote: if 'M' is your matrix, then try this: rowSums(M == 0) I hope it helps. Best, Dimitris onyourmark wrote: Hi. I have an n x m matrix M some of who's entries are zeros. I want to know how many zeros there are in each row-perhaps stored in a 1 x n vector which lists the number of zeros for each row of M. Before I had a vector V and I was able to get the number of zeros in V by doing length(V[ V==0]) but when I try something similar for M, like M[ M==0] it creates a vector not a matrix and so this does not work. Does anyone have a solution to this? Thank you. -- Dimitris Rizopoulos Assistant Professor Department of Biostatistics Erasmus University Medical Center Address: PO Box 2040, 3000 CA Rotterdam, the Netherlands Tel: +31/(0)10/7043478 Fax: +31/(0)10/7043014 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- View this message in context: http://www.nabble.com/number-of-zeros-in-a-matrix--row-by-row-tp22893147p22901010.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] line argument in mtext for axis ?
Is there a similar argument for axis that controls the position of labels via line argument in mtext ? __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.