[R] Repeated Measures ANOVA and Missing Values in the data set
I am doing Repeated Measures ANOVA with missing values. When i run my model i get this error message. *aov.out = aov(values ~ time + Error(subject/time), data=mydata2)Warning message:In aov(values ~ time + Error(subject/time), data = mydata2) : Error() model is singular* The missing Values are not a error of my instrument. They mean the element of my analysis is absent and i want to consider this. thanks in advance these are my data: subject <- c(1,2,3,4,5,6,7,8,9,10) time1 <- c(5040,3637,6384,5309,5420,3549,NA,5140,3890,3910) time2 <- c(5067, 3668, NA, 6489, NA, 3922, 3408, 6613, 4063, 3937) time3 <- c( 3278, 3814, 8745, 4760, 4911, 5716, 5547, 5844, 4914, 4390) time4 <- c( 0, 2971,0, 2776, 2128, 1208, 2935, 2739, 3054, 3363) time5 <- c(4161, 3483, 6728, 5008, 5562, 4380, 4006, 7536, 3805, 3923) time6 <- c( 3604, 3411, 2523, 3264, 3578, 2941, 2939, NA, 3612, 3604) mydata <- data.frame(time1, time2, time3, time4, time5, time6) mydata2 = stack(mydata) subject = factor(rep(subject,6)) mydata2[3] = subject colnames(mydata2) = c("values", "time", "subject") aov.out = aov(values ~ time + Error(subject/time), data=mydata2) summary(aov.out) model.tables(aov.out,"means") [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/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] Repeated Measures ANOVA and the Bonferroni post hoc test different results of significantly
Hi All and thanks for Help I am doing an Repeated Measures ANOVA and the Bonferroni post hoc test for my data using R project. The ANOVA gives a significantly difference between the data but not the Bonferroni post hoc test. > anova(aov2) numDF denDF F-value p-value (Intercept) 1 1366 110.51125 <.0001 time5 1366 9.84684 <.0001 while > pairwise.t.test(x=table.metric2$value, g=table.metric2$time, p.adj="bonf") Pairwise comparisons using t tests with pooled SD data: table.metric2$value and table.metric2$time Su1 Su2 Su3 Su4 Su5 Su2 1 - - - - Su3 1 1 - - - Su4 1 1 1 - - Su5 1 1 1 1 - Su6 1 1 1 1 1 P value adjustment method: bonferroni These are my data with the code used plot <- c(rep(1:275)) Su1 <- c(13.5584,0.,2.0710,0.4826,1.2761,1.6690,3.5188, 13.7578,0.,0.,0.0004,0.,0.,0.,4.4634,3.0151,2.1719,5.2861,4.9651,0.7908,0.,0.,0.,0., 0.,5.2749,5.4706,4.4416,3.2166,0.,0.1929,0.,0.,0.,0.,0.,4.6765,1.7761,4.3579,0.,0.,0.,0.,0.,0.,0.,0.,0., 0.,6.1794,2.4194,1.4319,0.,0.0327,0.,0.4633,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.0018,0.,0., 0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0., 0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.5425,0.5274,2.0883,0.6024,0.,0.0018,0.,0.,0.0015,0.0012,0.,0.,0., 0.4560,2.1256,0.0295,0.0328,0.,1.6447,0.0428,0.0067,0.0058,0.,0.,0.,0.0001,0.3317,0.2898,3.5134,0.1539,0.0199,0.,0.0494,0.1159,2.0976,0.0644,0.9730, 0.0010,0.5074,0.0003,0.,0.1188,1.8818,0.,0.1213,0.3585,7.3932,0.5492,0.0045,0.9879,0.0010,0.7625,0.1695,0.1211,0.3164,2.6750,3.8926,0.,3.4626,0.,5.8339, 6.7315,0.0244,4.8770,2.6237,2.3700,0.5338,0.0215,3.2196,1.9811,3.3825,3.3929,1.5426,0.9165, 10.6561,3.2154,4.1531,5.3381,3.9432,4.8675,0.0047,0.0026,0.2058,1.8509,0.3697, 0.3131,0.0707,4.7908,6.4087, 10.3670,5.7662,4.0460,3.2571,9.1767,0.0116,0.0908,0.0053,0.1480,0.9063,5.4331,5.7945, 14.4097,6.9635,7.0637,0.1064,9.9095, 11.8432, 10.0234,0., 0.0491,5.0472,5.3094,5.1657, 14.3944,7.6244,0.0034,1.4953, 14.7658,6.1775,7.1567,0.0296,0.0911,3.5552,4.9543,3.1200,1.9774,0.,4.1663,0.,2.3672,0.0638,1.8952,4.1948, 6.4229, 10.7573,0.0008,1.3818,6.0011,3.6791,9.7816,1.5203,0.8616,1.5483,5.4174,2.7070,2.1627,0.,1.7360,3.7292,2.4638,7.4498,4.2343,6.8263,3.2410,0.0001,0.0001,9.7424, 4.2861,2.9912,0.4316, 11.6082,2.0138,0.0002,1.8783,0.9934,0.2983,1.4013,0.1429) Su2 <- c(10.4361,0.,2.3346,0.5769,1.3392,1.5908,3.5759, 13.3183,0.0005,0.,0.0019,0.,0.,0.,4.4862,3.0418,2.3991,5.4263,4.9456,0.6907,0.,0.0007,0.,0., 0.,5.1065,5.1748,4.2888,3.3633,0.,0.1791,0.,0.,0.,0.,0.,5.3485,1.9216,4.2880,0.,0.0001,0.0001,0.,0.,0.,0.,0.,0., 0.,6.3664,2.2888,1.4636,0.,0.0282,0.,0.5838,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.0013,0.,0., 0.,0.,0.,0.,0.,0.,0.,0.,0.,0.0002,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0., 0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.6240,0.6401,1.9324,0.6263,0.,0.0020,0.,0.0001,0.0039,0.0018,0.,0.,0., 0.4321,1.9815,0.0528,0.0350,0.,1.6519,0.0328,0.0153,0.0171,0.,0.,0.,0.,0.3966,0.2957,3.4653,0.1473,0.0038,0.,0.0774,0.1180,2.2780,0.0611,0.9490, 0.0024,0.8337,0.,0.,0.1531,1.9778,0.,0.1171,0.3485,7.2378,0.7476,0.0028,1.2072,0.0015,0.7425,0.1352,0.1908,0.3092,2.3735,3.9768,0.,3.8786,0.,6.5420, 6.8490,0.0256,5.8268,2.7856,2.4640,0.6399,0.0215,5.3411,1.7939,3.3401,3.2942,1.4990,0.9264, 10.6864,3.3749,4.7978,5.2929,3.8639,5.3890,0.0027,0.0486,0.2105,1.8514,0.3526, 0.3146,0.0950,4.8061,6.2244, 10.1131,5.8538,3.8861,4.4240,9.5952,0.0255,0.0533,0.0085,0.1742,1.0188,7.7153,5.4663,14.6060,6.8725,6.7284,0.0841, 10.2016, 11.3384,9.8938,0., 0.0749,5.0774,7.2876,5.2040, 14.7609,7.7862,0.0005,1.4099, 14.8639,6.8801,7.2587,0.0125,0.0513,3.1338,7.1539,3.1733,1.9729,0.,4.5081,0.,2.1572,0.0452,1.8866,8.0198, 9.7868, 10.8746,0.,1.5011,6.0825,3.6705,9.9171,1.7091,0.8267,1.4186,6.4235,2.9303,1.9019,0.,1.5869,4.4028,2.4186,7.7739,4.6728,8.7722,3.9859,0.0001,0.0001, 10.1583, 5.4758,2.8977,0.9716, 11.5342,2.6111,0.,1.7497,1.0041,0.3273,1.1953,0.2289) Su3 <- c(12.5808,0.,2.7244,0.1119,1.2633,0.9702,2.4513,2.4419,0.0014,0.,0.,0.,0.0019,0.0105,3.5392,2.4033,2.1847,5.5912,4.8628,0
[R] the right reference for the R Stats package for a scientific journal
Dear Members list, I am writing a paper for a research where i used "the R Stats package". No one knows the right reference for this package? Thanks in Advance Gianni [[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] total number of citations for R project
Dear Member list, is there a weblink or a paper where the total number of citations for R project is report? Thanks in advance Gianni [[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] function to calculate a SMAPE (Symmetric mean absolute percentage error) to avoid the possibility of an inflation caused by zero values in the series
Dear Reseacher, i find this modified version of SMAPE at pag 13 formula "msSMAPE" to to avoid the possibility of an inï¬ation caused by zero values in the series using a Si component in the denominator of the symmetric MAPE http://www.stat.iastate.edu/preprint/articles/2004-10.pdf I wrote a function but I wish eventually correct error, change or improve with your suggestions. Please feel free to add comments or remarks because I am not a mathematicians and probably i wrong to write this function Thanks in advance Gianni SMAPEper <- function(x,y){mean((200*(abs(x-y)))/(x+y))} obs <- c(10,20,0,40) pred <-c(5,6,0,8) SMAPEper(obs,pred) msSMAPE <- function(obs,pred){ M <- numeric(length(obs)) for (i in 2:length(obs)) { n <- i-1 M[i] <- mean(obs[1:n]) } sum1 <- (abs(obs[1]-pred[1]))/(0.5*(abs(pred[1])+abs(obs[1]))) for (i in 2:length(obs)) { n <- i-1 sum2 <- 0 for (k in 1:n) {sum2 <- sum2+abs(obs[k]-M[k])} S <- sum2/n sum1 <- sum1+(abs(obs[i]-pred[i]))/(0.5*(abs(pred[i])+abs(obs[i]))+S) } my.SMAPE <- sum1/length(obs) return(my.SMAPE) } msSMAPE(obs,pred) [[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 to find or fix a SMAPE (Symmetric mean absolute percentage error ) with correction to avoid zero value in the series
Dear Researchers, I am looking for a library or a function to calculate SMAPE with same doneminator to avoid a problem to the presence of 0 values in the series. I find this variotion of SMAPE called mSMAPE mSMAPE page 13 http://www.stat.iastate.edu/preprint/articles/2004-10.pdf yesterday nigth i tried to create a function but without a good result. For this reason I wish to know if some reserchers develop a solution to avoid the possibility of an inï¬ation of sMAPE caused by zero values in the series. Thanks in advance Gianni # this is my mSMAPE function but there is some problem SMAPEper <- function(x,y){mean((200*(abs(x-y)))/(x+y))} obs <- c(10,20,30,40) pred <-c(5,6,7,8) mSMAPE <- function(obs,pred){ M <- numeric(length(obs)) for (i in 2:length(obs)) { n <- i-1 M[i] <- mean(obs[1:n]) } sum1 <- (abs(obs[1]-pred[1]))/(0.5*(abs(pred[1])+abs(obs[1]))) for (i in 2:length(obs)) { n <- i-1 sum2 <- 0 for (k in 1:n) {sum2 <- sum2+abs(obs[k]-M[k])} S <- sum2/n sum1 <- sum1+(abs(obs[i]-pred[i]))/(0.5*(abs(pred[i])+abs(obs[i]))+S) } my.SMAPE <- sum1/length(obs) return(my.SMAPE) } mSMAPE(obs,pred)*100 SMAPEper(obs,pred) [[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] some help to improve "hist to plot relative frequencies"
Dear Researches, sorry for disturb. I wish to improve my figure in R plotting the relative frequencies of my data set. library(lattice) a <- c(0,0,0,1,1,2,4,5,6,7,7,7,7,7,8,8,8,8,9,9,9,9,10,10,11) histogram(a, xlab="myData") what i wish to do is: 1) invert the order of X and Y (eg: Precent of Total on X-axis and "MyData" on X-axis) 2) plot not the bar of histogram but a line (i tried with "lines(density(a))" but the result is not what i wish) Great Thanks for any helps Gianni [[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] special simbol (±) in a legend
dear reserachers, I am looking for a expression about a special symbol (±) in a Legend and add on front of "Mean ± SD" the symbol of bracket sorry if the example is not great disptest <- matrix(rnorm(10),nrow=10) disptest.means<- rowMeans(disptest) plot(1:10,disptest.means) dispersion(1:10,disptest.means,(disptest.means+disptest ),(disptest.means-disptest ),intervals=FALSE,arrow.cap=0.001,arrow.gap=0) legend("topleft", legend=c("Mean","Mean +/- SD"), pch=c(1,NA), bty="n",cex=1.3) thanks for all suggestions Gianni [[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] error to convert a Compute A^-1 B from Matlab to R using solve(A, B)
Dear Researchers, I need to convert the following equation in R from Matlab a = [x y ones(size(x))]; b = [-(x.^2+y.^2)]; a\b ans = -9.9981 -16.4966 -7.6646 my solution in R is: a = cbind(x,y,rep(1,length(x))) b = cbind(-(x^2+y^2)) > head(a) xy [1,] 14.45319 5.065726 1 [2,] 14.99478 5.173893 1 [3,] 14.64158 5.616916 1 [4,] 14.61803 6.624069 1 [5,] 14.19997 6.794587 1 [6,] 15.08174 8.224843 1 > head(b) [,1] [1,] -234.5564 [2,] -251.6125 [3,] -245.9255 [4,] -257.5652 [5,] -247.8057 [6,] -295.1068 following MATLAB/ R Reference http://cran.r-project.org/doc/contrib/Hiebeler-matlabR.pdf the a\b could be converted by solve(a,b) but i get the following error: > solve(a,b) Error in solve.default(a, b) : 'b' must be compatible with 'a' thanks for any help Gianni [[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] Fit circle with R
Dear Researchers, I wrote two function to fit a circle using noisy data. 1- the fitCircle() is derived from MATLAB code of * zhak Bucher* from the link http://www.mathworks.com/matlabcentral/fileexchange/5557-circle-fit/content/circfit.m 2- the CircleFitByPratt() from MATLAB code of *Nikolai Chernov *from the link http://www.mathworks.com/matlabcentral/fileexchange/22643-circle-fit-pratt-method/content/CircleFitByPratt.m, based on: *V. Pratt, "Direct least-squares fitting of algebraic surfaces", Computer Graphics, Vol. 21, pages 145-152 (1987)* I am looking for new methods to compare and improve my analysis because the error increase with decreasing of points used in the functions. Thanks for all suggestions Gianni Here the funtions with example # fitCircle, returns: # xf,yf = centre of the fitted circle # Rf = radius of the fitted circle # Cf = circumference of the fitted circle # Af = Area of the fitted circle fitCircle <- function(x,y){ a = qr.solve(cbind(x,y,rep(1,length(x))),cbind(-(x^2+y^2))) xf = -.5*a[1] yf = -.5*a[2] Rf = sqrt((a[1]^2+a[2]^2)/4-a[3]) Cf = 2*pi*Rf Af = pi*(Rf^2) m <- cbind(xf,yf,Rf,Cf,Af) return(m)} # CircleFitByPratt, returns: # [,1] and [,2] = centre of the fitted circle # [,3] = radius of fitted cirlce CircleFitByPratt <- function(x,y){ n <- length(x) centroid <- cbind(mean(x),mean(y)) Mxx=0; Myy=0; Mxy=0; Mxz=0; Myz=0; Mzz=0; for(i in 1:n){ Xi <- x[[i]] - centroid[1] Yi <- y[[i]] - centroid[2] Zi <- (Xi*Xi) + (Yi*Yi) Mxy = Mxy + Xi*Yi; Mxx = Mxx + Xi*Xi; Myy = Myy + Yi*Yi; Mxz = Mxz + Xi*Zi; Myz = Myz + Yi*Zi; Mzz = Mzz + Zi*Zi; } Mxx = Mxx/n Myy = Myy/n Mxy = Mxy/n Mxz = Mxz/n Myz = Myz/n Mzz = Mzz/n # computing the coefficients of the characteristic polynomial Mz = Mxx + Myy; Cov_xy = Mxx*Myy - Mxy*Mxy; Mxz2 = Mxz*Mxz; Myz2 = Myz*Myz; A2 = 4*Cov_xy - 3*Mz*Mz - Mzz; A1 = Mzz*Mz + 4*Cov_xy*Mz - Mxz2 - Myz2 - Mz*Mz*Mz; A0 = Mxz2*Myy + Myz2*Mxx - Mzz*Cov_xy - 2*Mxz*Myz*Mxy + Mz*Mz*Cov_xy; A22 = A2 + A2; epsilon=1e-12; ynew=1e+20; IterMax=20; xnew = 0; # Newton's method starting at x=0 epsilon=1e-12; ynew=1e+20; IterMax=20; xnew = 0; iter=1:IterMax for (i in 1:IterMax){ yold = ynew; ynew = A0 + xnew*(A1 + xnew*(A2 + 4.*xnew*xnew)); if (abs(ynew) > abs(yold)){ print('Newton-Pratt goes wrong direction: |ynew| > |yold|') xnew = 0; break } Dy = A1 + xnew*(A22 + 16*xnew*xnew); xold = xnew; xnew = xold - ynew/Dy; if (abs((xnew-xold)/xnew) < epsilon) {break} if(iter[[i]] >= IterMax){ print('Newton-Pratt will not converge'); xnew = 0; } if(xnew < 0.){ print('Newton-Pratt negative root: x=',xnew); } } DET = xnew*xnew - xnew*Mz + Cov_xy; Center = cbind(Mxz*(Myy-xnew)-Myz*Mxy , Myz*(Mxx-xnew)-Mxz*Mxy)/DET/2; #computing the circle parameters DET = xnew*xnew - xnew*Mz + Cov_xy; Center = cbind(Mxz*(Myy-xnew)-Myz*Mxy , Myz*(Mxx-xnew)-Mxz*Mxy)/DET/2; Par = cbind(Center+centroid , sqrt(Center[2]*Center[2]+Mz+2*xnew)); return(Par) } #EXAMPLE library(plotrix) # Create a Circle of radius=10, centre=5,5 R = 10; x_c = 5; y_c = 5; thetas = seq(0,pi,(pi/64)); xs = x_c + R*cos(thetas) ys = y_c + R*sin(thetas) # Now add some random noise mult = 0.5; xs = xs+mult*rnorm(rnorm(xs)); ys = ys+mult*rnorm(rnorm(ys)); plot(xs,ys,pch=19,cex=0.5,col="red",xlim=c(-10,20),ylim=c(-10,20),asp=1) # real circle draw.circle(x_c,y_c,radius=10,border="black") points(x_c,y_c,,pch=4,col="black") CPrat <- CircleFitByPratt(xs,ys) draw.circle(CPrat[1],CPrat[2],radius=CPrat[3],border="blue") points(CPrat[1],CPrat[2],pch=4,col="blue") MyC <- fitCircle(xs,ys) draw.circle(MyC[1],MyC[2],radius=MyC[3],border="green") points(MyC[1],MyC[2],pch=4,col="green") # Select less points points(xs[20:49],ys[20:49]) MyC1 <- fitCircle(xs[20:49],ys[20:49]) draw.circle(MyC1[1],MyC1[2],radius=MyC1[3],border="blue",lty=2,lwd=2) CPrat1 <- CircleFitByPratt(xs[20:49],ys[20:49]) draw.circle(CPrat1[1],CPrat1[2],radius=CPrat1[3],border="green",lty=2,lwd=2) points(CPrat[1],CPrat[2],pch=4,col="red") [[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] Based-Density cluster algorithms, Help where find some algorithms
Dear Researches, Did anyone know where to find the following Based-Density cluster algorithms in R or function implemented? OPTICS (Ordering Points To Identify the Clustering Structure) DECLUN (DENsity CLUstering) Thanks in Advance Gianni [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] problem to tunning RandomForest, an unexpected result
Dear Researches, I am using RF (in regression way) for analize several metrics extract from image. I am tuning RF setting a loop using different range of mtry, tree and nodesize using the lower value of MSE-OOB mtry from 1 to 5 nodesize from1 to 10 tree from 1 to 500 using this paper as refery Palmer, D. S., O'Boyle, N. M., Glen, R. C., & Mitchell, J. B. O. (2007). Random Forest Models To Predict Aqueous Solubility. Journal of Chemical Information and Modeling, 47, 150-158. my problem is the following using data(airquality) : the tunning parameters with the lower value is: > print(result.mtry.df[result. mtry.df$RMSE == min(result.mtry.df$RMSE),]) *RMSE = 15.44751 MSE = 238.6257 mtry = 3 nodesize = 5 tree = 35* the numer of tree is very low, different respect how i can read in several pubblications And the second value lower is a tunning parameters with *tree = 1* print(head(result.mtry.df[with(result.mtry.df, order(MSE)), ])) RMSE MSE mtry nodesize tree 12035 15.44751 238.625735 35 *18001 15.44861 238.6595471 *7018 16.02354 256.753925 18 20031 16.02536 256.812151 31 11037 16.02862 256.916533 37 11612 16.05162 257.654434 112 i am wondering if i wrong in the setting or there are some aspects i don't conseder. thanks for attention and thanks in advance for suggestions and help Gianni require(randomForest) data(airquality) set.seed(131) MyOzone <- data.frame(na.omit(airquality)) str(MyOzone) #all data My.mtry=c(1,2,3,4,5) My.nodesize=c(seq(1,10,by=1)) My.tree=c(seq(1,500,by=1)) # {} result.mtry <- list() for(i in 1:length(My.mtry)){ result.nodesize <- list() for(m in 1:length(My.nodesize)){ result.tree <- list() for(l in 1:length(My.tree)){ ozone.rf <- randomForest(MyOzone[,-c(1)],MyOzone[,c(1)], mtry=My.mtry[[i]], nodesize=My.nodesize[[m]], ntree=My.tree[[l]], importance=TRUE) result.tree[[l]] <- data.frame(RMSE=sqrt(ozone.rf$mse[ozone.rf$ntree]),MSE=ozone.rf$mse[ozone.rf$ntree],mtry=My.mtry[[i]],nodesize=My.nodesize[[m]],tree=My.tree[[l]]) } result.tree.df <- do.call(rbind,result.tree) result.nodesize[[m]] <- result.tree.df } result.nodesize.df <- do.call(rbind,result.nodesize) result.mtry[[i]] <- result.nodesize.df } result.mtry.df <- do.call(rbind,result.mtry) print(result.mtry.df[result.mtry.df$RMSE == min(result.mtry.df$RMSE),]) RMSE MSE mtry nodesize tree 12035 15.44751 238.625735 35 head(result.mtry.df[with(result.mtry.df, order(RMSE)), ]) RMSE MSE mtry nodesize tree 12035 15.44751 238.625735 35 18001 15.44861 238.6595471 7018 16.02354 256.753925 18 20031 16.02536 256.812151 31 11037 16.02862 256.916533 37 11612 16.05162 257.654434 112 [[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] tuning random forest. An unexpected result
Dear Researches, I am using RF (in regression way) for analize several metrics extract from image. I am tuning RF setting a loop using different range of mtry, tree and nodesize using the lower value of MSE-OOB mtry from 1 to 5 nodesize from1 to 10 tree from 1 to 500 using this paper as refery Palmer, D. S., O'Boyle, N. M., Glen, R. C., & Mitchell, J. B. O. (2007). Random Forest Models To Predict Aqueous Solubility. Journal of Chemical Information and Modeling, 47, 150-158. my problem is the following using data(airquality) : the tunning parameters with the lower value is: > print(result.mtry.df[result.mtry.df$RMSE == min(result.mtry.df$RMSE),]) *RMSE = 15.44751 MSE = 238.6257 mtry = 3 nodesize = 5 tree = 35* the numer of tree is very low, different respect how i can read in several pubblications And the second value lower is a tunning parameters with *tree = 1* print(head(result.mtry.df[ with(result.mtry.df, order(MSE)), ])) RMSE MSE mtry nodesize tree 12035 15.44751 238.625735 35 *18001 15.44861 238.6595471 *7018 16.02354 256.753925 18 20031 16.02536 256.812151 31 11037 16.02862 256.916533 37 11612 16.05162 257.654434 112 i am wondering if i wrong in the setting or there are some aspects i don't conseder. thanks for attention and thanks in advance for suggestions and help Gianni [[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 to setting a multiple (linear) regression model with a 5% significance level (threshold) for the inclusion of the model variables.
Dear Researchers, someone know the right syntax to chose a 5% signiï¬cance level (threshold) for the inclusion of the model variables in a multiple (linear) regression in backward way? I set the formula in this way, but I don't know to choose the 5% significance? lmodelV <- step(lm(formula=MyFormula[[1]],data=LR.train),direction="backward") thanks in advance gianni [[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] explanation why RandomForest don't require a transformations (e.g. logarithmic) of variables
Dear Researches, sorry for the easy and common question. I am trying to justify the idea of RandomForest don't require a transformations (e.g. logarithmic) of variables, comparing this non parametrics method with e.g. the linear regressions. In leteruature to study my phenomena i need to apply a logarithmic trasformation to describe my model, but i found RF don't required this approach. Some people could suggest me text or bibliography to study? thanks in advance Gianni [[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] explanation why RandomForest don't require a transformations (e.g. logarithmic) of variables
about the " because they only use the ranks of the variables". Using a leave-one-out, in each interaction the the predictor variable ranks change slightly every time RF builds the model, especially for the variables with low importance. Is It correct to justify this because there are random splitting? Thanks in advance Gianni On Mon, Dec 5, 2011 at 7:59 PM, Liaw, Andy wrote: > Tree based models (such as RF) are invriant to monotonic transformations > in the predictor (x) variables, because they only use the ranks of the > variables, not their actual values. More specifically, they look for > splits that are at the mid-points of unique values. Thus the resulting > trees are basically identical regardless of how you transform the x > variables. > > Of course, the only, probably minor, differences is, e.g., mid-points can > be different between the original and transformed data. While this doesn't > impact the training data, it can impact the prediction on test data > (although difference should be slight). > > Transformation of the response variable is quite another thing. RF needs > it just as much as others if the situation calls for it. > > Cheers, > Andy > > > > -Original Message- > > From: r-help-boun...@r-project.org > > [mailto:r-help-boun...@r-project.org] On Behalf Of gianni lavaredo > > Sent: Monday, December 05, 2011 1:41 PM > > To: r-help@r-project.org > > Subject: [R] explanation why RandomForest don't require a > > transformations (e.g. logarithmic) of variables > > > > Dear Researches, > > > > sorry for the easy and common question. I am trying to > > justify the idea of > > RandomForest don't require a transformations (e.g. logarithmic) of > > variables, comparing this non parametrics method with e.g. the linear > > regressions. In leteruature to study my phenomena i need to apply a > > logarithmic trasformation to describe my model, but i found RF don't > > required this approach. Some people could suggest me text or > > bibliography > > to study? > > > > thanks in advance > > > > Gianni > > > > [[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. > > > Notice: This e-mail message, together with any attach...{{dropped:16}} __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 to split a ID column in a data.frame and create a new ID column
Dear Researchers, I have a data.frame with 2 columns like this: mydf <- data.frame(value=c(1,2,3,4,5),ID=c("Area_1","Area_2","Area_3","Area_4","Area_5")) > mydf value ID 1 1 Area_1 2 2 Area_2 3 3 Area_3 4 4 Area_4 5 5 Area_5 I need to convert the *ID *in the following version > mydf value ID newID 1 1 Area_1 AreaSample1 2 2 Area_2 AreaSample2 3 3 Area_3 AreaSample3 4 4 Area_4 AreaSample4 5 5 Area_5 AreaSample5 some people know the right function to split ID and create a new column Thanks in advance Gianni [[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] plot a BARPLOT with sd deviation bar up and down
dear Researchers, i am looking for a function to plot a barplot for each mean value and the related standard deviation, and i can close my week. This is an example of my data set. really Thanks in advance for any help or suggestions Gianni My.mean <- data.frame(Mean=c(0.4108926,0.3949009,0.4520346, 0.4091665,0.4664066,0.3048296,0.4297226,0.4056383, 0.4127453,0.3568891,0.3933964,0.3892999,0.4052982, 0.377359,0.3831106,0.4248397,0.4403693,0.9389882)) My.SD <- data.frame(SD = c(0.3225084,0.3756248,0.3708947, 0.2899242,0.394396,0.4920173,0.2674820,0.3233239,0.2913170, 0.4542726,0.4031899,0.2893581,0.403938,0.3686252,0.4014624, 0.4105261,0.2811270,0.4088456,0.4889143,0.3949252,1.338804)) [[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] Different result with "kruskal.test" and post-hoc analysis with Nemenyi-Damico-Wolfe-Dunn test implemented in the help page for oneway_test in the coin package that uses multcomp
Dear Researchers, Sorry for this email but I am not a statistician, and for this I have this problem to understand. Thanks in Advance for help and suggestions. Gianni I have 21 classes (00, 01, 02, 04, ,020) with different length. I did a kruskal wall test in R with the following code kruskal.test(m.class.l, m.class.length.lf) Kruskal-Wallis rank sum test data: m.class.l and m.class.length.lf Kruskal-Wallis chi-squared = 34.2904, df = 20, p-value = 0.02423 the p-value result to be < of 0.05, and i read in different books and forum that when the medium are significative difference to apply a post-hoc test to check which classes are different. i find this code in a R help: http://stackoverflow.com/questions/2478272/kruskal-wallis-test-with-details-on-pairwise-comparisons class <- m.class.length.lf var <- m.class.l dft <- data.frame(class,var) NDWD <- oneway_test(var ~ class, data = dft, ytrafo = function(data) trafo(data, numeric_trafo = rank), xtrafo = function(data) trafo(data, factor_trafo = function(x) model.matrix(~x - 1) %*% t(contrMat(table(x), "Tukey"))), teststat = "max", distribution = approximate(B=1000)) ### global p-value print(pvalue(NDWD)) [1] 0.074 99 percent confidence interval: 0.05425181 0.09791886 ### sites (I = II) != (III = IV) at alpha = 0.01 (page 244) print(pvalue(NDWD, method = "single-step")) when I apply the code for my data for see details for each pair comparison, no class comparision is significative > prova <- pvalue(NDWD, method = "single-step") > prova.df <- data.frame(prova) > prova.df$sing <- Signif(prova.df$V1) > prova.df V1 sing 01 - 00 1.000 ns 010 - 00 1.000 ns 011 - 00 1.000 ns 012 - 00 1.000 ns 013 - 00 1.000 ns 014 - 00 1.000 ns 015 - 00 1.000 ns 016 - 00 1.000 ns 017 - 00 0.996 ns 018 - 00 0.870 ns 019 - 00 1.000 ns 02 - 00 1.000 ns 020 - 00 0.642 ns 03 - 00 1.000 ns 04 - 00 1.000 ns 05 - 00 1.000 ns 06 - 00 1.000 ns 07 - 00 1.000 ns 08 - 00 1.000 ns 09 - 00 1.000 ns 010 - 01 1.000 ns 011 - 01 1.000 ns 012 - 01 1.000 ns 013 - 01 1.000 ns 014 - 01 1.000 ns 015 - 01 1.000 ns 016 - 01 1.000 ns 017 - 01 0.981 ns 018 - 01 0.737 ns 019 - 01 1.000 ns 02 - 01 1.000 ns 020 - 01 0.441 ns 03 - 01 1.000 ns 04 - 01 1.000 ns 05 - 01 1.000 ns 06 - 01 1.000 ns 07 - 01 1.000 ns 08 - 01 1.000 ns 09 - 01 1.000 ns 011 - 010 1.000 ns 012 - 010 1.000 ns 013 - 010 1.000 ns 014 - 010 1.000 ns 015 - 010 1.000 ns 016 - 010 1.000 ns 017 - 010 0.990 ns 018 - 010 0.774 ns 019 - 010 1.000 ns 02 - 010 1.000 ns 020 - 010 0.490 ns 03 - 010 1.000 ns 04 - 010 1.000 ns 05 - 010 1.000 ns 06 - 010 1.000 ns 07 - 010 1.000 ns 08 - 010 1.000 ns 09 - 010 1.000 ns 012 - 011 0.963 ns 013 - 011 1.000 ns 014 - 011 1.000 ns 015 - 011 1.000 ns 016 - 011 0.993 ns 017 - 011 0.673 ns 018 - 011 0.226 ns 019 - 011 0.869 ns 02 - 011 1.000 ns 020 - 011 0.074 ns 03 - 011 1.000 ns 04 - 011 1.000 ns 05 - 011 1.000 ns 06 - 011 1.000 ns 07 - 011 1.000 ns 08 - 011 1.000 ns 09 - 011 1.000 ns 013 - 012 1.000 ns 014 - 012 1.000 ns 015 - 012 1.000 ns 016 - 012 1.000 ns 017 - 012 1.000 ns 018 - 012 1.000 ns 019 - 012 1.000 ns 02 - 012 1.000 ns 020 - 012 1.000 ns 03 - 012 1.000 ns 04 - 012 0.991 ns 05 - 012 1.000 ns 06 - 012 1.000 ns 07 - 012 1.000 ns 08 - 012 0.995 ns 09 - 012 1.000 ns 014 - 013 1.000 ns 015 - 013 1.000 ns 016 - 013 1.000 ns 017 - 013 0.976 ns 018 - 013 0.707 ns 019 - 013 1.000 ns 02 - 013 1.000 ns 020 - 013 0.395 ns 03 - 013 1.000 ns 04 - 013 1.000 ns 05 - 013 1.000 ns 06 - 013 1.000 ns 07 - 013 1.000 ns 08 - 013 1.000 ns 09 - 013 1.000 ns 015 - 014 1.000 ns 016 - 014 1.000 ns 017 - 014 0.995 ns 018 - 014 0.839 ns 019 - 014 1.000 ns 02 - 014 1.000 ns 020 - 014 0.560 ns 03 - 014 1.000 ns 04 - 014 1.000 ns 05 - 014 1.000 ns 06 - 014 1.000 ns 07 - 014 1.000 ns 08 - 014 1.000 ns 09 - 014 1.000 ns 016 - 015 1.000 ns 017 - 015 0.978 ns 018 - 015 0.704 ns 019 - 015 1.000 ns 02 - 015 1.000 ns 020 - 015 0.383 ns 03 - 015 1.000 ns 04 - 015 1.000 ns 05 - 015 1.000 ns 06 - 015 1.000 ns 07 - 015 1.000 ns 08 - 015 1.000 ns 09 - 015 1.000 ns 017 - 016 1.000 ns 018 - 016 1.000 ns 019 - 016 1.000 ns 02 - 016 1.000 ns 020 - 016 0.992 ns 03 - 016 1.000 ns 04 - 016 1.000 ns 05 - 016 1.000 ns 06 - 016 1.000 ns 07 - 016 1.000 ns 08 - 016 1.000 ns 09 - 016 1.000 ns 018 - 017 1.000 ns 019 - 017 1.000 ns 02 - 017 1.000 ns 020 - 017 1.000 ns 03 - 017 0.916 ns 04 - 017 0.785 ns 05 - 017 0.999 ns 06 - 017 0.994 ns 07 - 017 1.000 ns 08 - 017 0.852 ns 09 - 017 0.994 ns 019 - 018 1.000 ns 02 - 018 0.905 ns 020 - 018 1.000
[R] help to slip a file name using "strsplit" function
Dear Researchers, I have several files as this example: Myfile_MyArea1_sample1.txt i wish to split in "Myfile", "MyArea1", "sample1", and "txt", becasue i need to use "sample1" label. I try to use "strsplit" but I am able just to split as "Myfile_MyArea1_sample1" and "txt" OR "Myfile", "MyArea1", "sample1.txt" using strsplit(as.character("Myfile_MyArea1_sample1.txt"), ".", fixed = TRUE) strsplit(as.character("Myfile_MyArea1_sample1.txt"), "_", fixed = TRUE) Thanks in advance Gianni [[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 with Box plot
Dear researchers I wish to plot a box plot without the mean line (the black line) and the i wish a full line for the standard deviation This is an example mytest <- c(2.1,2.6,2.7,3.2,4.1,4.3,5.2,5.1,4.8,1.8,1.4,2.5,2.7,3.1,2.6,2.8) boxplot(mytest) really thanks Gianni [[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 boxplot to add mean, standard error and/or stadard deviation
Dear researchers I wish to plot a box plot without the mean line (the black line) and plot only the mean (red square). Futhermore, is it possible to add standard error and/or stadard deviation? This is an example mytest <- c(2.1,2.6,2.7,3.2,4.1,4.3,5.2,5.1,4.8,1.8,1.4,2.5,2.7,3.1,2.6,2.8) boxplot(mytest,lty = "solid") means <- mean(mytest,na.rm=TRUE) points(means, pch = 22, col = "red") thanks in advance Gianni [[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 to improve bwplot plot (lattice)
Dear Researchers, I wish to plot mean, standard deviation, and standard error and I am using bwplot(). I have the following problems and sorry if maybe there are simple questions: 1- use color black for standarddevuation line and add horizontal end bar (as the commun graphic in scientific papers) 2- change the label in x axes (es; A and B, respect 1 and 2) 3- add the standard error around the mean with a shape of a rectangular (high equal to standard error and length possibile to be visible in a graphic) this is an exmple to explain better A <- data.frame(class=rep(1,25)) B <- data.frame(class=rep(2,25)) class <- data.frame(rbind(A,B)) mydata <- data.frame(H=runif(50),class) bwplot(H~class, data=mydata, horizontal=FALSE, panel=panel.meansdplot,mean.col = "black") really thanks and good week-end Gianni [[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] plot with ylim with regural interval
Dear Researchers, sorry for the easy question but Is it possible to plot with an interval of 1 or .5 in a plot using ylim? Thanks gianni x = 0:10; y = 0:10; plot(x~y,ylim=c(0,10),las=1) [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] R citation for the 2012
Dear Reasearchers, I am writing a report and i need (and wish) cite R. somebody know the citation of R for the 2012? or the more actual? thanks in advance Gianni [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] R citation for the 2012
Sorry I don't know the citation inside the report is corret The Analysis was done using a script written in the statistical computing environment of R (R Development Core Team, 2011) is It correct the form "statistical computing environment" ? Gianni On Wed, Feb 15, 2012 at 5:38 PM, Sarah Goslee wrote: > Typing > citation() > at an R prompt will provide you with complete citation information for > the version of > R you are using. Same goes for packages, with citation("pkgname"). > > Sarah > > On Wed, Feb 15, 2012 at 11:29 AM, gianni lavaredo > wrote: > > Dear Reasearchers, > > > > I am writing a report and i need (and wish) cite R. somebody know the > > citation of R for the 2012? or the more actual? > > > > thanks in advance > > Gianni > > > > -- > Sarah Goslee > http://www.functionaldiversity.org > [[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] count how many row i have in a txt file in a directory
Dear Researchers, I have a large TXT (X,Y,MyValue) file in a directory and I wish to import row by row the txt in a loop to save only the data they are inside a buffer (using inside.owin of spatstat) and delete the rest. The first step before to create a loop row-by-row is to know how many rows there are in the txt file without load in R to save memory problem. some people know the specific function? thanks in advance for all suggestions Gianni [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Problem to resolve a step for reading a large TXT and split in several file
Dear Researchs, It's the first time I am trying to resolve this problem. I have a TXT file with 1408452 rows. I wish to split file-by-file where each file has 1,000,000 rows with the following procedure: # split in two file one with 1,000,000 of rows and one with 408,452 of rows file <- "09G001_72975_7575_25_4025.txt" fileNames <- strsplit(as.character(file), ".", fixed = TRUE) fileNames.temp.1 <- unique(as.vector(do.call("rbind", fileNames)[, 1])) con <- file(file, open = "r") # n is the number of row n <- 100 i <- 0 while (length(readLines(con, n=n)) > 0 ) { i <- i + 1 pv <- read.table(con,header=F,sep="\t", nrow=n) write.table(pv, file = paste(fileNames.temp.1,"_",i,".txt",sep = ""), sep = "\t") } close(con) when I use 1,000,000 I have in the directory only "09G001_72975_7575_25_4025_1.txt" (with 100 of rows) and not "09G001_72975_7575_25_4025_2.txt" (with 408,452). I din't understand where is my bug Furthermore when i wish for example split in 3 files (where n is 469484 = 1408452/3) i have this message: *Error in read.table(con, header = F, sep = "\t", nrow = n) : no lines available in input* Thanks for all help and sorry for the disturb [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] how disable the Error massage in read.table() " no lines available in input"
Dear Researchers, I am looking a way to disable the Error massage in read.table() as warn = TRUE in readLines(), when the lines are empty Error in read.table(con, header = F, sep = " ", nrow = n) : no lines available in input thanks for all suggestions Gianni [[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] question how to add Standard Deviation as "Whiskers" in a simple plot
Dear Researchers, sorry for this simple question. I have a point plot with mean values and i wish to plot line with Standard Deviation as "Whiskers". I calculate the mean+sd and mean-sd, but i can not figure out the way to add the line. mydata <- data.frame(mean=c(0.42,0.41,0.41,0.43,0.45,0.43,0.43,0.42,0.44,0.45,0.45,0.45,0.46,0.43,0.42,0.37,0.44,0.46,0.46,0.39,0.40), sdUP=c(0.58,0.56,0.55,0.57,0.61,0.55,0.57,0.59,0.61,0.60,0.57,0.60,0.62,0.57,0.59,0.56,0.57,0.61,0.61,0.56,0.54), sdDOWN=c(0.26,0.26,0.28,0.29,0.30,0.30,0.29,0.26,0.28,0.31,0.34,0.30,0.31,0.30,0.25,0.19,0.31,0.31,0.31,0.22,0.25)) plot(mydata$mean, type="o", ylab="mean", xlab="class") thanks in advance and sorry for any disturb Gianni [[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 how to add Standard Deviation as "Whiskers" in a simple plot
The function i am looking is a bars from the mean points of the plot in boxplot style. I tryed several forum but I have no clear the way to create these bars. Gianni On Mon, May 28, 2012 at 7:13 PM, David Winsemius wrote: > > On May 28, 2012, at 9:55 AM, gianni lavaredo wrote: > > Dear Researchers, >> >> sorry for this simple question. I have a point plot with mean values and >> i >> wish to plot line with Standard Deviation as "Whiskers". I calculate the >> mean+sd and mean-sd, but i can not figure out the way to add the line. >> >> mydata <- >> data.frame(mean=c(0.42,0.41,0.**41,0.43,0.45,0.43,0.43,0.42,0.** >> 44,0.45,0.45,0.45,0.46,0.43,0.**42,0.37,0.44,0.46,0.46,0.39,0.**40), >> >> sdUP=c(0.58,0.56,0.55,0.57,0.**61,0.55,0.57,0.59,0.61,0.60,0.** >> 57,0.60,0.62,0.57,0.59,0.56,0.**57,0.61,0.61,0.56,0.54), >> >> sdDOWN=c(0.26,0.26,0.28,0.29,**0.30,0.30,0.29,0.26,0.28,0.31,** >> 0.34,0.30,0.31,0.30,0.25,0.19,**0.31,0.31,0.31,0.22,0.25)) >> >> plot(mydata$mean, >> type="o", >> ylab="mean", >> xlab="class") >> >> thanks in advance and sorry for any disturb >> > > If you tried lines() usin either sdUP or sdDOWN you saw nothing because > the ylim was set by default using ony hte information in the mean vector. > > Try: > > > plot(mydata$mean, > +type="o", > +ylab="mean", > +xlab="class", ylim=range(c(mydata$sdUP, mydata$sdDOWN))) > > lines(mydata$sdDOWN, col="blue", lty=3) > > lines(mydata$sdUP, col="blue", lty=3) > > > -- > David. > [[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.