Re: [R] incomplete final line found warning
On Sun, 11 Dec 2011, David Winsemius wrote: On Dec 10, 2011, at 10:01 PM, Xiaobo Gu wrote: without following the posting guide in several respects and hence leaving us guessing Hi, I saved the following as a UTF-8 encoded file named amberutil.r BTW, it is hard to know how you know that ASCII is encoded as UTF-8, and on Windows (which from the file path it appears to be) it would not have worked had it been UTF-8 encoded. Let's hope this did not mean what Windows calls 'Unicode', that is UTF-16LE. as.factor.loop - function(df, cols){ if (!is.null(df) !is.null(cols) length(cols) 0) { for(col in cols) { df[[col]] - as.factor(df[[col]]) } } df } And got this warning message, source('D:/ambertuil.r') Warning message: In readLines(file) : incomplete final line found on 'D:/ambertuil.r' Can you help with this? Help with what? You got a warning. And it had information that should tell you how to edit the file if the warning bothers you. Also, we were not told the version of R. Updating (as requested by the posting guide prior to posting) would most likely remove the harmless warning (AFAIK it occurs only in 2.14.0 and not in R-patched) if this were an ASCII file. -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] incomplete final line found warning
On Sun, Dec 11, 2011 at 5:03 PM, Prof Brian Ripley rip...@stats.ox.ac.uk wrote: On Sun, 11 Dec 2011, David Winsemius wrote: On Dec 10, 2011, at 10:01 PM, Xiaobo Gu wrote: without following the posting guide in several respects and hence leaving us guessing Hi, I saved the following as a UTF-8 encoded file named amberutil.r BTW, it is hard to know how you know that ASCII is encoded as UTF-8, and on Windows (which from the file path it appears to be) it would not have worked had it been UTF-8 encoded. Let's hope this did not mean what Windows calls 'Unicode', that is UTF-16LE. I use RStudio to edit the source file, there is as save as encoding option, and I chose UTF-8 as.factor.loop - function(df, cols){ if (!is.null(df) !is.null(cols) length(cols) 0) { for(col in cols) { df[[col]] - as.factor(df[[col]]) } } df } And got this warning message, source('D:/ambertuil.r') Warning message: In readLines(file) : incomplete final line found on 'D:/ambertuil.r' Can you help with this? Help with what? You got a warning. And it had information that should tell you how to edit the file if the warning bothers you. Can you help finding the reason about this warning. Also, we were not told the version of R. Updating (as requested by the posting guide prior to posting) would most likely remove the harmless warning (AFAIK it occurs only in 2.14.0 and not in R-patched) if this were an ASCII file. I am using R 2.14.0 64 bit on Windows. -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UK Fax: +44 1865 272595 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Any packages or ways available providing un-/a-/non-symmetric centered distributions?
Dear R folks, I would like to use un-/a-/non-symmetric centered, i. e. the expectation vanishes/is zero, distributions. Searching for that is not that easy since all three words un-, a-, non-symmetric seem to be common. Do I need to create such a distribution myself by for example composing a distribution? I would use `sample()` and for integer-valued the following would be one example. The expectation vanishes too, 1/4 * (-4 - 1 + 2 + 3) = 0. x = sample(c(-4L, -1L, 2L, 3L)) y = sample(c(-4L, -1L, 2L, 3L), 1e7, replace=TRUE); mean(y) [1] 0.0006821 # replicate will produce a matrix, but in this case it should not matter for the mean x = replicate(2, sample(c(-4L, -1L, 2L, 3L), replace=TRUE)); x [,1] [,2] [1,] -12 [2,]3 -1 [3,]3 -1 [4,] -13 mean(x) [1] 0.875 # Passing the length explicitly to `sample()` as a work around. x = replicate(1e1, sample(c(-4L, -1L, 2L, 3L), 1, replace=TRUE)); x [1] -4 2 2 3 -4 -4 2 2 2 2 For continues (without 0) distributions I came up with the following. # -rexp(1, 2) has expectation -1. # sqrt(pi/2) * abs(rnorm(1)) # http://en.wikipedia.org/wiki/Half-normal_distribution # X ~ N(0, 1) ⇒ E(|X|) = √(2/π) x = sample(c(sqrt(pi/2) * abs(rnorm(1)), -rexp(1, 1)), 1); x [1] 0.2730692 Without passing `, 1` as the length to `sample` I would get two values. x = sample(c(sqrt(pi/2) * abs(rnorm(1)), -rexp(1, 1))); x [1] -0.4055199 1.8830235 Wanting to use `replicate()` therefore the explicit length for `sample()` should be passed here too. x = replicate(1e5, sample(c(sqrt(pi/2) * abs(rnorm(1)), -rexp(1, 1)), 1)); mean(x) [1] -0.001099127 Are there other sophisticated ways to create such distributions easily? Thanks, Paul signature.asc Description: This is a digitally signed message part __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] incomplete final line found warning
Xiaobo.Gu wrote Hi, I saved the following as a UTF-8 encoded file named amberutil.r as.factor.loop - function(df, cols){ if (!is.null(df) !is.null(cols) length(cols) 0) { for(col in cols) { df[[col]] - as.factor(df[[col]]) } } df } And got this warning message, source('D:/ambertuil.r') Warning message: In readLines(file) : incomplete final line found on 'D:/ambertuil.r' Can you help with this? A warning message such as this could not be clearer. It means that the last line of the file does not end with a newline sequence == the final line of the file is incomplete. In an editor go to the end of that line and press Enter or Return And save. Alternatively configure your editor to always terminate the last line of a file with a newline sequence. /Berend -- View this message in context: http://r.789695.n4.nabble.com/incomplete-final-line-found-warning-tp4181838p4182612.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] multiple comparison of interaction of ANCOVA
Hi there, The following data is obtained from a long-term experiments. mydata - read.table(textConnection( +y year Trt + 9.37 1993 A + 8.21 1995 A + 8.11 1999 A + 7.22 2007 A + 7.81 2010 A +10.85 1993 B +12.83 1995 B +13.21 1999 B +13.70 2007 B +15.15 2010 B + 5.69 1993 C + 5.76 1995 C + 6.39 1999 C + 5.73 2007 C + 5.55 2010 C), header = TRUE) closeAllConnections() The experiments is designed without replication, thus I have to use ANCOVA or linear mixed effect model to analyze the data. In the model, variable year is coded as a continuous variable, and Trt is factor variable. mydata.aov - aov(y~Trt*year, mydata) anova(mydata.aov) Analysis of Variance Table Response: y Df Sum Sq Mean Sq F valuePr(F) Trt2 140.106 70.053 197.9581 3.639e-08 *** year 1 0.610 0.610 1.7246 0.221600 Trt:year 2 8.804 4.402 12.4387 0.002567 ** Residuals 9 3.185 0.354 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 As you have seen, the interaction effect is significant. I hope to use TukeyHSD() or glht() to do multiple comparison on Trt:year. However, for variable year is not a factor, they all give error messages. I try to follow the demo(MMC.WoodEnergy) in HH package, as follwoing: library(HH) mca.1993 - mcalinfct(mydata.aov, Trt) non.zero - mca.1993[,5:6] != 0 mca.1993[,5:6][non.zero] - 1993 * sign(mca.1993[,5:6][non.zero]) summary(glht(mydata.aov, linfct=mca.1993)) Simultaneous Tests for General Linear Hypotheses Fit: aov(formula = y ~ Trt * year, data = mydata) Linear Hypotheses: Estimate Std. Error t value Pr(|t|) B - A == 0 2.8779 0.5801 4.961 0.00215 ** C - A == 0 -2.8845 0.5801 -4.972 0.00191 ** C - B == 0 -5.7624 0.5801 -9.933 0.001 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 (Adjusted p values reported -- single-step method) It can give comparison between levels of Trt within one year, e.g., 1993. Is it possible to do multiple comparison for the following pairs: A.1995 - A.1993 B.1995 - A.1993 C.1995 - A.1993 Any suggestions or comments will be really appreciated. Thanks in advance! Regards, Jinsong __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Titelpage pdf-manual, packaging
Am 10.12.2011 um 16:50 schrieb Uwe Ligges: On 09.12.2011 14:15, Johannes Radinger wrote: Hi, I am trying to get write my first own package. I followed the instructions in Creating R Packages: A Tutorial and Writing R Extensions. So far everything works really fine, the script works and even the man-pages don't show any problems during the check process. During check there is also the package-manual.pdf created. When I compare my manual to those at the CRAN rep, i looks different: I seems the other manuals have a kind of titlepage and even the document title is different (package 'xxx' vs. R documentation of ‘/Users/...’ etc.). Is that behavior simply a matter of being published at CRAN and mine is just a local package for myself? Or do I miss a man-page which is responsible for that first page and the title? Which R version are you using? I typically get package 'xxx' when running R CMD check but R documentation of ‘...’ when runnign R CMD Rd2pdf on some Rd files with R-release. Hi, I am using R 2.14.0, so I don't know what is the reason for that behavior... Best, Uwe Ligges Best regards, Johannes -- __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Bioconductor. MA plot for qPCR array
Dear all, Is there anyway too generate MA plot for 2 qPCR assays (an array of 2x 400). -- View this message in context: http://r.789695.n4.nabble.com/Bioconductor-MA-plot-for-qPCR-array-tp4182805p4182805.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] no non-missing arguments to max; returning -Inf
Hi, I just wanted to object of class xts, zoo and I got the error no non-missing arguments to max; returning -Inf After all the values were printed. Tried searching on this forum, but no luck. Anyone has any idea ? Thanks. George -- View this message in context: http://r.789695.n4.nabble.com/no-non-missing-arguments-to-max-returning-Inf-tp4182296p4182296.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] Titelpage pdf-manual, packaging
On 11/12/2011 12:34, Johannes Radinger wrote: Am 10.12.2011 um 16:50 schrieb Uwe Ligges: On 09.12.2011 14:15, Johannes Radinger wrote: Hi, I am trying to get write my first own package. I followed the instructions in Creating R Packages: A Tutorial and Writing R Extensions. So far everything works really fine, the script works and even the man-pages don't show any problems during the check process. During check there is also the package-manual.pdf created. When I compare my manual to those at the CRAN rep, i looks different: I seems the other manuals have a kind of titlepage and even the document title is different (package 'xxx' vs. R documentation of ‘/Users/...’ etc.). Is that behavior simply a matter of being published at CRAN and mine is just a local package for myself? Or do I miss a man-page which is responsible for that first page and the title? Which R version are you using? I typically get package 'xxx' when running R CMD check but R documentation of ‘...’ when runnign R CMD Rd2pdf on some Rd files with R-release. Hi, I am using R 2.14.0, so I don't know what is the reason for that behavior... And nor can we without the reproducible example requested at the foot of this and every R-help message. But one quick point: the way to make a package manual is R CMD Rd2pdf top package source directory and 'R CMD check' is not intended to be the way to produce that manual, merely a way to check that CRAN could produce it. Most likely there is diagnostic information in the R CMD check output, but we aren't being shown it Best, Uwe Ligges Best regards, Johannes -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] no non-missing arguments to max; returning -Inf
Yes - E-mail us reproducible code... Contact Details:--- Contact me: tal.gal...@gmail.com | 972-52-7275845 Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) | www.r-statistics.com (English) -- On Sun, Dec 11, 2011 at 10:25 AM, George Kumar grgkum...@gmail.com wrote: Hi, I just wanted to object of class xts, zoo and I got the error no non-missing arguments to max; returning -Inf After all the values were printed. Tried searching on this forum, but no luck. Anyone has any idea ? Thanks. George -- View this message in context: http://r.789695.n4.nabble.com/no-non-missing-arguments-to-max-returning-Inf-tp4182296p4182296.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.
Re: [R] Overlaying density plot on forest plot
At 07:16 10/12/2011, Frank Peter wrote: Dear R User, Please, I am new to R. I want to overlay density plot for predictive interval pooled result in meta-analysis. http://addictedtor.free.fr/graphiques/graphcode.php?graph=114 It is hard to be sure from your rather brief question but does the addcred parameter to forest.rma in package metafor do what you want? Regards Frank Peter Michael Dewey i...@aghmed.fsnet.co.uk http://www.aghmed.fsnet.co.uk/home.html __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] multiple comparison of interaction of ANCOVA
Inline below. -- Bert On Sun, Dec 11, 2011 at 4:15 AM, Jinsong Zhao jsz...@yeah.net wrote: Hi there, The following data is obtained from a long-term experiments. mydata - read.table(textConnection( + y year Trt + 9.37 1993 A + 8.21 1995 A + 8.11 1999 A + 7.22 2007 A + 7.81 2010 A + 10.85 1993 B + 12.83 1995 B + 13.21 1999 B + 13.70 2007 B + 15.15 2010 B + 5.69 1993 C + 5.76 1995 C + 6.39 1999 C + 5.73 2007 C + 5.55 2010 C), header = TRUE) closeAllConnections() The experiments is designed without replication, thus I have to use ANCOVA or linear mixed effect model to analyze the data. In the model, variable year is coded as a continuous variable, and Trt is factor variable. mydata.aov - aov(y~Trt*year, mydata) anova(mydata.aov) Analysis of Variance Table Response: y Df Sum Sq Mean Sq F value Pr(F) Trt 2 140.106 70.053 197.9581 3.639e-08 *** year 1 0.610 0.610 1.7246 0.221600 Trt:year 2 8.804 4.402 12.4387 0.002567 ** Residuals 9 3.185 0.354 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 As you have seen, the interaction effect is significant. I hope to use TukeyHSD() or glht() to do multiple comparison on Trt:year. However, for variable year is not a factor, they all give error messages. I try to follow the demo(MMC.WoodEnergy) in HH package, as follwoing: library(HH) mca.1993 - mcalinfct(mydata.aov, Trt) non.zero - mca.1993[,5:6] != 0 mca.1993[,5:6][non.zero] - 1993 * sign(mca.1993[,5:6][non.zero]) summary(glht(mydata.aov, linfct=mca.1993)) Simultaneous Tests for General Linear Hypotheses Fit: aov(formula = y ~ Trt * year, data = mydata) Linear Hypotheses: Estimate Std. Error t value Pr(|t|) B - A == 0 2.8779 0.5801 4.961 0.00215 ** C - A == 0 -2.8845 0.5801 -4.972 0.00191 ** C - B == 0 -5.7624 0.5801 -9.933 0.001 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 (Adjusted p values reported -- single-step method) It can give comparison between levels of Trt within one year, e.g., 1993. Is it possible to do multiple comparison for the following pairs: A.1995 - A.1993 B.1995 - A.1993 C.1995 - A.1993 Any suggestions or comments will be really appreciated. Thanks in advance! Graph the data sensibly to figure out what's going on. Statistical machinationsand anova tables with P values alone are not sufficient and can be opaque or misleading. If you do not know what sensibly is (or even if you do), consult: http://addictedtor.free.fr/graphiques/ Cheers, Bert Regards, Jinsong __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Bert Gunter Genentech Nonclinical Biostatistics Internal Contact Info: Phone: 467-7374 Website: http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] how to assign a value?
Hi there, I hope to modify values in a vector or matrix in the following code: for (i in 1:9) { assign(paste(a., i, sep = ), 1:i) get(paste(a., i, sep = ))[i] - i+50 } I get the following error message: Error in get(paste(a., i, sep = ))[i] - i + 50 : target of assignment expands to non-language object I have read the FAQ How can I turn a string into a variable?, however, I don't find a way to deal with: get(paste(a., i, sep = ))[i] - i+50 Any suggestions? Thanks in advance! Regards, Jinsong __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 assign a value?
On Dec 11, 2011, at 10:27 AM, Jinsong Zhao wrote: Hi there, I hope to modify values in a vector or matrix in the following code: for (i in 1:9) { assign(paste(a., i, sep = ), 1:i) get(paste(a., i, sep = ))[i] - i+50 } Just one matrix? Then you seem to have inappropriately borrowed using . as an indexing operation. In R that is just another character when used as an object name. a.1 is notgoing to evaulate to a[1]. Look at what you would have had after for (i in 1:9) { + assign(paste(a., i, sep = ), 1:i) + } ls() [1] aa.1 a.2 [4] a.3 a.4 a.5 [7] a.6 a.7 a.8 [10] a.9 a.1 [1] 1 a.2 [1] 1 2 Each of those assign() operations created a single vector of length i. I doubt that was what you intended, Better would be to describe your objects and your intentions, rather than expecting us to understand your goals by just looking at code that doesn't achieve thos goals. (There is no `get-` function which was the source of the error.) I get the following error message: Error in get(paste(a., i, sep = ))[i] - i + 50 : target of assignment expands to non-language object I have read the FAQ How can I turn a string into a variable?, however, I don't find a way to deal with: get(paste(a., i, sep = ))[i] - i+50 Any suggestions? Thanks in advance! Regards, Jinsong __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to assign a value?
You are basically in R Inferno Circle 8.1.40. http://www.burns-stat.com/pages/Tutor/R_inferno.pdf On 11/12/2011 15:27, Jinsong Zhao wrote: Hi there, I hope to modify values in a vector or matrix in the following code: for (i in 1:9) { assign(paste(a., i, sep = ), 1:i) get(paste(a., i, sep = ))[i] - i+50 } I get the following error message: Error in get(paste(a., i, sep = ))[i] - i + 50 : target of assignment expands to non-language object I have read the FAQ How can I turn a string into a variable?, however, I don't find a way to deal with: get(paste(a., i, sep = ))[i] - i+50 Any suggestions? Thanks in advance! Regards, Jinsong __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 Burns pbu...@pburns.seanet.com twitter: @portfolioprobe http://www.portfolioprobe.com/blog http://www.burns-stat.com (home of 'Some hints for the R beginner' and 'The R Inferno') __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] efficiently finding the integrals of a sequence of functions
Thanks,Hans! I agree that this is a good way of solving this problem. Here is another way. Instead of defining a vector of uni-dimensional functions and trying to integrating each component (a uni-dimensional function), we can do something below my.integrand-function(x,k) { return(f[x,k]) ## use matrix to represent the sequence of uni-dimensional functions } my.integral-function(k) ## k=1,...,5000 denotes the function labels { return(integrate(my.integrand,lower=...,upper=...,k)$value) } When calculating the integrals, just perform sapply(1:5000, my.integral) This is a way of avoiding loops but the computing time needs to be carefully examined. Jeff -- View this message in context: http://r.789695.n4.nabble.com/efficiently-finding-the-integrals-of-a-sequence-of-functions-tp4179452p4183323.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] some CRAN mirrors not accessible
Dear all, Since some weeks, look like the following CRAN mirrors are no longer accessible for package update: http://cran.univ-lyon1.fr http://mirror.ibcp.fr/pub/CRAN/ update.packages(ask='graphics',checkBuilt=TRUE) Warning: unable to access index for repository http://cran.univ-lyon1.fr/bin/windows/contrib/2.14 update.packages(ask='graphics',checkBuilt=TRUE) Warning: unable to access index for repository http://mirror.ibcp.fr/pub/CRAN/bin/windows/contrib/2.14 Furthermore, attempting to connect via a web browser gives However, I don't know how to get in touch with the webmasters in charge. Any idea about hos to signal the trouble ? PG __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] multiple comparison of interaction of ANCOVA
Thank you for you use of HH. I think the right graph for this data is the much simpler ancova function library(HH) ancova(y ~ year * Trt, data=mydata) where we see that the three treatments have totally different slopes. The WoodEnergy example doesn't apply here. The WoodEnergy example illustrates a way of finding differences among treatments for a fixed value of the covariate when the slopes are similar. Rich On Sun, Dec 11, 2011 at 7:15 AM, Jinsong Zhao jsz...@yeah.net wrote: Hi there, The following data is obtained from a long-term experiments. mydata - read.table(textConnection( +y year Trt + 9.37 1993 A + 8.21 1995 A + 8.11 1999 A + 7.22 2007 A + 7.81 2010 A +10.85 1993 B +12.83 1995 B +13.21 1999 B +13.70 2007 B +15.15 2010 B + 5.69 1993 C + 5.76 1995 C + 6.39 1999 C + 5.73 2007 C + 5.55 2010 C), header = TRUE) closeAllConnections() The experiments is designed without replication, thus I have to use ANCOVA or linear mixed effect model to analyze the data. In the model, variable year is coded as a continuous variable, and Trt is factor variable. mydata.aov - aov(y~Trt*year, mydata) anova(mydata.aov) Analysis of Variance Table Response: y Df Sum Sq Mean Sq F valuePr(F) Trt2 140.106 70.053 197.9581 3.639e-08 *** year 1 0.610 0.610 1.7246 0.221600 Trt:year 2 8.804 4.402 12.4387 0.002567 ** Residuals 9 3.185 0.354 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 As you have seen, the interaction effect is significant. I hope to use TukeyHSD() or glht() to do multiple comparison on Trt:year. However, for variable year is not a factor, they all give error messages. I try to follow the demo(MMC.WoodEnergy) in HH package, as follwoing: library(HH) mca.1993 - mcalinfct(mydata.aov, Trt) non.zero - mca.1993[,5:6] != 0 mca.1993[,5:6][non.zero] - 1993 * sign(mca.1993[,5:6][non.zero]) summary(glht(mydata.aov, linfct=mca.1993)) Simultaneous Tests for General Linear Hypotheses Fit: aov(formula = y ~ Trt * year, data = mydata) Linear Hypotheses: Estimate Std. Error t value Pr(|t|) B - A == 0 2.8779 0.5801 4.961 0.00215 ** C - A == 0 -2.8845 0.5801 -4.972 0.00191 ** C - B == 0 -5.7624 0.5801 -9.933 0.001 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 (Adjusted p values reported -- single-step method) It can give comparison between levels of Trt within one year, e.g., 1993. Is it possible to do multiple comparison for the following pairs: A.1995 - A.1993 B.1995 - A.1993 C.1995 - A.1993 Any suggestions or comments will be really appreciated. Thanks in advance! Regards, Jinsong __** R-help@r-project.org mailing list https://stat.ethz.ch/mailman/**listinfo/r-helphttps://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/** posting-guide.html http://www.r-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] some CRAN mirrors not accessible
First, see the status links in the first para of http://cran.r-project.org/mirrors.html . So it seems that the issue is not that the mirror is not accessible, but the part you are looking for is not current/available. Second, the mirror list is http://cran.r-project.org/CRAN_mirrors.csv and that lists the maintainers. It is also in your R distribution, in directory 'doc' (but that is as old as your distribution, and mirrors do change). On 11/12/2011 17:43, Patrick Giraudoux wrote: Dear all, Since some weeks, look like the following CRAN mirrors are no longer accessible for package update: http://cran.univ-lyon1.fr http://mirror.ibcp.fr/pub/CRAN/ update.packages(ask='graphics',checkBuilt=TRUE) Warning: unable to access index for repository http://cran.univ-lyon1.fr/bin/windows/contrib/2.14 update.packages(ask='graphics',checkBuilt=TRUE) Warning: unable to access index for repository http://mirror.ibcp.fr/pub/CRAN/bin/windows/contrib/2.14 Furthermore, attempting to connect via a web browser gives However, I don't know how to get in touch with the webmasters in charge. Any idea about hos to signal the trouble ? PG __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] some CRAN mirrors not accessible
Le 11/12/2011 18:57, Prof Brian Ripley a écrit : First, see the status links in the first para of http://cran.r-project.org/mirrors.html . So it seems that the issue is not that the mirror is not accessible, but the part you are looking for is not current/available. yes indeed. For one R for Windows contrib deny access, and for the other has no 2.14/ folder Second, the mirror list is http://cran.r-project.org/CRAN_mirrors.csv and that lists the maintainers. It is also in your R distribution, in directory 'doc' (but that is as old as your distribution, and mirrors do change). Got it ! Thanks. Also fubbling in the chooseCRANmirror() doc I also discovered that getCRANmirrors() returns a data.frame with the maintainer address; looking into the function it either straight reads http://cran.r-project.org/CRAN_mirrors.csv, and if no connection reads de csv file stored in the local directory 'doc'. Now the info about accessibility/availability has been conveyed to the maintainers. Thanks again, On 11/12/2011 17:43, Patrick Giraudoux wrote: Dear all, Since some weeks, look like the following CRAN mirrors are no longer accessible for package update: http://cran.univ-lyon1.fr http://mirror.ibcp.fr/pub/CRAN/ update.packages(ask='graphics',checkBuilt=TRUE) Warning: unable to access index for repository http://cran.univ-lyon1.fr/bin/windows/contrib/2.14 update.packages(ask='graphics',checkBuilt=TRUE) Warning: unable to access index for repository http://mirror.ibcp.fr/pub/CRAN/bin/windows/contrib/2.14 Furthermore, attempting to connect via a web browser gives However, I don't know how to get in touch with the webmasters in charge. Any idea about hos to signal the trouble ? PG __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to assign a value?
I find that get() and assign() are awkward to use and that the syntax is easier if you put your objects into a list or environment. To me, it also makes it clearer what the code is doing and keeps the output of objects() shorter and easier to manage. E.g., nResults - 9 results - vector(list, nResults) # or results - new.env() for(i in 1:nResults) { resultName - paste(a., i, sep=) results[[resultName]] - 1:i results[[resultName]][i] - i+50 } Bill Dunlap Spotfire, TIBCO Software wdunlap tibco.com -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Jinsong Zhao Sent: Sunday, December 11, 2011 7:28 AM To: r-help@r-project.org Subject: [R] how to assign a value? Hi there, I hope to modify values in a vector or matrix in the following code: for (i in 1:9) { assign(paste(a., i, sep = ), 1:i) get(paste(a., i, sep = ))[i] - i+50 } I get the following error message: Error in get(paste(a., i, sep = ))[i] - i + 50 : target of assignment expands to non-language object I have read the FAQ How can I turn a string into a variable?, however, I don't find a way to deal with: get(paste(a., i, sep = ))[i] - i+50 Any suggestions? Thanks in advance! Regards, Jinsong __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] nls start values
I'm using nls to fit periodic gene-expression data to sine waves. I need to set the upper and lower boundaries, because I do not want any negative phase and amplitude solutions. This means that I have to use the port algorithm. The problem is, that depending on what start value I choose for phase, the fit works for some cases, but not for others. In the example below, the fit works using phase=pi, but not using phase=0. But there are many examples which fit just fine using 0. Is there a comparable alternative to nls that is not so extremely influenced by the start values? # Data for example fit lowervals - list(phase=0, amp=0) uppervals - list(phase=2*pi, amp=2) afreq - 1 / (24 / 2 / pi) gene_expression - c(1.551383, 1.671742, 1.549499, 1.694480, 1.632436, 1.471568, 1.623381, 1.579361, 1.809394, 1.753223, 1.685918, 1.754968, 1.963069, 1.820690, 1.985159, 2.205064, 2.160308, 2.120189, 2.194758, 2.165993, 2.189981, 2.098671, 2.122207, 2.012621, 1.963610, 1.884184, 1.955160, 1.801175, 1.829686, 1.773260, 1.588768, 1.563774, 1.559192) tpoints - c(0,0,0,2,2,2,4,4,4,6,6,6,8,8,8,12,12,12,14,14,14,16,16,16,18,18,18,20,20,20,24,24,24) shift=mean(gene_expression) # y-axis (expression) shift # Perfect fit startvals - list(phase=pi, amp=0.5) sine_nls - nls(gene_expression ~ sin(tpoints * afreq + phase) * amp + shift, start=startvals, algorithm=port, lower=lowervals, upper=uppervals) # Convergence failure startvals - list(phase=0, amp=0.5) sine_nls - nls(gene_expression ~ sin(tpoints * afreq + phase) * amp + shift, start=startvals, algorithm=port, lower=lowervals, upper=uppervals) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Boston/Cambridge -- Statistical Programming Language Technology Breakthroughs
If you are a statistician or researcher working in Boston/Cambridge, and you have a strong interest in breakthroughs in statistical programming language technology, contact me. Robert __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Kaplan-Meier survfit problem
2011/12/11 Esteban Cervetto estebancs...@gmail.com I am working with uncensored data. I have duration of workers compensation. Then I have for each the number of days thet it doesn't work. This sample is not censored at right because I query only work accidents with date of return of work (saned) That is because I have only one vector: the number of days that the worker doesn't work. Reading works that uses this library, noticed that it need's a vector to mark the type of termination. That is because I did a formula Y(x) = 1 as.numeric(T.**201110))~1 Your words were inspired me to do this: I solvet it putting to T.201110 this proper vector and a vector of ones. T.201110$censor - apply(T.201110,1,function(row) 1)##is there a best method to do that? it takes much time I believe that this vector is superfluous, because the result of the formula is ever 1 km1-survfit(Surv(T.201110$dias,T.201110$censor)~1) 2011/12/10 David Winsemius dwinsem...@comcast.net On Dec 10, 2011, at 6:39 PM, capitantyler wrote: done it, again, i have the next problem my traduction: The object (list) cannot be corced as double Original: *km1 - survfit(Surv(as.numeric(T.**201110))~1)* Error en Surv(as.numeric(T.201110)) : el objeto (list) no puede ser coercionado a 'double' I do not read that language, but I am surprised to see a single vector being used as an argument to Surv(). When I use Surv(. , .) it is with two vectors, an interval and a censor variable. note that need it convert to numeric class, otherwise: *km1 - survfit(Surv((T.201110))~1)* Error en Surv((T.201110)) : Time variable is not numeric -- View this message in context: http://r.789695.n4.nabble.com/** Kaplan-Meier-survfit-problem-**tp2015369p4181476.htmlhttp://r.789695.n4.nabble.com/Kaplan-Meier-survfit-problem-tp2015369p4181476.html Context? Yes. We do want context, but we don't want to go to no steenking Nabble. Sent from the R help mailing list archive at Nabble.com. Arrrgh. The rhelp mailing list is NOT at Nabble. Nabbel si a commercial mirror of the real thing. And it isn't really an archive. either, since they start discarding posts after a year or two. -- David Winsemius, MD West Hartford, CT [[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] using dcast to reshape a DF from long to wide with multiple measured variables per obs
I have data in the following format: person- c(1,1,1,1,2,2,2,2,2,3,3,3,3,3,3) v2- c(2011-01-01, 2011-02-01, 2011-03-01, 2011-04-01, 2011-01-01, 2011-02-01, 2011-03-01, 2011-04-01, 2011-05-01, 2011-01-01, 2011-02-01, 2011-03-01, 2011-04-01, 2011-05-01, 2011-06-01) v3 - rep(30, 15) DF -data.frame(person, v2, v3) DF$fillno - with(DF, ave(person, person, FUN = seq)) DF$v2 - as.Date(DF$v2) DF$v2 - as.numeric(DF$v2) It represents data in long format that needs to be reshaped by person, with each fill number+v2 and fill number+v3 as a column, with missing values in cases where a person does not have a fill number. Essentially, I would like the DF to eventually assume the following strucutre: DFnew - structure(list(person = 1:3, fillno1_v2 = c(14975L, 14975L, 14975L ), fillno1_v3 = c(30L, 30L, 30L), fillno2_v2 = c(15006L, 15006L, 15006L), fillno2_v3 = c(30L, 30L, 30L), fillno3_v2 = c(15034L, 15034L, 15034L), fillno3_v3 = c(30L, 30L, 30L), fillno4_v2 = c(15065L, 15065L, 15065L), fillno4_v3 = c(30L, 30L, 30L), fillno5_v2 = c(NA, 15095L, 15095L), fillno5_v3 = c(NA, 30L, 30L), fillno6_v2 = c(NA, NA, 15126L), fillno6_v3 = c(NA, NA, 30L)), .Names = c(person, fillno1_v2, fillno1_v3, fillno2_v2, fillno2_v3, fillno3_v2, fillno3_v3, fillno4_v2, fillno4_v3, fillno5_v2, fillno5_v3, fillno6_v2, fillno6_v3), class = data.frame, row.names = c(NA, -3L)) I have tried melt and dcast, but can only get reshape to organize by one measured variable at a time (and not two, as I would like). For example, the following organizes my data by fill-number and v2: dcast(DF, person ~ fillno, value.var=c(v2)) likewise the following works for v3 dcast(DF, person ~ fillno, value.var=c(v3)) But the following does not work: dcast(DF, person ~ fillno, value.var=c(v2, v3)) and produces the following error: Error in .subset2(x, i, exact = exact) : subscript out of bounds Thank you in advance for any help you may be able to provide! Chris [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Labelling R plots using Greek letters
Dear R Users, Please I have the following query. I want to label one of the axes of my graph with the follwing latex expression- \beta^{\prime}x, i.e I have the transpose of beta. How do I go abt this. My second query is similar but it has to do with the conditioning symbol. The other axis is to be labelled as Y|X + \beta^{\prime}x. Thank you very much as I await assistance on this. John [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Kaplan-Meier survfit problem
On Dec 11, 2011, at 5:48 PM, Esteban Cervetto wrote: 2011/12/11 Esteban Cervetto estebancs...@gmail.com I am working with uncensored data. I have duration of workers compensation. Then I have for each the number of days thet it doesn't work. This sample is not censored at right because I query only work accidents with date of return of work (saned) That is because I have only one vector: the number of days that the worker doesn't work. Reading works that uses this library, noticed that it need's a vector to mark the type of termination. That is because I did a formula Y(x) = 1 as.numeric(T.**201110))~1 Your words were inspired me to do this: I solvet it putting to T. 201110 this proper vector and a vector of ones. T.201110$censor - apply(T.201110,1,function(row) 1)##is there a best method to do that? it takes much time I believe that this vector is superfluous, because the result of the formula is ever 1 You may be correct. The help page for Surv says Although unusual, the event indicator can be omitted, in which case all subjects are assumed to have an event. km1-survfit(Surv(T.201110$dias,T.201110$censor)~1) 2011/12/10 David Winsemius dwinsem...@comcast.net On Dec 10, 2011, at 6:39 PM, capitantyler wrote: done it, again, i have the next problem my traduction: The object (list) cannot be corced as double Original: *km1 - survfit(Surv(as.numeric(T.**201110))~1)* If you are adding the ** for emphasis, it is certainly confusing my understanding of what your original code was, which i'm now wondering you ever provided. Perhaps it was the lack of a data argument. Hard to tell. This works: fit - survfit(Surv(time) ~ 1, data = aml[aml$status==1, ]) plot(fit) -- David Error en Surv(as.numeric(T.201110)) : el objeto (list) no puede ser coercionado a 'double' I do not read that language, but I am surprised to see a single vector being used as an argument to Surv(). When I use Surv(. , .) it is with two vectors, an interval and a censor variable. note that need it convert to numeric class, otherwise: *km1 - survfit(Surv((T.201110))~1)* Error en Surv((T.201110)) : Time variable is not numeric -- View this message in context: http://r.789695.n4.nabble.com/** Kaplan-Meier-survfit-problem-**tp2015369p4181476.htmlhttp://r.789695.n4.nabble.com/Kaplan-Meier-survfit-problem-tp2015369p4181476.html Context? Yes. We do want context, but we don't want to go to no steenking Nabble. Sent from the R help mailing list archive at Nabble.com. Arrrgh. The rhelp mailing list is NOT at Nabble. Nabbel si a commercial mirror of the real thing. And it isn't really an archive. either, since they start discarding posts after a year or two. -- David Winsemius, MD West Hartford, CT [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Labelling R plots using Greek letters
On Dec 11, 2011, at 6:38 PM, john james wrote: Dear R Users, Please I have the following query. I want to label one of the axes of my graph with the follwing latex expression- \beta^{\prime}x, i.e I have the transpose of beta. How do I go abt this. My second query is similar but it has to do with the conditioning symbol. The other axis is to be labelled as Y|X + \beta^{\prime}x. ?plotmath plot(1,1, xlab=expression(beta^`*X), ylab=expression(Y*|*X +beta^`*X) ) (Unfortunately the prime on the y lab is noe on the printed page in my device ( mac). so this hack seems to work better) plot(1,1, xlab=expression(beta^'*X), ylab=expression(over(phantom(), Y*|*X+beta^'*X) )) Thank you very much as I await assistance on this. John [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] limit ranges in hexbin
Hello everybody, I hope you can give me some help with limiting the ranges in x, y, and z for a hexbin plot. All I have found on the net is an unanswered message to this list from last year, so I hope my problem is not too stupid. I would like to plot some data using hexagonal binning. Currently, I am doing the following: library(hexbin) data-scan(input.dat, what = list('numeric', 'numeric')) x-unlist(data[1]) y-unlist(data[2]) hbin-hexbin(x, y) plot(hbin) That does get me a hexgonal bin plot, so far, so good. My problem now is twofold: 1) Some of my data points have really high values in x and y, and I would like to exclude them. First idea was to set xbnds and ybnds parameters for hexbin like this. hbin-hexbin(x, y, xbins = 50, xbnds = c(0, 5000)) However, hexbin does not like it: Error in hexbin(x, y, xbins = 50, xbnds = c(0, 5000)) : 'xbnds' must encompass range(x). Is there some other way to enforce limits in x and y? Otherwise, I could, of course, prefilter my data. 2) The second problem is limiting the count values. My counts are strongly dominated by one bin. So I end up with that bin being black and some others very light gray. I would like to set the black value to some lower number. Of course, some bins will be saturated then, but for that I will get better contrast in the range I am interested in. Thanks a lot for your help, Lutz __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 assign a value?
On 2011-12-12 0:00, David Winsemius wrote: On Dec 11, 2011, at 10:27 AM, Jinsong Zhao wrote: Hi there, I hope to modify values in a vector or matrix in the following code: for (i in 1:9) { assign(paste(a., i, sep = ), 1:i) get(paste(a., i, sep = ))[i] - i+50 } Just one matrix? Then you seem to have inappropriately borrowed using . as an indexing operation. In R that is just another character when used as an object name. a.1 is notgoing to evaulate to a[1]. Look at what you would have had after for (i in 1:9) { + assign(paste(a., i, sep = ), 1:i) + } ls() [1] a a.1 a.2 [4] a.3 a.4 a.5 [7] a.6 a.7 a.8 [10] a.9 a.1 [1] 1 a.2 [1] 1 2 Each of those assign() operations created a single vector of length i. I doubt that was what you intended, yes, it was what I intended. Better would be to describe your objects and your intentions, rather than expecting us to understand your goals by just looking at code that doesn't achieve thos goals. (There is no `get-` function which was the source of the error.) The question is why get(paste(a., i, sep = ))[i] - i+50 give the following error message: Error in get(paste(a., i, sep = ))[i] - i + 50 : target of assignment expands to non-language object The a.1 to a.9 was created in the previous step. if only get(paste(a., i, sep = ))[i] can give correct output. Why I cannot assign values to it? Regards, Jinsong __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Lagged values problem requiring short solution time
Hello, I am hoping someone can help tackle the problem below, for which I require a fast solution. It feels like there should be an elegant approach, but I am drawing blanks. Take a vector 'x' with random values 0: x = runif(10,1,5) Assume some reasonably small positive value 'delta': delta = 0.75 The task is to find a solution vector 'y' of same length as 'x' such that: The absolute difference between y[i] and y[i-1] is = 'delta' and y[i] = x[i] and that the sum of y[i] - x[i] be as small as possible -- i.e. minimize sum(y-x). The real-world application that (loosely) inspires this problem is the case of thermal power plants that face limits ('delta') in the speed with which they can ramp output up (or down) in response to changing demand. The period-to-period difference in output cannot exceed the absolute value of 'delta'. The other constraints I've imposed are specific to my application, but also provide a more neatly defined problem. A real-world problem would not have random starting values for 'x', but I figure the random values will present a particular difficulty in terms of solution time. SPEED IS CRITICAL here, as this example must handle 'x' with length=10,000 in practice and is located within an optimization routine that requires it be iterated over different 'x' vectors many times. My Neanderthal-ish solution (below) may or may not give the theoretically optimal solution, but, regardless, is too slow when 'x' becomes lengthy due to its reliance on loops. Hope you can help! x = runif(10,1,5) delta = 0.75 chg = diff(c(x,x[1])) y = x while (any(abs(chg)delta)) { temp = sign(chg)*chg - delta temp1=temp temp1[chg=(-delta)] = 0 temp1 = c(temp1[length(temp1)],temp1[-length(temp1)]) temp2 = temp temp2[chg=delta] = 0 y = y+temp1+temp2 chg = diff(c(y,y[1])) } #Solution vector: y Thank you, Kevin Kevin Ummel CARMA Project Manager Center for Global Development [[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 about fitting seasonal ARIMA in R?
Hi all, I just couldn't find a R function which can fit multiple seasonal patters... i.e. in the following code: *arima(x = data, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), ... *** * there can be only one period, am I right? What if the data seem to have three different seasonality cycles, 5, 12, 21? Thanks a lot! * [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to assign a value?
On 2011-12-12 1:16, Patrick Burns wrote: You are basically in R Inferno Circle 8.1.40. http://www.burns-stat.com/pages/Tutor/R_inferno.pdf Thanks. R_inferno is a good material for me. In this document, there is several sections titled string not the name. I try try to change get(paste(a., i, sep = ))[i] - i+50 to assign(get(paste(a., i, sep = ))[i], i+50) however, error message: Error in assign(get(paste(a., i, sep = ))[i], i + 50) : invalid first argument I don't know why... Regards, Jinsong On 11/12/2011 15:27, Jinsong Zhao wrote: Hi there, I hope to modify values in a vector or matrix in the following code: for (i in 1:9) { assign(paste(a., i, sep = ), 1:i) get(paste(a., i, sep = ))[i] - i+50 } I get the following error message: Error in get(paste(a., i, sep = ))[i] - i + 50 : target of assignment expands to non-language object I have read the FAQ How can I turn a string into a variable?, however, I don't find a way to deal with: get(paste(a., i, sep = ))[i] - i+50 Any suggestions? Thanks in advance! Regards, Jinsong __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 about fitting seasonal ARIMA in R?
sorry for the re-post... On Sun, Dec 11, 2011 at 8:22 PM, Michael comtech@gmail.com wrote: Hi all, I just couldn't find a R function which can fit multiple seasonal patters... i.e. in the following code: *arima(x = data, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), ... *** * there can be only one period, am I right? What if the data seem to have three different seasonality cycles, 5, 12, 21? Thanks a lot! * [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how to assign a value?
On Dec 11, 2011, at 9:07 PM, Jinsong Zhao wrote: On 2011-12-12 0:00, David Winsemius wrote: On Dec 11, 2011, at 10:27 AM, Jinsong Zhao wrote: Hi there, I hope to modify values in a vector or matrix in the following code: for (i in 1:9) { assign(paste(a., i, sep = ), 1:i) get(paste(a., i, sep = ))[i] - i+50 } Just one matrix? Then you seem to have inappropriately borrowed using . as an indexing operation. In R that is just another character when used as an object name. a.1 is notgoing to evaulate to a[1]. Look at what you would have had after for (i in 1:9) { + assign(paste(a., i, sep = ), 1:i) + } ls() [1] a a.1 a.2 [4] a.3 a.4 a.5 [7] a.6 a.7 a.8 [10] a.9 a.1 [1] 1 a.2 [1] 1 2 Each of those assign() operations created a single vector of length i. I doubt that was what you intended, yes, it was what I intended. Then you are free to continue banging your head against a wall. Better would be to describe your objects and your intentions, rather than expecting us to understand your goals by just looking at code that doesn't achieve thos goals. (There is no `get-` function which was the source of the error.) The question is why get(paste(a., i, sep = ))[i] - i+50 give the following error message: What part of THERE IS NO get- function (much less a `get[-` function) don't you understand? Error in get(paste(a., i, sep = ))[i] - i + 50 : target of assignment expands to non-language object The a.1 to a.9 was created in the previous step. if only get(paste(a., i, sep = ))[i] can give correct output. Right. They are there and can even be indexed: get(paste(a, 9, sep=.))[9] [1] 9 You could assign the value of get(paste(a, 9, sep=.)) to an intermediate object, which you could then reference using [ and then use `assign` to push that object's value back to an object named a. 1, , a.2, etc. Very clumsy and not an idiom that people want to promote. x - get(paste(a, 9, sep=.)) x[9] - x[9]+50 assign(paste(a, 9, sep=.), x) a.9 [1] 1 2 3 4 5 6 7 8 59 Why I cannot assign values to it? Using get, you mean? Because that is not the way R is designed. get() returns a value. `assign` is used... wait for it ... assignment. get(paste(a, 1, sep=.)) [1] 1 Not a.1 but rather a.1's value. You cannot assign something else to the number 1. You are free to complain about the fact that R is is not languageX as much as you like, but it won't create new capabilities for functions. You've been given advice about how to get to the goal you desire by both Dunlap and Burns. The counter-question is why you have such trouble accepting advice. Regards, Jinsong David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Labelling R plots using Greek letters
A short answer is you can consider the tikzDevice package; then a long (really really long) answer is this: http://yihui.github.com/knitr/demo/graphics/ Regards, Yihui -- Yihui Xie xieyi...@gmail.com Phone: 515-294-2465 Web: http://yihui.name Department of Statistics, Iowa State University 2215 Snedecor Hall, Ames, IA On Sun, Dec 11, 2011 at 7:17 PM, David Winsemius dwinsem...@comcast.net wrote: On Dec 11, 2011, at 6:38 PM, john james wrote: Dear R Users, Please I have the following query. I want to label one of the axes of my graph with the follwing latex expression- \beta^{\prime}x, i.e I have the transpose of beta. How do I go abt this. My second query is similar but it has to do with the conditioning symbol. The other axis is to be labelled as Y|X + \beta^{\prime}x. ?plotmath plot(1,1, xlab=expression(beta^`*X), ylab=expression(Y*|*X+beta^`*X) ) (Unfortunately the prime on the y lab is noe on the printed page in my device ( mac). so this hack seems to work better) plot(1,1, xlab=expression(beta^'*X), ylab=expression(over(phantom(), Y*|*X+beta^'*X) )) Thank you very much as I await assistance on this. John [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. David Winsemius, MD West Hartford, CT __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 assign a value?
On 2011-12-12 10:58, David Winsemius wrote: On Dec 11, 2011, at 9:07 PM, Jinsong Zhao wrote: On 2011-12-12 0:00, David Winsemius wrote: On Dec 11, 2011, at 10:27 AM, Jinsong Zhao wrote: Hi there, I hope to modify values in a vector or matrix in the following code: for (i in 1:9) { assign(paste(a., i, sep = ), 1:i) get(paste(a., i, sep = ))[i] - i+50 } Just one matrix? Then you seem to have inappropriately borrowed using . as an indexing operation. In R that is just another character when used as an object name. a.1 is notgoing to evaulate to a[1]. Look at what you would have had after for (i in 1:9) { + assign(paste(a., i, sep = ), 1:i) + } ls() [1] a a.1 a.2 [4] a.3 a.4 a.5 [7] a.6 a.7 a.8 [10] a.9 a.1 [1] 1 a.2 [1] 1 2 Each of those assign() operations created a single vector of length i. I doubt that was what you intended, yes, it was what I intended. Then you are free to continue banging your head against a wall. Better would be to describe your objects and your intentions, rather than expecting us to understand your goals by just looking at code that doesn't achieve thos goals. (There is no `get-` function which was the source of the error.) The question is why get(paste(a., i, sep = ))[i] - i+50 give the following error message: What part of THERE IS NO get- function (much less a `get[-` function) don't you understand? Sorry, I didn't understand it in the previous post. Now, it seems clear... Error in get(paste(a., i, sep = ))[i] - i + 50 : target of assignment expands to non-language object The a.1 to a.9 was created in the previous step. if only get(paste(a., i, sep = ))[i] can give correct output. Right. They are there and can even be indexed: get(paste(a, 9, sep=.))[9] [1] 9 You could assign the value of get(paste(a, 9, sep=.)) to an intermediate object, which you could then reference using [ and then use `assign` to push that object's value back to an object named a.1, , a.2, etc. Very clumsy and not an idiom that people want to promote. x - get(paste(a, 9, sep=.)) x[9] - x[9]+50 assign(paste(a, 9, sep=.), x) a.9 [1] 1 2 3 4 5 6 7 8 59 Yes, the intermediate object could be used to archive my goal: for (i in 1:9) { assign(paste(a, i, sep = .), 1:i) x - get(paste(a, i, sep = .)) x[[i]] - x[[i]] + 50 assign(paste(a, i, sep = .), x) } Why I cannot assign values to it? Using get, you mean? Because that is not the way R is designed. get() returns a value. `assign` is used... wait for it ... assignment. get(paste(a, 1, sep=.)) [1] 1 Not a.1 but rather a.1's value. You cannot assign something else to the number 1. You are free to complain about the fact that R is is not languageX as much as you like, but it won't create new capabilities for functions. You've been given advice about how to get to the goal you desire by both Dunlap and Burns. The counter-question is why you have such trouble accepting advice. I don't have trouble accepting advice. I am just curious about the error. Thank you very much for your patience. Regards, Jinsong David Winsemius, MD West Hartford, CT Regards, Jinsong __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] color2D.maplot fixed color range
Hello, Were you ever able to find a solution to your problem? I have been trying on and off to solve this problem for several months now. I would love to hear if you found a solution. -- View this message in context: http://r.789695.n4.nabble.com/color2D-maplot-fixed-color-range-tp3298679p4184667.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] Color2D.matplot uniform color range
Hi all. I am having a problem using color2D.matplot I am trying to visualize several different matrices under two color ranges. One color range corresponds to values less than 1 and the second color range for values greater than 1. However, the minimum value of each matrix differs and is automatically set to have the smallest value in the color range. Similarly, the maximum value of each matrix differs leading to non-uniform color scales. Below is the function that is causing me trouble. It is designed to take as input 2 csv files with row and column headers describing a 20X20 matrix. The goal would be to modify the code below as follows. (1) if an entry is in range [0,1] the color scale should always be uniform with the minimum color being 0.5 and the maximum color being 1. (2) if an entry is in the range [1,infinity] the color scale should have minimum color value 1 and maximum color value 3. (I used infinity simply because there is no way of knowing how high the value will be beforehand). Any help would be greatly appreciated. Thanks, Jav library(plotrix) make_S_figure-function(filename,alias){ h0 - read.csv(file=filename,head=TRUE,sep=,,row.names=1) d =data.matrix(h0) m - 1:20 cellcolors-matrix(NA,nrow=20,ncol=20) cellcolors2--matrix(NA,nrow=20,ncol=20) for(i in 1:length(m)){ cellcolors[d = 1]-color.scale(d[d=1],cs1=c(1,1),cs2=c(1,0),cs3=c(1,0)) cellcolors[d1]-color.scale(d[d1],cs1=c(0,1),cs2=c(0,1),cs3=c(1,1)) color2D.matplot(d,cellcolors=cellcolors,show.values=F,na.color=white,axes=FALSE,main=alias, xlab=,ylab=) axis(1,at=0.5:19.5,labels=colnames(d)) axis(2,at=0.5:19.5,labels=rev(rownames(d))) } -- View this message in context: http://r.789695.n4.nabble.com/Color2D-matplot-uniform-color-range-tp4184701p4184701.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] Is there a way to print branch distances for hclust function?
The R function hclust is used to do cluster analysis, but based on R help I see no way to print the actual fusion distances (that is, the vertical distances for each connected branch pairs seen in the cluster dendrogram). Any ideas? I'd like to use them test for significant differences from the mean fusion distance (i.e. The Best Cut Test). To perform a cluster analysis I'm using: x - dist(mydata, method = euclidean) # distance matrix y - hclust(x, method=ward) #clustering (i.e. fusion) method plot(y) # display dendogram Thanks, kbrownk __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 about spaces in r
Hi First of all, it's R, not r, and on this mailing list people care about this kind of thing. Second, you will need to provide more information in order to get better help. Please read the posting guide. There are a number of introductory level documents available via CRAN, please pick one and study the basics. That said, here are a few basics: When R reads data that is a mixture of letters and digits, as yours is, it will interpret all of them as characters. (Possibly converting to something called a factor, depending on exactly how the data is being input.) read.csv() produces data frames, not matrices. In a matrix, all values must be the same type, numeric or character. In a data frame, different columns can be different types, but each column must be all the same type. Your column has a mixture of tyeps. Your input column apparently has a mixture of types, violating the data frame rule, so read.csv() treats the numbers as characters. You can use as.numeric() on this column to convert to numeric. The Not if it is factor (which you can check by str(your.object). With factor you need to do as.numeric(as.character(factor.variable)) But it is always better to read data properly instead of changing them later after reading. To get more help you shall follow Don§s advice to provide some more info about what you have, what you did and what you get. Usefull commands are ?dput ?str And probably spending some of your time reading R-intro document. Regards Petr elements that look like numbers will become numeric, and the elements that have letters will be converted to NA for missing. Just guessing, but from what you've shown in looks like maybe your input data (outside R) is structured in a way that is not conducive to reading into R. This is quite common when the input data is a spreadsheet. -Don -- Don MacQueen Lawrence Livermore National Laboratory 7000 East Ave., L-627 Livermore, CA 94550 925-423-1062 On 12/9/11 7:44 AM, Matt Spitzer matthewjspit...@gmail.com wrote: Hello, I would like to please ask if someone would explain how r reads characters and numbers differently. Using read.csv, I had a matrix that resembled the following, only with many more ids and data: ID Visit variable 2 1 5 2 1 3 2 3 4 2 41 1 2 42 34 2 5 54 2 9 1 2 10 3 2 12 5 5 1 54 5 2 9 5 3 3 5 41 54 5 41 2 5 5 235 5 9 4 5 10 2 5 12 2 I then tried to subset for Visit==3. However, subset == was not working properly. This gave me zero rows. I printed the matrix/dataframe and found that this was because r viewed the 3 as 3 (space three). So, I had to type subset == 3 to select for the data instead. I think this has to do with character, number and string properties, but I am quite a novice. Would anyone be able to instruct me how one tells a dafaframe/matrix to convert numbers such as 3 to 3 so that I do not get confused in the future? I guess another problem I have is that I am still learning the differences between matrices and dataframes. Thanks so much, Matt [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Colours for sunflowerplot
Dear fellow R users, I would like to draw a sunflowerplot because I have data (decade by month) that plots multiple times on the same x-y co-ordinates. Further I would like to colour each of the points/sunflower leaves on the plot according to the group they belong to (i.e. which type of event each represents within that decade and month). I thought that this would be relatively straight forward - I have a series of x and y co-ords and a colour associated with each and I just use the column with the colours in it to colour-code the points using col = column.name when calling the plot. However, the sunflower plot uses the function xy.coords to calculate the number of times that multiple points are plotted at the same coordinates and in so doing creates a new dataset (x.coords, y.coords and the number of times each is repeated) and this dataset no longer has the colours associated with the individual points. The colours in a resultant plot merely plot in the original order. I would like to know whether there is a way of associating the correct colours with individual points that are now represented by the xy.coords output. If you are interested in this example, an example of the code that I have been using and the input data used can be found at the end of this message. Many thanks in advance for your input! Best wishes, Nicola Code: For sunflowerplot par(las = 1) sunflowerplot(extremes.decade$Decade, extremes.decade$Month, col = extremes.decade$Extreme, cex = 1.2, xlab = Decade, ylab = Month, axes = F , pch = 16, seg.col = extremes.decade$Extreme, size = 0.2, seg.lwd =2) axis(side = 1, at = c(6:11), labels = c(,1960s,1970s,1980s,1990s,2000s), tick = TRUE, line = NA) axis(side = 2, at = c(0:12), labels = c(,Jan,Feb,Mar,Apr,May,Jun,Jul,Aug,Sep,Oct,Nov, Dec), tick = TRUE, line = NA) legend(7, 12.5, legend = c(Grp1, Grp2, Grp3, Grp4, Grp5, Grp6, Grp7, Grp8), col = c(blue,brown,orange,green,red,lightblue,purple,pink), pch = 19, bty = n) For comparison of where colours should lie - normal plot: #add random number to make all points different so that they dont overplot a - runif(109,1,1.001) extremes.decade$Year.a - a*extremes.decade$Year plot(extremes.decade$Year.a, extremes.decade$Month, col = extremes.decade$Extreme, cex = 1.2, xlab = Decade, ylab = Month, yaxt =n , pch = 16) #axis(side = 1, at = c(6:11), labels = c(,1960s,1970s,1980s,1990s,2000s), tick = TRUE, line = NA) axis(side = 2, at = c(0:12), labels = c(,Jan,Feb,Mar,Apr,May,Jun,Jul,Aug,Sep,Oct,Nov, Dec), tick = TRUE, line = NA) legend(1950, 12.5, legend = c(Grp1, Grp2, Grp3, Grp4, Grp5, Grp6, Grp7, Grp8), col = c(blue,brown,orange,green,red,lightblue,purple,pink), pch = 19, bty = n) Data: Decade Year Month Extreme 7 1968 1 brown 7 1966 2 green 7 1966 2 lightblue 7 1968 2 green 7 1963 3 lightblue 7 1966 4 green 7 1961 5 lightblue 7 1967 5 blue 7 1964 6 blue 7 1965 6 green 7 1968 6 blue 7 1968 6 lightblue 7 1965 7 orange 7 1963 9 red 7 1964 9 lightblue 7 1968 9 blue 7 1961 10 green 7 1966 10 brown 7 1965 12 green 8 1973 1 brown 8 1976 1 red 8 1977 1 brown 8 1977 2 blue 8 1977 2 orange 8 1974 3 blue 8 1970 4 brown 8 1974 4 brown 8 1975 4 brown 8 1975 5 green 8 1973 6 red 8 1979 7 blue 8 1979 8 blue 8 1974 9 lightblue 8 1979 9 green 8 1974 10 brown 8 1975 10 green 8 1976 10 blue 8 1976 11 green 8 1970 12 lightblue 8 1975 12 lightblue 9 1980 1 pink 9 1982 1 pink 9 1980 2 pink 9 1986 2 pink 9 1981 3 green 9 1982 5 green 9 1987 5 red 9 1985 6 pink 9 1980 7 red 9 1983 7 blue 9 1980 8 pink 9 1981 8 lightblue 9 1985 8 red 9 1981 9 green 9 1982 9 pink 9 1985 10 orange 9 1987 10 brown 9 1983 11 pink 10 1990 1 red 10 1999 1 brown 10 1997 2 brown 10 1999 2 brown 10 1999 3 orange 10 1991 4 brown 10 1991 4 red 10 1997 4 lightblue 10 1997 6 blue 10 1997 6 lightblue 10 1994 7 green 10 1996 7 lightblue 10 1999 7 orange 10 1993 9 blue 10 1994 10 pink 10 1997 10 red 10 1994 11 brown 10 1996 11 blue 10 1996 11 lightblue 10 1993 12 blue 10 1997 12 purple 11 2004 1 orange 11 2001 2 brown 11 2006 2 orange 11 2006 2 red 11 2000 3 lightblue 11 2000 3 orange 11 2002 3 red 11 2004 3 purple 11 2003 4 orange 11 2004 4 purple 11 2001 5 purple 11 2003 5 purple 11 2002 7 purple 11 2005 7 red 11 2002 8 blue 11 2002 8 orange 11 2002 8 purple 11 2006 8 blue 11 2001 10 orange 11 2002 10 brown 11 2002 10 purple 11 2005 10 red 11 2002 11 purple 11 2004 11 orange 11 2004 11 red 11 2002 12 orange 11 2004 12 orange
Re: [R] Is there a way to print branch distances for hclust function?
On Sun, Dec 11, 2011 at 8:43 PM, kbrownk kbro...@gmail.com wrote: The R function hclust is used to do cluster analysis, but based on R help I see no way to print the actual fusion distances (that is, the vertical distances for each connected branch pairs seen in the cluster dendrogram). Any ideas? I'd like to use them test for significant differences from the mean fusion distance (i.e. The Best Cut Test). To perform a cluster analysis I'm using: x - dist(mydata, method = euclidean) # distance matrix y - hclust(x, method=ward) #clustering (i.e. fusion) method plot(y) # display dendogram Thanks, kbrownk You need to dig a bit deeper in the help file :) The return value is a list that contains, among others, components 'merge' and 'height'. The 'merge' component tells you which objects were merged at each particular step, and the 'height' component tells you what the merging height at that step was. The (slightly) tricky part is to relate the merge component to actual objects - AFAIK there is no function for that. The function cutree() using the argument k and varying it between 2 and n should basically do it for you but you need to match it to the entries in 'merge'. Maybe someone else knows a better way to do this. HTH, Peter __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.