Re: [R] sin(pi)?
Umm. pi has been rounded to 6 decimal places in the second example. So it isn't surprising that the results differ. sin(pi) is not zero, as it also has been rounded, and you can't represent irrational numbers exactly in a numerical form anyway. R agrees with Octave: octave:1 sin(pi) ans = 1.2246e-16 octave:2 sin(3.141593) ans = -3.4641e-07 octave:3 To paraphrase someone else on this list: I think it is strange that you think it is strange. Simon. As someone On Mon, 2007-09-03 at 16:43 +1000, Nguyen Dinh Nguyen wrote: Dear all, I found something strange when calculating sin of pi value sin(pi) [1] 1.224606e-16 pi [1] 3.141593 sin(3.141593) [1] -3.464102e-07 Any help and comment should be appreciated. Regards Nguyen Nguyen Dinh Nguyen Garvan Institute of Medical Research Sydney, Australia __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320 Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au Policies: 1. I will NOT analyse your data for you. 2. Your deadline is your problem. The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Interpreting the eigen value of a population matrix (2nd try)
To get a confidence interval on lambda, you need to have measures of variability in the elements of the transition matrix. If you have that, you can use a parametric bootstrap to get approximate confidence intervals. I have done this, and it seems to work. Alternatively, you could calculate a Bayesian posterior density for lambda using the Bayesian melding methods developed by Adrian Raftery et al., and calculate an HPD interval from that. I've done that too. It's slightly more difficult, however. Simon. Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au Policies: 1. I will NOT analyse your data for you. 2. Your deadline is your problem. The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. -Original Message- From: [EMAIL PROTECTED] on behalf of Anouk Simard Sent: Wed 29/08/2007 1:17 AM To: r-help@stat.math.ethz.ch Subject: [R] Interpreting the eigen value of a population matrix (2nd try) Thanks for telling me that you could not get my message, I hope this work better... so my question was: I built a population matrix to which I applied the fonction eigen in order to find the main parameters about my population. I know that the first eigen value correspond to lambda or exponential growth rate of my population. My problem is that I want to have the 95% confidence interval of the specific lambda (1.056 in the case). Is there a way to do that? Are the other eigen value shown in the output could help me doing it. I would very appreciate any help. Thanks for your time $values [1] 1.0561867+0.000i 0.0749653+0.5249157i 0.0749653-0.5249157i [4] 0.4498348+0.0795373i 0.4498348-0.0795373i -0.3357868+0.000i $vectors [1,] -0.72849129+0i -0.11058308+0.3293511i -0.11058308-0.3293511i 0.00244042+0.03012017i 0.00244042-0.03012017i [2,] -0.41384232+0i 0.35124594+0.1765638i 0.35124594-0.1765638i 0.01004458+0.03839895i 0.01004458-0.03839895i [3,] -0.27427879+0i 0.29630718-0.4260863i 0.29630718+0.4260863i 0.02540181+0.05526223i 0.02540181-0.05526223i [4,] -0.34274458+0i -0.62502691+0.000i -0.62502691+0.000i 0.55688585-0.17705587i 0.55688585+0.17705587i [5,] -0.31754610+0i 0.19351247+0.1625154i 0.19351247-0.1625154i -0.73460380+0.i -0.73460380+0.i [6,] -0.06705781+0i -0.00340804-0.0295753i -0.00340804+0.0295753i 0.30711075+0.13557984i 0.30711075-0.13557984i __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Interpreting the eigen value of a population matrix (2nd try)
[I've cc'ed this to the list because I think that it may be of value to other useRs, at least for archiving purposes.] Sorry, I did that work a long time ago in XLispStat (before I switched to R). The delta method and the bootstrap method for standard errors for eigenvalues are described in Caswell (2001) Matrix Population Models, 2nd Ed., Chapter 12. Simon. On Wed, 2007-08-29 at 13:31 -0400, Anouk Simard wrote: Thanks a lot for your help! Following you advice I have a couple more questions. I you could help me with it I would really appreciate although I would understand if you don't have time for it. I do have SE estimate on the survival and reproductive rate at each age class but I don't know how to incorporate this in a matrix model in order to do the bootstrap. So I was just wondering if you have any example of codes for it. I am not very familiar with R coding for this and I would very appeciate any help! Thank for your time again, Anouk rom: Simon Blomberg s.blomberg1 at uq.edu.au Subject: Re: Interpreting the eigen value of a population matrix (2nd try) Newsgroups: gmane.comp.lang.r.general Date: 2007-08-29 13:10:59 GMT (3 hours and 31 minutes ago) To get a confidence interval on lambda, you need to have measures of variability in the elements of the transition matrix. If you have that, you can use a parametric bootstrap to get approximate confidence intervals. I have done this, and it seems to work. Alternatively, you could calculate a Bayesian posterior density for lambda using the Bayesian melding methods developed by Adrian Raftery et al., and calculate an HPD interval from that. I've done that too. It's slightly more difficult, however. Simon. Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320 Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au Policies: 1. I will NOT analyse your data for you. 2. Your deadline is your problem. The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Sphericity test in R for repeated measures ANOVA
?mauchly.test also ?anova.mlm Another approach would be to do a direct test of compound symmetry by fitting models with and without that structure, and compare them via a likelihood ratio test. Simon. On Mon, 2007-08-20 at 19:08 -1000, Orou Gaoue wrote: Hi, Is there a way to do a sphericity test in R for repeated measures ANOVA (using aov or lme)? I can't find anything about it in the help. Thanks Orou [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320 Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au Policies: 1. I will NOT analyse your data for you. 2. Your deadline is your problem. The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] LSD, HSD,...
If you have a priori planned comparisons, you can just test those using linear contrasts, with no need to correct for multiple testing. If you do not, and you are relying on looking at the data and analysis to tell you which treatment means to compare, and you are considering several tests, then you should consider correcting for multiple testing. There is a large literature on the properties of the various tests. (Tukey HSD usually works pretty well for me.) rant Why do people design experiments with a priori hypotheses in mind, yet test them using post hoc comparison procedures? It's as if they are afraid to admit that they had hypotheses to begin with! Far better to test what you had planned to test using the more powerful methods for planned comparisons, and leave it at that. /rant On Mon, 2007-07-16 at 09:52 +0200, Adrian J. Montero Calvo wrote: Hi, I'm designing a experiment in order to compare the growing of several clones of a tree specie. It will be a complete randomized block design. How can I decide what model of mean comparision to choose? LSD, HSD,TukeyHSD, Duncan,...? Thanks in advance __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320 Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au Policies: 1. I will NOT analyse your data for you. 2. Your deadline is your problem. The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] type III ANOVA for a nested linear model
I second the nomination! Simon. On Tue, 2007-07-10 at 10:02 -0600, Greg Snow wrote: I nominate the following 2 pieces from Bill's reply for fortunes (probably 2 separate fortunes): All this becomes even more glaring if you take the unusal step of plotting the data. and What sort of editor would overlook this clear and demonstrable message leaping out from the data in favour of some arcane argument about types of sums of squares? Several answers come to mind: A power freak, a SAS afficianado, an idiot. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320 Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au Policies: 1. I will NOT analyse your data for you. 2. Your deadline is your problem. The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] EM algorithm for Missing Data.
Sure! Read this: MAXIMUM LIKELIHOOD FROM INCOMPLETE DATA VIA EM ALGORITHM Author(s): DEMPSTER AP, LAIRD NM, RUBIN DB Source: JOURNAL OF THE ROYAL STATISTICAL SOCIETY SERIES B-METHODOLOGICAL 39 (1): 1-38 1977 then read the posting guide. Simon. On Sun, 2007-07-08 at 23:20 -0300, Marcus Vinicius wrote: Dear all, I need to use the EM algorithm where data are missing. Example: x- c(60.87, NA, 61.53, 72.20, 68.96, NA, 68.35, 68.11, NA, 71.38) May anyone help me? Thanks. Marcus Vinicius [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Fees to use R
Of course, if your company would like to donate some money to the R foundation, I'm sure it would be greatly appreciated. http://www.r-project.org/foundation/donations.html Cheers, Simon. Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. -Original Message- From: [EMAIL PROTECTED] on behalf of Gavin Simpson Sent: Sat 7/7/2007 6:36 AM To: [EMAIL PROTECTED] Cc: r-help Subject: Re: [R] Fees to use R On Fri, 2007-07-06 at 09:31 +0200, [EMAIL PROTECTED] wrote: Good morning to all, I work for a bank in Italy, I want to know if i can install R and relative add on like Rbloomberg for free or my company has to pay some fee. tanks to all. Stefano Colucci R is released under the GNU GPL licence version 2. You can read the licence online here: http://www.gnu.org/copyleft/gpl.html As such R is free (as in beer) and you can install it without paying a fee. The source code is also free (as in speech) and is available from www.r-project.org as are pre-compiled binaries for various systems. You are however bound by the GPL licence and you should evaluate the implication of the GPL for the use you/your employer has in mind. HTH G __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] ?replace characters within vector data
sub(xxx, aaa, vectorx) or maybe gsub, depending on your application. Cheers, Simon. On Fri, 2007-07-06 at 12:40 +1000, [EMAIL PROTECTED] wrote: Hi List, I want replace characters within a vector. Outside R I could use sed, but I'd like to automate it in R. For example vectorx xxxyyz xxxyyza xxxyyzzb I want to change to: vectorx aaayyz aaayyza aaayyzzb The obvious replace command only deals with whole data entries? Any hints would be appreciated. Thanks Herry __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Parsimony Analysis of Encemism in R?
Short answer: No. Longer answer: R currently has some useful but relatively limited functions for phylogenetics. See the Statistical genetics task view on CRAN. The ape package is the most full-featured. Simon. On Tue, 2007-07-03 at 07:59 -0700, Milton Cezar Ribeiro wrote: Hi R-gurus, Is there a package for Parsimony Analysis of Endemism (Cladist) in R? Kind regards, Miltinho Brazil http://yahoo.com.br/oqueeuganhocomisso [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Extracting sums for individual factors in data frames
Does this do what you want? with(dat, tapply(BA, Species, sum)) ACSA AEGL Dead FRAM VIPR 565.172518 11.780972 157.393792 122.993352 3.926991 Cheers, Simon. On Sun, 2007-07-01 at 23:15 -0400, James R. Milks wrote: I have a data frame with two columns, one of which is a factor (Species) and the other is numeric (BA, which stands for basal area). Here's a sample: Species BA ACSA 55.7632696 FRAM 122.9933524 ACSA 67.54424205 ACSA 89.22123136 ACSA 82.46680716 ACSA 22.46238747 ACSA 19.94911335 ACSA 20.42035225 ACSA 19.00663555 ACSA 21.67698931 ACSA 57.80530483 ACSA 30.31636911 Dead 43.98229715 Dead 40.21238597 Dead 16.49336143 Dead 40.21238597 Dead 16.49336143 ACSA 78.53981634 VIPR 3.926990817 AEGL 11.78097245 AEGL 0 AEGL 0 ACSA 0 ACSA 0 ACSA 0 VIPR 0 I would like to calculate relative basal area for each species in this plot. For that, I need to divide the total basal area per species by the total basal area in the plot. Getting the total basal area in the plot is easy. However, I'm mystified on how to get the total basal area per species. Is there a way to extract and/or sum the total basal area per species? Thank you in advance. Jim Milks Graduate Student Environmental Sciences Ph.D. Program Wright State University 3640 Colonel Glenn Hwy Dayton, OH 45435 [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] why this doesn't work for qqnorm
You need to specify a whole column. try qqnorm(table[,1])
On Thu, 2007-06-28 at 18:50 -0700, Jiong Zhang, PhD wrote:
I want to qqnorm every column in a table. When I try the first column
using
qqnorm(table$column1), it worked.
But when I use
qqnorm(table[1]), it tells me Error in stripchart(x1, ...) : invalid
plotting method.
What happen? How can I make a function that qqnorms every column?
Try this:
dat - data.frame(x=rnorm(20), y= rnorm(20), z =rnorm(20))
par(mfrow=c(3,1))
sapply(dat, function (x) {qqnorm (x) ; qqline (x) } )
Cheers,
Simon
thanks a lot.
-jiong
The email message (and any attachments) is for the sole use of the intended
recipient(s) and may contain confidential information. Any unauthorized
review, use, disclosure or distribution is prohibited. If you are not the
intended recipient, please contact the sender by reply email and destroy all
copies of the original message (and any attachments). Thank You.
[[alternative HTML version deleted]]
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Simon Blomberg, BSc (Hons), PhD, MAppStat.
Lecturer and Consultant Statistician
Faculty of Biological and Chemical Sciences
The University of Queensland
St. Lucia Queensland 4072
Australia
Room 320, Goddard Building (8)
T: +61 7 3365 2506
email: S.Blomberg1_at_uq.edu.au
The combination of some data and an aching desire for
an answer does not ensure that a reasonable answer can
be extracted from a given body of data. - John Tukey.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] unequal variance assumption for lme (mixed effect model)
The default settings for lme do assume equal variances within groups. You can change that by using the various varClasses. see ?varClasses. A simple example would be to allow unequal variances across groups. So if your call to lme was: lme(...,random=~1|group,...) then to allow each group to have its own variance, use: lme(...,random=~1|group, weights=varIdent(form=~1|group),...) You really really should read Pinheiro Bates (2000). It's all there. HTH, Simon. , On Wed, 2007-06-27 at 21:55 -0400, shirley zhang wrote: Dear Douglas and R-help, Does lme assume normal distribution AND equal variance among groups like anova() does? If it does, is there any method like unequal variance T-test (Welch T) in lme when each group has unequal variance in my data? Thanks, Shirley __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] ANOVA non-sphericity test and corrections (eg, Greenhouse-Geisser)
On Mon, 2007-06-25 at 17:53 +1000, Simon Blomberg wrote: If you use lme, you can fit a general correlation structure to the within-subject data, and compare the fit to a model assuming uncorrelated within-subjects errors. That should tell you whether your data are ... correlated... (damn email gremlins.) Aren't the G-G and H-F corrections only approximate fixes? Surely it is better to work with a model that actually fits your data, rather than using ad hoc adjustments towards a model that doesn't quite fit. But I'm no psychologist. :-) Cheers, Simon. On Mon, 2007-06-25 at 08:22 +0200, Peter Dalgaard wrote: DarrenWeber wrote: I'm an experimental psychologist and when I run ANOVA analysis in SPSS, I normally ask for a test of non-sphericity (Box's M-test). I also ask for output of the corrections for non-sphericity, such as Greenhouse-Geisser and Huhn-Feldt. These tests and correction factors are commonly used in the journals for experimental and other psychology reports. I have been switching from SPSS to R for over a year now, but I realize now that I don't have the non-sphericity test and correction factors. This can be done using anova.mlm() and mauchly.test() which work on mlm objects, i.e., lm() output where the response is a matrix. There is no theory, to my knowledge, to support it for general aov() models, the catch being that you need to have a within-subject covariance matrix. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] ANOVA non-sphericity test and corrections (eg, Greenhouse-Geisser)
If you use lme, you can fit a general correlation structure to the within-subject data, and compare the fit to a model assuming uncorrelated within-subjects errors. That should tell you whether your data are Aren't the G-G and H-F corrections only approximate fixes? Surely it is better to work with a model that actually fits your data, rather than using ad hoc adjustments towards a model that doesn't quite fit. But I'm no psychologist. :-) Cheers, Simon. On Mon, 2007-06-25 at 08:22 +0200, Peter Dalgaard wrote: DarrenWeber wrote: I'm an experimental psychologist and when I run ANOVA analysis in SPSS, I normally ask for a test of non-sphericity (Box's M-test). I also ask for output of the corrections for non-sphericity, such as Greenhouse-Geisser and Huhn-Feldt. These tests and correction factors are commonly used in the journals for experimental and other psychology reports. I have been switching from SPSS to R for over a year now, but I realize now that I don't have the non-sphericity test and correction factors. This can be done using anova.mlm() and mauchly.test() which work on mlm objects, i.e., lm() output where the response is a matrix. There is no theory, to my knowledge, to support it for general aov() models, the catch being that you need to have a within-subject covariance matrix. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] R vs. Splus in Pharma/Devices Industry
Furthermore, one day we will all be doing our statistics in C. So higher level programming is predicted to become more difficult! Cheers, Simon. Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. -Original Message- From: [EMAIL PROTECTED] on behalf of Ted Harding Sent: Sat 6/16/2007 9:33 AM To: r-help@stat.math.ethz.ch Subject: Re: [R] R vs. Splus in Pharma/Devices Industry On 15-Jun-07 19:30:45, Douglas Bates wrote: On 6/15/07, Nicholas Lewin-Koh [EMAIL PROTECTED] wrote: Hi, Probably when the statistical community is using Z big pharma will be ready to use R. %P I think you mean A, not Z. First there was S, then there was R. And, furthermore, A is further from R than Z is, which seems fitting. Ted. E-Mail: (Ted Harding) [EMAIL PROTECTED] Fax-to-email: +44 (0)870 094 0861 Date: 15-Jun-07 Time: 23:57:13 -- XFMail -- __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] if statement
My solutions are usually too baroque, but does this do what you want? x - rnorm(100) quants - quantile(x, c(.3, .7)) Case - rep(2, length(x)) # 2 lies in the middle of the distribution Case[x = quants[1]] - 0 Case[x = quants[2]] - 1 Case Cheers, Simon. On Mon, 2007-06-11 at 15:14 -0700, Jiong Zhang, PhD wrote: Hi all, I have a rather naive question. I have the height of 100 individuals in a table and I want to assign the tallest 30% as Case=1 and the bottom 30% as Case=0. How do I do that? thanks. jiong The email message (and any attachments) is for the sole use of the intended recipient(s) and may contain confidential information. Any unauthorized review, use, disclosure or distribution is prohibited. If you are not the intended recipient, please contact the sender by reply email and destroy all copies of the original message (and any attachments). Thank You. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] PCA for Binary data
You might try (detrended) correspondence analysis, which is designed for count data, if it makes sense to treat your binary data that way. I've used ade4 and also vegan, and they are both good packages for these types of ordinations. You could also look at non-metric multidimensional scaling. There seems to be 2 schools of ordination. The Europeans like eigenanalysis methods (PCA, correspondence analysis, multiple correspondence analysis, coinertia analysis etc.). The Americans seem to prefer MDS. Cheers, Simon. This is On Tue, 2007-06-12 at 20:17 -0700, Spencer Graves wrote: The problem with applying prcomp to binary data is that it's not clear what problem you are solving. The standard principal components and factor analysis models assume that the observations are linear combinations of unobserved common factors (shared variability), normally distributed, plus normal noise, independent between observations and variables. Those assumptions are clearly violated for binary data. RSiteSearch(PCA for binary data) produced references to 'ade4' and 'FactoMineR'. Have you considered these? I have not used them, but FactoMineR included functions for 'Multiple Factor Analysis for Mixed [quantitative and qualitative] Data' Hope this helps. Spencer Graves Josh Gilbert wrote: I don't understand, what's wrong with using prcomp in this situation? On Sunday 10 June 2007 12:50 pm, Ranga Chandra Gudivada wrote: Hi, I was wondering whether there is any package implementing Principal Component Analysis for Binary data Thanks chandra - [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Expand duplicated observations
Does this do what you want?
dat - c(NA, 0, 3.2, 4)
fn - function (x) {
z - round(x)
if (is.na(x) | x = 1) z else rep(z, each=z)
}
unlist(sapply(dat, fn))
[1] NA 0 3 3 3 4 4 4 4
HTH,
Simon.
On Wed, 2007-06-06 at 01:54 +0100, M. P. Papadatos wrote:
Dear all,
I am trying to expand duplicated observations. I need to replace
each observation in the dataset with n copies of the observation,
where n is equal to the required expression rounded to the nearest
integer. If the expression is less than 1 or equal to missing, it is
interpreted as if it were 1, and the observation is retained but not
duplicated.
Example
From
c(1,2,3)
To
c(1,2,2,3,3,3)
Thank you in advance.
Best wishes,
Martin
--Apple-Mail-4-920612661--
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Simon Blomberg, BSc (Hons), PhD, MAppStat.
Lecturer and Consultant Statistician
Faculty of Biological and Chemical Sciences
The University of Queensland
St. Lucia Queensland 4072
Australia
Room 320, Goddard Building (8)
T: +61 7 3365 2506
email: S.Blomberg1_at_uq.edu.au
The combination of some data and an aching desire for
an answer does not ensure that a reasonable answer can
be extracted from a given body of data. - John Tukey.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Expand duplicated observations
D'Oh! yet again my first inclination is to write something complicated when a little thought shows a short, neat solution. Ah, well, I live and learn. Cheers, Simon. On Wed, 2007-06-06 at 09:59 +0800, Jared O'Connell wrote: Also, to handle NAs and non-integers: x = c(1:3,9.4,NA) tmp = round(x) tmp[is.na(tmp)]=1 rep(x,tmp) [1] 1.0 2.0 2.0 3.0 3.0 3.0 9.4 9.4 9.4 9.4 9.4 9.4 9.4 9.4 9.4 NA On 6/6/07, Francisco J. Zagmutt [EMAIL PROTECTED] wrote: I think this will do what you want x=c(1,2,3) rep(x,x) [1] 1 2 2 3 3 3 Regards Francisco M. P. Papadatos wrote: Dear all, I am trying to expand duplicated observations. I need to replace each observation in the dataset with n copies of the observation, where n is equal to the required expression rounded to the nearest integer. If the expression is less than 1 or equal to missing, it is interpreted as if it were 1, and the observation is retained but not duplicated. Example From c(1,2,3) To c(1,2,2,3,3,3) Thank you in advance. Best wishes, Martin __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] fixed effects anova in lme lmer
?gls in package nlme. It's like lme but with no random effects. But you can still model the variance-covariance properties of the data. Simon. On Tue, 2007-06-05 at 19:11 -0700, [EMAIL PROTECTED] wrote: Can lme or lmer fit a plain regular fixed effects anova? Ie a model without a random effect, or have there be at least one random effect in order for these functions to work? Trying to run such, (1) without specifying a random effect produces an error, (2) specifying that there is no random effect does not produce the same output as an anova run in lm(); (2b) specifying that there is no random effect in lmer crashed R (division by zero, I think). Just trying to see the connection of fixed and random effects anova in R. STATA gives me same results for both models up to the point where they differ. Best Toby dt1 = as.data.frame(cbind(c(28,35,27,21,21,36,25,18,26,38,27,17,16,25,22,18),c(1,1,1,1,2,2,2,2,3,3,3,3,4,4,4,4))) summary(a1 - lm(V1~factor(V2)-1, dt1)) anova(a1) summary(a1 - lm(V1~factor(V2), dt1)) anova(a1) dt1$f = factor(dt1$V2) summary(a2 - lme(V1~f, dt1)) #1a summary(a2 - lme(V1~f, dt1, ~-1|f)) #2a anova(a2) lmer(V1~f, dt1) #1b lmer(V1~f+(-1|f), dt1) #2b __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] ANOVA results in R conflicting with results in other software packages
R uses treatment contrasts for factors (ie 0/1 coding) by default. Systat is using sum (ie sum to zero) contrasts: Try this: options(contrasts=c(contr.sum, contr.poly) lm(maladapt~host*increase*size2)-fm Anova(fm, type=III) I won't discuss the dangers of types of sums of squares and different contrast codings. That would be tempting the wrath of the gods. See section 7.18 in the R FAQ. John Fox's Companion book also has a brief discussion (p. 140). Cheers, Simon. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Reduced Error Logistic Regression, and R?
From what I've read (which isn't much), the idea is to estimate a utility (preference) function for discrete categories, using logistic regression, under the assumption that the residuals of the linear predictor of the utilities are ~ Type I Gumbel. This implies the independence of irrelevant alternatives in economic jargon. ie the utility of choice a versus choice b is independent of the introduction of a third choice c. It also implies homoscedasticity of the errors. The model can be generalized in various ways. If you are willing to introduce extra parameters into the model, such as the parameters of the Gumbel distribution, you may get more precision in the estimates of the utility function. An alternative (without the independence of irrelevant alternatives assumption) is to model the errors as multivariate normal (ie use probit regression), which is computationally much more difficult. Whether it makes substantive sense to use these models outside of discrete choice experiments is another question. Patenting these methods is worrying. There have been a lot of people working on discrete choice experiments over the years. It's hard to believe that a single company could have ownership over an idea that is the result of a collaborative effort such as this. Cheers, Simon. On Thu, 2007-04-26 at 12:29 +1000, Tim Churches wrote: This news item in a data mining newsletter makes various claims for a technique called Reduced Error Logistic Regression: http://www.kdnuggets.com/news/2007/n08/12i.html In brief, are these (ambitious) claims justified and if so, has this technique been implemented in R (or does anyone have any plans to do so)? Tim C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Suggestions for statistical computing course
On Fri, 2007-04-20 at 12:13 -0400, Fred Bacon wrote: Ideally, it would work like this: The free VMware player is installed on each of the lab computers. The lab manager uses a licensed copy of VMware Workstation to create a clean image of a computer. You can use the open source QEMU program to create VMware machines. http://fabrice.bellard.free.fr/qemu/ After installing QEMU, the following command creates a machine with 20 Gb disk space, onto which you can load a (licensed!) copy of Windows (or better, Linux :-) ): qemu-img.exe create -f vmdk VMmachine.vmdk 20G The instructor makes a copy of the clean image and installs the necessary software and instructional materials. The instructor can use either the free player or the paid workstation version to do this. After the virtual machine is completed, the image is sent back to the lab where it is made available to the lab computers. If you use the paid workstation version rather than the free player version on the lab computers, then you can use the Snapshot feature to create a consistent image for every student. Every time the virtual machine is shutdown, the system can revert back to the snapshot for the next student. It all depends on your budget. Again, you can do this for free with QEMU, using the -snapshot option. How you handle the OS licensing issue for the guest operating system is up to you. I personally would recommend using Linux, but some of our customers are terrified of anything that doesn't look like a Microsoft OS. The only caveat is the disk space utilization. Having a complete OS image for every student for every class could eat up terabytes of space. But heck, terabyte RAID arrays are readily available these days. Fred -- Simon Blomberg, BSc (Hons), PhD, MAppStat. Lecturer and Consultant Statistician Faculty of Biological and Chemical Sciences The University of Queensland St. Lucia Queensland 4072 Australia Room 320, Goddard Building (8) T: +61 7 3365 2506 email: S.Blomberg1_at_uq.edu.au The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] novice/beginner's reading list for non-programmers learning R?
You will find a lot of information on CRAN, from simple cheat sheets, through to advanced books. Perhaps start with the freely-downloadable documents until you find your feet? Cheers, Simon. new ruser wrote: Can someone please recommend a novice/beginner's reading list for non-programmers learning R? - 8:00? 8:25? 8:40? Find a flick in no time [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Using boot.ci to find lower confident bound
Start by reading the references in ?boot.ci, to understand the different methods. Then you pays your money and you takes your choice. (Hint: R is free.) Cheers, Simon. Nguyen Manh The wrote: Dear R-users, I am trying to find lower confident bound (one side CI) for a variable of interest using boot.ci command in bootstrap method. In the boot.ci, it provides 5 methods for finding CI (two-side CI). Anyone can help me with that? Or any suggestions? Do I have to use a different commands? Many thanks, The Nguyen Manh The Department of Statistics University of Glasgow Don't pick lemons. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Help on variable ranking
Before you do that, you might try reading this paper: Bring, J. 1995. Variable importance by partitioning R^2. Quality and Quantity 29:173-189. Cheers, Simon. Andrew Robinson wrote: Rupendra, depending on the nature of your data (which you haven't mentioned), you might try hierarchical partitioning, as found in the hier.part package on CRAN. Cheers Andrew On Wed, Jan 17, 2007 at 11:07:18AM +0545, Rupendra Chulyadyo wrote: Hello all, I want to assign relative score to the predictor variables on the basis of its influence on the dependent variable. But I could not find any standard statistical approach appropriate for this purpose. Please suggest the possible approaches. Thanks in advance, Rupendra Chulyadyo Institute of Engineering, Tribhuvan University, Nepal [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] SNA Matrix
?source
Salvaj, Erica wrote:
Hello
I export a one mode network from Pajek to R, and the former made an .r file
called PajekR.r, that is actually an script to be run in R
The problem is that what the file actually does is to set a 0 martrix and
then assing to each pair of nodes the corresponding values. Since the matrix
is huge (10.000 nodes) a copy/paste procedure is very time consuming, and I
guess pretty inefficient. So the question is: are there any way to read
directly the .r file in R, without copy and paste the sentences of the script
in the R console?
Erica H. Salvaj
PhD Candidate
IESE Business School
[EMAIL PROTECTED]
This message has been scanned for viruses by TRENDMICRO,\ an...{{dropped}}
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat.
Centre for Resource and Environmental Studies
The Australian National University
Canberra ACT 0200
Australia
T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au
F: +61 2 6125 0757
CRICOS Provider # 00120C
The combination of some data and an aching desire for
an answer does not ensure that a reasonable answer
can be extracted from a given body of data.
- John Tukey.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Geometric Brownian Process
This function samples from the solution to the stochastic differential
equation of the geometric wiener process:
rgwiener - function (n=1, t=1, S0=1000, mu=0, sigma=1) {
S0 * exp((mu - 1/2 * sigma^2) * t + sigma * sqrt(t) * rnorm(n))
}
To test this, note that the geometric wiener process has log(S_t/S0) ~
N((mu -1/2*sigma^2)* t, sigma^2 * t)
So if we let mu = 0, sigma = 1, t =1, S0 = 1000, log(S_t/S0) should be ~
Normal(-0.5, 1).
samp - log(rgwiener(10)/1000)
hist(samp)
library(MASS)
fitdistr(samp, normal)
meansd
-0.5005531891.32814
( 0.003162381) ( 0.002236141)
which looks good to me.
Because the geometric wiener process is a transformation of the ordinary
wiener process, we can simulate it easily enough (apologies to the
authors of package e1071, whose function rwiener I hacked):
gwiener - function (mu=0, sigma=1, S0=1000, frequency=1000) {
z - S0 * exp((mu - 1/2 * sigma^2) * seq(0, 1, length=frequency) +
sigma * cumsum(rnorm(frequency)/sqrt(frequency)))
ts(z, start = 1/frequency, frequency=frequency)
}
plot(gwiener(), type=l)
Looks ok.
Cheers,
Simon.
Michael Graber wrote:
Dear R People,
Consider I have 3 realizations of an Geometric Brownian process.
Now i want to overlay the plot of these realizations. In a future point
in time a probability density curve in this specific point of time
should overlay this plot ( view rotated 90°).
I am sorry for not providing any source code. Can anybody point mo to an
package, or has anybody an idea how to simulate an geometric brownian
process in R?
Thanks in advance,
Michael Graber
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat.
Centre for Resource and Environmental Studies
The Australian National University
Canberra ACT 0200
Australia
T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au
F: +61 2 6125 0757
CRICOS Provider # 00120C
The combination of some data and an aching desire for
an answer does not ensure that a reasonable answer
can be extracted from a given body of data.
- John Tukey.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to sum one column in a data frame keyed on other columns
You could look at the reshape package, using sum as the aggregate function.
HTH,
Simon.
George Nachman wrote:
I have a data frame that looks like this:
url time somethingirrelevant visits
www.foo.com 1:00 xxx 100
www.foo.com 1:00 yyy 50
www.foo.com 2:00 xyz 25
www.bar.com 1:00 xxx 200
www.bar.com 1:00 zzz 200
www.foo.com 2:00 xxx 500
I'd like to write some code that takes this as input and outputs
something like this:
url time total_vists
www.foo.com 1:00 150
www.foo.com 2:00 525
www.bar.com 1:00 400
In other words, I need to calculate the sum of visits for each unique
tuple of (url,time).
I can do it with this code, but it's very slow, and doesn't seem like
the right approach:
keys = list()
getkey = function(m,cols,index) { paste(m[index,cols],collapse=,) }
for (i in 1:nrow(data)) { keys[[getkey(data,1:2,i)]] = 0 }
for (i in 1:nrow(data)) { keys[[getkey(data,1:2,i)]] =
keys[[getkey(data,1:2,i)]] + data[i,4] }
I'm sure there's a more functional-programming approach to this
problem! Any ideas?
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat.
Centre for Resource and Environmental Studies
The Australian National University
Canberra ACT 0200
Australia
T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au
F: +61 2 6125 0757
CRICOS Provider # 00120C
The combination of some data and an aching desire for
an answer does not ensure that a reasonable answer
can be extracted from a given body of data.
- John Tukey.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] lmer, p-values and all that
Try using mcmcsamp() to sample from the posterior distribution of the parameter estimates. You can calculate a p-value from that, if that is your desire. Instructions are in the R wiki: http://wiki.r-project.org/rwiki/doku.php?id=guides:lmer-tests HTH, Simon. Dan Bebber wrote: Hello, I've just located the illuminating explanation by Douglas Bates on degrees of freedom in mixed models. The take-home message appears to be: don't trust the p-values from lme. Questions: Should I give up hypothesis testing for fixed effects terms in mixed models? Has my time spent reading Pinheiro Bates been in vain? Is there a publication on this issue? Thanks, Dan Bebber Department of Plant Sciences University of Oxford __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Adding terms to a function
How about this:
form - formula(paste(y ~, paste(x, 1:n, sep=, collapse= + )))
model - lm(form)
HTH,
Simon.
Brooke LaFlamme wrote:
Hi all,
I am running R version 2.4.0 on Windows XP. I am new and have the following
question:
I have a dataset of columns named x1, x2, x3...xn. I would like to write a
linear regression using lm that looks like this:
lm(y~x1+x2+x3+...+xn)
If I try to use the following code, I only get the model for y~x1+xn:
n-ncol(dataset)
model-lm(y~x1)
for(i in 1:n) {
model.new-update(model, .~.+dataset[,i])
}
The purpose of this is so I can use stepAIC with model.new as the upper scope
and model as the lower.
I know there must be a simple way to do this, but I am not yet familiar with
much syntax. Any help appreciated!
--
Brooke LaFlamme
Cornell University
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat.
Centre for Resource and Environmental Studies
The Australian National University
Canberra ACT 0200
Australia
T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au
F: +61 2 6125 0757
CRICOS Provider # 00120C
The combination of some data and an aching desire for
an answer does not ensure that a reasonable answer
can be extracted from a given body of data.
- John Tukey.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Nested for loop
Please supply some code that reproduces the problem! cheers, Simon. Natalie Zayats wrote: Hi all, Does anybody know whether one can nest IF statement in multiple FOR loops ? According to results of my code, only the first iteration of each loop is performed... Thanks, Regards, Natalie [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Chi-Square Goodness-of-Fit test
See pearson.test in the nortest package. Also, read the notes section in ?pearson.test. You may not really want to do this test. HTH, Simon. Ethan Johnsons wrote: Do you know/have a function that takes a vector x and provides a returned p-value that uses the Chi-Square Goodness-of-Fit test to test the goodness of fit of a standard normal distribution. Awaiting your positive reply. Thx ej __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] lmer and a response that is a proportion
, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] comments in scan
From ?scan : If |comment.char| occurs (except inside a quoted character field), it signals that the rest of the line should be regarded as a comment and be discarded. Lines beginning with a comment character (possibly after white space with the default separator) are treated as blank lines. So, set comment.char=# and you should be fine. Cheers, Simon. Jarrett Byrnes wrote: I had a question about scan in R. For better code readability, I would like to have lines in the block of data to be scanned that are commented - not just lines that have a comment at the end. For example #age, weight, height 33,128,65 34,56,155 instead of having to do something like 33,128,65 #age, weight, height 34,56,155 Is this at all possible? __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] comments in scan
You may want to set blank.lines.skip=TRUE too. Simon. From ?scan : If |comment.char| occurs (except inside a quoted character field), it signals that the rest of the line should be regarded as a comment and be discarded. Lines beginning with a comment character (possibly after white space with the default separator) are treated as blank lines. So, set comment.char=# and you should be fine. Cheers, Simon. Jarrett Byrnes wrote: I had a question about scan in R. For better code readability, I would like to have lines in the block of data to be scanned that are commented - not just lines that have a comment at the end. For example #age, weight, height 33,128,65 34,56,155 instead of having to do something like 33,128,65 #age, weight, height 34,56,155 Is this at all possible? __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] ANOVA in Randomized-complete blocks design
You have mis-transcribed the data:
SeriesGenotypeWeight
1pb0.986
3bb0.829
anova(aov(weight ~ series + genotype, data=dat))
gives the correct results when compared to SR.
Cheers,
Simon.
January Weiner wrote:
Dear all,
I am trying to repeat an example from Sokal and Rohlfs Biometry --
Box 11.4, example of a randomized-complete-blocks experiment.
The data is fairly simple:
series genotype weight
1 pp 0.958
1 pb 0.985
1 bb 0.925
2 pp 0.971
2 pb 1.051
2 bb 0.952
3 pp 0.927
3 pb 0.891
3 bb 0.892
4 pp 0.971
4 pb 1.010
4 bb 0.955
The model used for ANOVA in the book is
Y_{ij} = \mu + \alpha_i + B_i + [(\alpha B)_{ij}] + \epsilon_{ij}
(I am not quite confident how to represent this model in R, see below)
The ANOVA table from SR looks like this:
MSBseries 3 0.021 0.007 10.23 **
MSAgenotypes 2 0.010 0.005 6.97 *
MSEerror 6 0.004 0.001
In R, I am using the following model (is this correct?)
aov.ob = aov( genotype ~ genotype + series )
However, my results are
DfSum Sq Mean Sq F value Pr(F)
series 3 0.0135867 0.0045289 6.6360 0.02469 *
genotype 2 0.0056732 0.0028366 4.1563 0.07367 .
Residuals6 0.0040948 0.0006825
What am I doing wrong?
Regards,
January
--
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat.
Centre for Resource and Environmental Studies
The Australian National University
Canberra ACT 0200
Australia
T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au
F: +61 2 6125 0757
CRICOS Provider # 00120C
The combination of some data and an aching desire for
an answer does not ensure that a reasonable answer
can be extracted from a given body of data.
- John Tukey.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] ANOVA in Randomized-complete blocks design
January Weiner wrote: Ops, sorry. OK, as you see, I am going through SR and doing the examples using R. I have a further question, on the nested ANOVA example from Box 10.1 (mosquito female wing lengths). I made sure that the data is correct :-) (data below). I am not sure how to create in R a Model II nested ANOVA. aov( wing ~ cage / female ) which is, I believe, a fixed effect nested ANOVA where females are nested within cages produced the correct sum of squares, but computates a different F value for the cage effect: Df Sum Sq Mean Sq F valuePr(F) cage 2 665.68 332.84 255.70 1.452e-10 *** cage:female 9 1720.68 191.19 146.88 6.981e-11 *** Residuals 12 15.621.30 (whereas in SR, among cages F = 1.741; the rest is same). The F value for the cage:female effect is the same as in SR. Why do I get a higly significant cage effect? In SR, the significance of the cage effect is done by F = MSamong / MSsubgr = 332.84 / 191.19 = 1.74. In the R model, F = 255.70, and I do not understand where this value comes from (at first I thought that it is MSamong / MSerror = 332.84 / 1.3 = 256.0308). I am confused... aov(wing ~ cage/female) gives you a test of cage against the residual variance and cage:female against the residual variance. As you noticed, the residual variance is the wrong error stratum for the cage effect. To get the correct test of the cage effect, try: summary(aov(wing ~ cage + Error(female))) You will get the correct F test from that. To get the variance components, it is easier to use the lme function in the nlme package: fit - lme(wing ~ cage, random=~1|female) summary(fit) where cage is fixed and female is random. Compare: fit - lme(wing ~ 1, random=~1|cage/female) summary(fit) This treats both cage and female as random effects. VarCorr(fit) will give you the variance components, or you can read off the standard deviations (sqrt of variance components) from summary(fit). Yet another way to analyse the problem is with the lmer function (package lme4). Cheers, Simon. How should I modify the model? Another question: is there a way to automatically estimate the variance components, or do I have to take the respective MS' and calculate it myself? Thanks, January OK, here is the data: cage female wing cage1 f1 58.5 cage1 f1 59.5 cage1 f2 77.8 cage1 f2 80.9 cage1 f3 84.0 cage1 f3 83.6 cage1 f4 70.1 cage1 f4 68.3 cage2 f1 69.8 cage2 f1 69.8 cage2 f2 56.0 cage2 f2 54.5 cage2 f3 50.7 cage2 f3 49.3 cage2 f4 63.8 cage2 f4 65.8 cage3 f1 56.6 cage3 f1 57.5 cage3 f2 77.8 cage3 f2 79.2 cage3 f3 69.9 cage3 f3 69.2 cage3 f4 62.1 cage3 f4 64.5 Cheers, January -- __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] New package Ryacas
Hi Gabor,
I'm running Quantian (Debian) inside a VMware virtual machine, on a
Windows XP host.
I installed the latest version of yacas from the source tarball. I
remembered to ./configure --enable-server to allow server connections.
make and make install worked ok, after some fiddling. I checked that the
yacas server option worked, by doing yacas --server , and then
telnet'ing to 127.0.0.1 to check. It worked fine. I installed
Ryacas. I then tried it out and got the following error:
library(Ryacas)
Loading required package: XML
yacas('Integrate(x)x;')
[1] Starting Yacas!
Error in socketConnection(host = 127.0.0.1, port = 9734, server =
FALSE, :
unable to open connection
In addition: Warning message:
127.0.0.1:9734 cannot be opened
Accepting requests from port 9734
I tried again (stubborn, I guess):
yacas('Integrate(x)x;')
[1] Starting Yacas!
Accepting requests from port 9734
YacasServer Could not bind to the socket
: Address already in use
/usr/local/lib/R/site-library/Ryacas/yacdir/R.ys(1) : File not found
CommandLine(1) : Expecting ) closing bracket for sub-expression, but
got x instead
Any ideas where I may be going wrong? I don't know anything about
sockets. I've cross-posted to r-sig-debian. They may be interested.
Cheers,
Simon.
--
Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat.
Centre for Resource and Environmental Studies
The Australian National University
Canberra ACT 0200
Australia
T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au
F: +61 2 6125 0757
CRICOS Provider # 00120C
The combination of some data and an aching desire for
an answer does not ensure that a reasonable answer
can be extracted from a given body of data.
- John Tukey.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] remotely saving an R session
Does save.image do what you want? Gamal Azim wrote: Forgot to indicate that the remote system is Linux, accessed remotely by ssh. --- Gamal Azim [EMAIL PROTECTED] wrote: Is it possible to remotely save an R session then terminate R? Of course the destructive task after 'then' is rather straightforward by itself. Thanks __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Extracting overdispersion estimates from lmer amd glm objects
summary(modeltest)@sigma Toby Gardner wrote: Dear list, I am needing to extract the estimate of overdispersion (deviance / residual degrees of freedom or c-hat) from multiple model objects - so they can then be used to compare the extent of overdispersion among alternative models as well as calculate qausi-AIC values. I have been unable to do this, despite consulting a number of manuals and searching the R-help. I am imaging that in theory it should be possible with some call to attr(), but i have so far had no success. An example model output would be: modeltest-lmer(Coleodactylus_amazonicus_N~USD + (1|site),data=SFArray,family=poisson,method=Laplace,control=list(usePQL=FALSE, msVerbose=TRUE)) summary(modeltest) Generalized linear mixed model fit using Laplace Formula: Coleodactylus_amazonicus_N ~ USD + (1 | site) Data: SFArray Family: poisson(log link) AIC BIClogLik deviance 75.94996 81.68603 -34.97498 69.94996 Random effects: Groups NameVariance Std.Dev. site (Intercept) 2.6076 1.6148 number of obs: 50, groups: site, 5 Estimated scale (compare to 1) 1.080798 What I need is to extract this value (1.080798) from multiple lmer objects. Has anyone any recommendations? I also need to do this for glm objects although I suspect if someone was able to kindly point me in the right direction then the solution is likely to be similar. Very many thanks, Toby Gardner sessionInfo() Version 2.3.1 (2006-06-01) i386-pc-mingw32 attached base packages: [1] datasets graphics grDevices methods stats utils base other attached packages: JGR iplots JavaGD lme4 Matrixlattice MASS rJava 1.4-71.0-30.3-4 0.995-2 0.995-15 0.13-8 7.2-27.1 0.4-6 School of Environmental Sciences University of East Anglia Norwich, NR4 7TJ United Kingdom Email: [EMAIL PROTECTED] Website: www.uea.ac.uk/~e387495 [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] R in Nature
Hi all, We've just had a paper accepted for publication in Nature. We used R for 95% of our analyses (one of my co-authors sneaked in some GenStat when I wasn't looking.). The preprint is available from the Nature web site, in the open peer-review trial section. I searched Nature for previous references to R Development Core Team, and I received no hits. So I tentatively conclude that our paper is the first Nature paper to cite R. A great many thanks to the R Development Core Team for R, and Prof. Bates for lmer. Cheers, Simon. (I'm off to the pub to celebrate.) -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] R in Nature
AAh. Then my hypothesis has been rejected. Oh well! Cheers, Simon. Simon, Congratulations! It used to be that R Ihaka and R Gentleman R: A Language for Data Analysis and Graphics Journal of Computational and Graphical Statistics, 1996 was used to cite R. I see a handful of these in Nature from around 2003. Chuck On Fri, 25 Aug 2006, Simon Blomberg wrote: Hi all, We've just had a paper accepted for publication in Nature. We used R for 95% of our analyses (one of my co-authors sneaked in some GenStat when I wasn't looking.). The preprint is available from the Nature web site, in the open peer-review trial section. I searched Nature for previous references to R Development Core Team, and I received no hits. So I tentatively conclude that our paper is the first Nature paper to cite R. A great many thanks to the R Development Core Team for R, and Prof. Bates for lmer. Cheers, Simon. (I'm off to the pub to celebrate.) -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. - John Tukey. [ Part 3.91: Included Message ] Charles C. Berry(858) 534-2098 Dept of Family/Preventive Medicine E mailto:[EMAIL PROTECTED] UC San Diego http://biostat.ucsd.edu/~cberry/ La Jolla, San Diego 92093-0717 -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] merge 2 data frame based on more than 2 variables
Wensui Liu wrote: Dear Lister, I understand merge() can be used to join 2 data frames based on 1 variable. But how about merge based on more than 2 variables? Thank you so much! Just specify the 2 (or more) variable names in a column vector for by) merge(dat1, dat2, by= c(VarA, VarB)) assuming both data frames have columns VarA and VarB. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] merge 2 data frame based on more than 2 variables
Then instead of by, use by.x and by.y to specifiy the variable names separately for both data frames. See ?merge, especially the examples. Wensui Liu wrote: what if the names are different in 2 data frames? On 8/14/06, Simon Blomberg [EMAIL PROTECTED] wrote: Wensui Liu wrote: Dear Lister, I understand merge() can be used to join 2 data frames based on 1 variable. But how about merge based on more than 2 variables? Thank you so much! Just specify the 2 (or more) variable names in a column vector for by) merge(dat1, dat2, by= c(VarA, VarB)) assuming both data frames have columns VarA and VarB. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] basic question re lm()
You could look at using lm.fit instead of lm.
Alternatively, you can paste the names of the variables together using
the following approach. It's a bit baroque, but it works:
form.fn - function (dframe) {
nms - names(dframe)
formula(paste(nms[1], ~, paste(nms[2:length(nms)], collapse=+)))
}
Then do:
form - form.fn(data1)
lm(form, data=data1)
HTH,
Simon.
r user wrote:
I am using R in a Windows environment.
I have a basic question regarding lm().
I have a dataframe “data1” with ncol=w.
I know that my dependent variable is in column1.
Is there a way to write the regression formula so that
I can use columns 2 thru w as my independent
variables?
e.g. something like: “ lm(data1[,1] ~ data1[,2:w] ) “
Thanks
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat.
Centre for Resource and Environmental Studies
The Australian National University
Canberra ACT 0200
Australia
T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au
F: +61 2 6125 0757
CRICOS Provider # 00120C
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Data transformation
If x is your data, and is a vector, then rep(x, each=2) should do it for you. Cheers, Simon. jenny tan wrote: Hi there, Wonder if someone who is R-savvy can help me with the following task (see below) that I occasionally do work and gets quite tedious if I do it manually. Thanks in advance! jenny. - I have a column of data that looks like this: NA18501 NA18502 NA18504 NA18505 NA18507 NA18508 NA18516 NA18517 NA18522 NA18523 And I want to duplicate the values and sort of interweave them to look like this: NA18501 NA18501 NA18502 NA18502 NA18504 NA18504 NA18505 NA18505 NA18507 NA18507 NA18508 NA18508 NA18516 NA18516 NA18517 NA18517 NA18522 NA18522 NA18523 NA18523 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Help with matrix manipulation
apply is your friend:
fn - function (x, a, b) {
if (x a) return(a) else
if (x b) return(b) else
x
}
apply(mat, c(1,2), fn, -4, 3)
[,1] [,2] [,3] [,4]
[1,] -4 -13 -4
[2,] -403 -3
[3,] -313 -2
[4,] -22 -4 -1
HTH,
Simon.
Kartik Pappu wrote:
Hi all,
I have a square (a x a) matrix with values in a range. For example:
[,1] [,2] [,3] [,4]
[1,] -5 -13 -4
[2,] -404 -3
[3,] -315 -2
[4,] -22 -5 -1
I want to take any number smaller than -4 (in this example -5) and
replace it with -4 and similarly take any number greater than 3 (in
this case 4 and 5) and replace it with 3. The other numbers (and the
overall structure of the matrix should remain unchanged.
Seems like something that would use an if ab then c else d kind of
logic, but I cannot figure out how to manipulate the entire matrix.
Thanks much
Kartik
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat.
Centre for Resource and Environmental Studies
The Australian National University
Canberra ACT 0200
Australia
T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au
F: +61 2 6125 0757
CRICOS Provider # 00120C
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to obtain 95th percentile of a normal distribution of a continuous variable
?quantile jenny tan wrote: Hi, How do I get R to output the 95% cutoff from a distribution of a continous variable? summary() only displays a few statistics Thanks! __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] planned comparisons for ANOVA
Darren Weber wrote: [snip] Our planned comparisons are: 1. test the group mean difference for S1 vs S2 (in the absence of S3) 2. test the group mean difference for S2 vs S3 (in the absence of S1) This is the current form of the ANOVA specification for R: aov( Y ~ (Task*Hemisphere*Group) + Error( Subject/(Task*Hemisphere) ) How can we add planned comparisons to this specification? Can we add just one planned comparison matrix, with rows for 1 2 above, or do we need to run the model twice, once for each planned comparison? There are a couple of ways to do this in R. Perhaps the easiest is to use make.contrasts in the gmodels package. cmat - rbind(S1 v S2 = c(1, -1, 0), S2 v S3 = c(0, 1, -1)) library(gmodels) fit - aov( Y ~ Task*Hemisphere*Group + Error( Subject/(Task*Hemisphere), contrasts=list(Task=make.contrasts(cmat ))) summary(fit) Alternatively, is there a function to compute the planned comparisons after running the full ANOVA model? see ?fit.contrast in gmodels for this. HTH, Simon. Thanks, Darren PS, I was trained well on using SPSS, but I am trying to make a switch to R. You help with this would be really appreciated. We need to get our results revised and published soon. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Scalling/Centering the Data by an Index
If you want to centre the data, but not scale as well, turn off scale: unlist(tapply(x, group, scale, scale=FALSE)) HTH, Simon. Ashraf Chaudhary wrote: Dear All: I would like to center the data in 'x' by 'group'. The following code scale the data and I have not been able to figure out how to change it so I get the centered data. x - c(1, 2, 3, 4, 5, 6, 7, 8) group - c(1,1,1,2,2,2,2,2) unsplit(lapply(split(x,group),scale),group) I would appreciate your help. Ashraf -__ M. Ashraf Chaudhary, Ph.D. Associate Scientist/Biostatistician Department of International Health Johns Hopkins Bloomberg School of Public Health 615 N. Wolfe St. Room W5506 Baltimore MD 21205-2179 (410) 502-0741/fax: (410) 502-6733 mailto:[EMAIL PROTECTED] [EMAIL PROTECTED] Web:http://faculty.jhsph.edu/?F=MohammadL=Chaudhary [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Suggestion for ?split
Hi all, I noticed an undocumented feature for split. It sorts the resulting list according to the grouping factor. An example: test - data.frame(x=rnorm(48), f=letters[sample(1:8)]) split(test, test$f) I wasn't expecting this behaviour, although I was pleasantly surprised. I suggest that the help page for split be amended to include this feature. I know it's a small thing, but someone else may also find it useful to know. Cheers, Simon. -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Decimal series how to make.........
?seq anil kumar rohilla wrote: Hi List, I am new to this Rsoftware, i want to make a sereis for example which is having values like this, s- 0,0.1,0.2,0.3,0.4,1 i tryed this statement s-0:0.1:1 but this giving an error megssage. but by default increment 1 it is taking ,so what to do ,, i want to use this varible in for loop.. like for(j in s) thanks in advance ANIL KUMAR( METEOROLOGIST) LRF SECTION NATIONAL CLIMATE CENTER ADGM(RESEARCH) INDIA METEOROLOGICAL DEPARTMENT SHIVIJI NAGAR PUNE-411005 INDIA MOBILE +919422023277 [EMAIL PROTECTED] [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Zero-inflated Poisson Repeated Measures Data
I have had success with AD model builder: http://otter-rsch.ca/admbre/examples/glmmadmb/glmmADMB.html Cheers, Simon. [EMAIL PROTECTED] wrote: Does someone have code, or point to a source to get it, to model repeated measures zero-inflated poisson data. The data come from a replicated field trial comparing two treatments - a control and a test treatment. Thanks in advance Subhash Chandra, DSc Senior Biometrician Primary Industries Research Victoria Department of Primary Industries 1 Ferguson Road Tatura 3616, Victoria Australia Tel: (03) 5833 5397 Email: [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Help In Function
Use readline.
x - readline(Give the value of x? )
cat(The value of x is, x, \n)
Cheers,
Simon.
SUMANTA BASAK wrote:
Hi All,
I need a basic help from you. I've built a function like this,
windowlength-function(x)
{
z - rep(seq(0,331,by=x-1)+1, each=2)
zz - z[-c(1,length(z))]
ind - as.data.frame(matrix(zz, nr=2))
j-lapply(ind, function(x) mat[x[1]:x[2],])
cat(For,x/4,month number of windows is = ,length(ind),\n)
}
windowlength(x=12)
I need to know how can i give command in R so that instead of giving the
last line, i.e R will ask the user to give the value of x? I mean to say,
1) It will ask user Give the value of x
2) Then user inputs 12, and R gives the ultimate result.
Thanks,
Sumanta Basak.
-
What makes Sachin India's highest paid sports celebrity?, Share your
knowledge on Yahoo! India Answers
Send instant messages to your online friends - NOW
[[alternative HTML version deleted]]
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
--
Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat.
Centre for Resource and Environmental Studies
The Australian National University
Canberra ACT 0200
Australia
T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au
F: +61 2 6125 0757
CRICOS Provider # 00120C
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] ape comparative analysis query
Chris Knight wrote: This raises various questions: 1) Was I misleading myself that my independent contrasts were valid in the first place? I think partly. Independent contrasts were designed for data for tip taxa only, not ancestral state data. One of the steps in the algorithm requires a correction to account for the fact that ordinarily the values for traits at nodes are estimated and not known, basically lengthening the branch length to that node by a certain amount (See Felsenstein's original paper). If you have known ancestral states, then you should not make this correction. 2) What is it, if anything, about the root taxon that causes this issue, given that other taxa also have zero branch lengths? When converting a tree to a correlation matrix for use in GLS or GEE, the Brownian motion assumption implies that the covariances among taxa are represented by their shared branch length from the root, and the height of each terminal taxon is equal to the variance of the trait for that taxon. If the tree is ultrametric, then the variances are equal. In your case, the root has zero covariance with each taxon, and zero variance. Now, the error in the gls call is caused by a division by zero in the calculation of the correlation matrix, because the root has zero variance. try: vcv.phylo(tree) vcv.phylo(tree, cor=TRUE) 3) Is there any way of getting around this and including data on the root taxon, or am I better off just dropping it (ultimately I want to work with much larger trees (up to tens of thousands of taxa) where that one piece of information will become relatively less important) I'm surprised that you have known values for the root data. Are you sure that the root taxon is not actually an outgroup? If so, it will have some branch length (variance), and the model should fit OK. HTH, Simon. Any help very much appreciated, Chris I'm working with ape 1.8-2 in R 2.1.1 under ubuntu 'Breezy' linux (unfortunately 2.1.1 is the latest easily available in breezy) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] How to find row maximum of two columns
apply(tmax.m, 2, max) Dale Steele wrote: Hello -- Given the 10x2 matrix below, my goal is to create a vector which contains the maximum of column 1 and column 2, or the only value if there is an NA in one column. I experimented with max.col without success. Thanks. --Dale tmax.m tmaxhme tmaxer [1,] 101.0 99.8 [2,] 102.5 99.0 [3,] 100.6 98.4 [4,] NA 100.5 [5,] 101.0 99.4 [6,] NA 97.6 [7,] NA 99.0 [8,]99.0 98.4 [9,] NA 98.5 [10,] NA 99.1 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] predict.glmmPQL Problem
I have found a similar problem when constructing formulae and passing them to glmmPQL. My solution was to use do.call(glmmPQL,...). See ?do.call. HTH, Simon. pencer Graves wrote: Please try again after upgrading to the versions of R and MASS. If you still have a problem, PLEASE do read the posting guide! www.R-project.org/posting-guide.html before you submit another post. Most of the people who donate their time to answer questions on this listserve are more likely to answer a question if it is simple and completely self contained -- including a very brief toy example that they can copy from an email into R and reproduce the problem. If they can't do that, they are less likely to understand your question and therefore less likely to produce a useful answer -- and less likely to bother to even read carefully your question. hope this helps, spencer graves rsubcriber wrote: Dear all, for a cross-validation I have to use predict.glmmPQL() , where the formula of the corresponding glmmPQL call is not given explicitly, but constructed using as.formula. However, this does not work as expected: x1-rnorm(100); x2-rbinom(100,3,0.5); y-rpois(100,2) mydata-data.frame(x1,x2,y) library(MASS) # works as expected model1-glmmPQL(y~x1, ~1 | factor(x2), family=poisson, data=mydata) predict(model1, newdata=mydata, type=response) f-as.formula(paste(y, ~,x1)) # predict does not work: # Error in mCall[[fixed]][-2] : object is not subsettable model2-glmmPQL(f, ~1 | factor(x2), family=poisson, data=mydata) predict(model2, newdata=mydata, type=response) Has anyone an idea what goes wrong? I am using R 2.0.1 under Windows 2000. Thanks, Anneke __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] write.table command
The help file for load says that it is to be used for objects saved using the save command. Perhaps you should try read.table. Cheers, Simon. Abd Rahman Kassim wrote: Dear All, I'm trying to save a dataframe using write.table command. It works, but when I retrieved, there's an error message as shown below: write.table(soil.dat,file=C:/soil.rdata) load(C:/soil.rdata) Error: bad restore file magic number (file may be corrupted)-- no data loaded I can figure out the error message. Any assistance to solve the problem is very much appreciated. Regards. Abd. Rahman Kassim, PhD Forest Management Ecology Program Forestry Conservation Division Forest Research Institute Malaysia Kepong 52109 Selangor Malaysia Fax: 603-62729852 Tel: 603-62797179 * * [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] How do I tell it which directory to use?
I think you want getwd() and setwd(). HTH, Simon. Thomas L Jones wrote: From Tom: In R 2.2.0 under Windows, I want to be able to give it a filename such as myFile.txt without the quotes. But actually I mean: C:\Documents and Settings\Tom\My Documents\qpaper7\R Project Started 19 Dec 05\myFile.txt If I were to repeat this each time, my computer would get all bored and cranky and start to drop bits (only a joke, of course). I think I want to set the Home directory or the working directory or some directory or other to the above directory. I may or may not want to set some environmental variables. R 2.2.0; working directly from the console and copying and pasting code which I want to test into the console. Windows XP Home Edition. Administrator privileges are enabled. A curve ball: There are two accounts, Tom and Jones; the data are stored under Tom, whereas the computation is being done under the Jones account. I won't bore you with the details of why I am doing this. I was able to call Sys.getenv (R_USER) and get the home directory. I am a newbie to R and not familiar with the terminology. Tom Thomas L. Jones, Ph.D., Computer Science __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] stripchart-y axis labels
Try par(las=2) Cheers, Simon. [EMAIL PROTECTED] wrote: Hello, I am trying to create a stripchart for my data, where y axis labels are characters (ie,names of cities). I would like to change the orientation of the y - axis labels, ie perpendicular to y axis. Below is the code i am using: par(srt=90) with(ozone.ne.trim, stripchart(Median~City,main = stripchart(ozone))) The par option doesn't seem to working. Kindly help. Thanks mathangi. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] How to sum values across multiple variables using a wildcard?
data - data.frame(var1=c(1,2,3), var2=c(3,4,5), var3=c(4,5,6), foo = c(100,200,300)) # sum rows with var in their name rowSums(data[, grep(var, names(data))]) 1 2 3 8 11 14 [EMAIL PROTECTED] wrote: I have a dataframe called data with 5 records (in rows) each of which has been scored on each of many variables (in columns). Five of the variables are named var1, var2, var3, var4, var5 using headers. The other variables are named using other conventions. I can create a new variable called var6 with the value 15 for each record with this code: var6=var1+var2+var3+var4+var5 but this is tedious for my real dataset with dozens of variables. I would rather use a wildcard to add up all the variables that begin with Var like this pseudocode: Var6=sum(var*) Any suggestions for implementing this in R? Thanks! Mark __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Help! What does this R command mean?
RTFM! Especially the Introduction to R tutorial. Type help.start() at the prompt. It should open your web browser. Choose Introduction to R. Only then will you attain enlightenment. Simon. Michael wrote: Hi all, R is so difficult. I am so desperate. What does the : mean in the following statement? What does the [, -1] mean? # Leaps takes a design matrix as argument: throw away the intercept # column or leaps will complain X - model.matrix(lm(V ~ I + D + W +G:I + P + N, election.table))[,-1] Thanks a lot! [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Help! What does this R command mean?
?formula Michael wrote: After reading your pointers, I still don't understand that particular usage of :: btw, the data after read.table(): Year V I D W G P N 1 1916 0.5168 1 1 0 2.229 4.252 3 2 1920 0.3612 1 0 1 -11.463 16.535 5 3 1924 0.4176 -1 -1 0 -3.872 5.161 10 On 1/29/06, Gabor Grothendieck [EMAIL PROTECTED] wrote: ?: ?[ Also try ?? Also 1. read the An Introduction to R manual. On R Home page click on Manuals under Documentation in the left hand column. 2. go to: http://cran.r-project.org/other-docs.html and read your choice of documents. 3. print out this reference card and keep it handy: http://www.rpad.org/Rpad/R-refcard.pdf On 1/29/06, Michael [EMAIL PROTECTED] wrote: Hi all, R is so difficult. I am so desperate. What does the : mean in the following statement? What does the [, -1] mean? # Leaps takes a design matrix as argument: throw away the intercept # column or leaps will complain X - model.matrix(lm(V ~ I + D + W +G:I + P + N, election.table))[,-1] Thanks a lot! [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Extracting variance components in lme
v - VarCorr(fm1Rail.lme) str(v) # get an idea of how v is structured. This suggests: as.numeric(v[1, 2]) [1] 24.80547 There may be easier and better ways HTH, Simon. At 11:02 AM 3/11/2005, you wrote: Consider the output for the inroductory Rail example in Mixed Effects Models in S and S-PLUS by Pinheiro and Bates: summary(fm1Rail.lme) Linear mixed-effects model fit by REML Data: Rail AIC BIC logLik 128.177 130.6766 -61.0885 Random effects: Formula: ~1 | Rail (Intercept) Residual StdDev:24.80547 4.020779 Fixed effects: travel ~ 1 Value Std.Error DF t-value p-value (Intercept) 66.5 10.17104 12 6.538173 0 Standardized Within-Group Residuals: Min Q1 Med Q3 Max -1.61882658 -0.28217671 0.03569328 0.21955784 1.61437744 Number of Observations: 18 Number of Groups: 6 I want to extract the variance components sigma = 4.020779 (the within component) and sigma_b = 24.80547 (the between component). I can get sigma easily: sigma - fm1Rail.lme$sigma but how can I get sigma_b ? Cheers, Murray Jorgensen -- Dr Murray Jorgensen http://www.stats.waikato.ac.nz/Staff/maj.html Department of Statistics, University of Waikato, Hamilton, New Zealand Email: [EMAIL PROTECTED]Fax 7 838 4155 Phone +64 7 838 4773 wkHome +64 7 825 0441Mobile 021 1395 862 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R- exp(-1000) ? - how to get R to give me an actual answer ?
At 03:42 PM 26/10/2005, you wrote: We can double check this with yacas which is a free computer algebra system that supports arbitrary precision math: In N(Exp(-1000)) Out 0.5075958898e-434 Or in Axiom, which is open source too: numeric exp(-1000) 0.50759588975494567653E-434 Type: Float On 10/26/05, rramacha [EMAIL PROTECTED] wrote: Dear All, I am a novice user of the R software package. When I try and compute, exp(-1000) or exp(-2000), i get the answer as zero. Is there any way i can get R to compute the answer and give me an actual number ? ( by increasing the precision or any other method). If I cannot get R to give me a number, can anybody give me some advice on how to manually compute this number ? I would greatly appreciate your help on this issue. Thanks Ravi Ravichandran Ravi Ramachandran Graduate Teaching Assistant Department of Statistics The University of Tennessee,Knoxville USA E-mail : [EMAIL PROTECTED] Phone : Official - 865 - 974 - 2739 (Aconda Court) Residence - 865 - 946 - 5155 Webpage : http://web.utk.edu/~rramacha Space maybe the final frontier But its made in a Hollywood basement - The Red Hot Chilli Peppers - 'Californication' Space maybe the final frontier But it can be 'replicated' in a Hollywood basement - Me, after taking the Stat 573 class Thank God every day when you get up that you have something to do that day which must be done whether you like it or not. Being forced to work and forced to do your best will breed in you temperance and self-control, diligence and strength of will - Basil Carpenter __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Population Projections
Depends what you want to do. Forming population projection matrices and doing the eigenanalysis to work out population growth rates, stable stage distributions, elasticities etc. is simple enough. Do you really need a package to do that? See ?eigen. What would you like to see in a population projection package? Cheers, Simon. At 06:44 AM 27/10/2005, you wrote: --- Rodriguez, Orlando [EMAIL PROTECTED] wrote: Is there an R Package for Population Projections? Orlando J. Rodriguez The Center for Population Research University of Connecticut Unit 2068 344 Mansfield Road Storrs, CT 06269-2068 (860)486-9269 (860)486-6356 (f) http://popcenter.uconn.edu Visit Dept of Demography at UC Berkeley courses page: http://www.demog.berkeley.edu/courses/ They use R and provide some code. Tomas Tomas Aragon, MD, DrPH Tel: 510-847-9139 (mobile) URL: http://www.idready.org __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Extracting Variance Components
?VarCorr At 12:02 PM 27/10/2005, you wrote: Dear List, Is there a way to extract variance components from lmeObjects or summary.lme objects without using intervals()? For my purposes I don't need the confidence intervals which I'm obtaining using parametric bootstrap. Thanks, Mike [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] a silly question on index of a matrix
?which X - matrix(2:10, nrow=3) which(X == 10, arr.ind=TRUE) row col [1,] 3 3 HTH, simon. At 11:52 AM 26/10/2005, you wrote: Hi netters, This is probably a silly question,but I can't find the answer after searching the R-help archives online. ok, I have a matrix. I know there is a 10 somewhere in it. Now I want to know the index of the element 10 in this matrix. That is, if X[i,j]=10, I want to know i and j. Is there a R function to do this? Just like the find function in matlab. Thanks all! __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] newbie questions - looping through hierarchial datafille
21.3528 F 0 21.3528 SFNSW_DIC:P F 21.3528 100 SFNSW_DIC:P T 2 25 L 0 32 23.1 F 0 6.5 SFNSW_DIC:A F 6.5 23.1 SFNSW_DIC:C F 23.1 100 SFNSW_DIC:C T 3 25 L 0 39.5 22.2407 F 0 4.7 SFNSW_DIC:A F 4.7 6.7 SFNSW_DIC:C P 2 20.25 slope=13 SPP:P.RAD T 1 25 L 0 38 22.1474 F 0 1 SFNSW_DIC:G F 1 2.3 SFNSW_DIC:A T 1001 25 L 0 38 22.1474 F 0 1 SFNSW_DIC:G F 1 2.3 SFNSW_DIC:A T 2 25 L 0 32.5 21.7386 F 0 2 SFNSW_DIC:A F 2 3.3 SFNSW_DIC:G F 3.3 10.4 SFNSW_DIC:C [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Translating lme model call to lme4
There is a slight caveat in that lmer does not respect implicit nesting, so you need to make sure your nested groups have unique levels: http://finzi.psych.upenn.edu/R/Rhelp02a/archive/47413.html http://finzi.psych.upenn.edu/R/Rhelp02a/archive/47423.html Simon. At 11:42 AM 13/09/2005, Doran, Harold wrote: Only the random portion will differ as in: lmer(lognrms ~ Group*Rotation*muscle*side*support*arms + (1|Subject) + (1|Stratum) + (1|rep), Data) -Original Message- From: [EMAIL PROTECTED] on behalf of Ross Darnell Sent: Mon 9/12/2005 9:28 PM To: r-help@stat.math.ethz.ch Cc: Subject:[R] Translating lme model call to lme4 I would appreciate help translating the following lme model to an lmer function. lme(lognrms ~ Group*Rotation*muscle*side*support*arms, random=~1|Subject/Stratum2/rep, data=Data) Many thanks Ross Darnell [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to do multiple comparisons in R?
R has various methods for multiple comparison procedures. See package multcomp, or ?TukeyHSD or ?pairwise.t.test for example. An RSiteSearch(multiple comparison) returned 187 results. A priori contrasts can be constructed using the make.contrasts function in the gmodels package, for example. We try not to do anything like in SAS. Simon. At 11:12 AM 12/09/2005, Hongyu Sun wrote: Hi, Sorry I have to bother a question. Does R have the functions to do lsd, tukey, bonferonni, contrast etc. like in SAS? Many thanks, HS __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] barplots - text direction
See las under ?par. You can often pass par parameters to higher-level graphics functions, so barplot(h, las=2) (for example) works. HTH, Simon. At 11:27 AM 24/08/2005, Murray Pung wrote: If the variable names are too long to allow room for each to be displayed on a barplot, how can the direction of the text be changed? a - cbind(2.4,2.4,2.5,2.6,2.6,2.6,2.6,2.6,2.9,2.9,2.9) b - cbind(2.3,2.5,2.4,2.2,3.2,2.4,2.9,2.6,2.9,3.0,2.8) h - rbind(a,b) colnames(h) - c(one,two,three,four,five,six,seven,eight,nine,ten,eleven) rownames(h) - c(Pre-stage,Post-stage) barplot(h, beside = T, legend = colnames(g), horiz = T, xlim = c(0, 5)) Many Thanks Murray Murray Pung | Research Analyst AIM Research HR Consulting PO Box 328, Nth Sydney NSW 2060 P +61 (02) 9956 3951 F +61 (02) 9922 2210 www.aimsurveys.com.au __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] extract t-values from pairwise.t.test
I think the questioner was interested in pairwise.t.test. See ?pairwise.t.test. From looking at the source, pairwise.t.test calls t.test if sd's are not pooled, or calculates its own t.val if sds are pooled. It looks very easy to hack to return the t values instead of the p values. Simon. At 04:46 PM 8/08/2005, Petr Pikal wrote: Hallo My output lists more than p-values ttt-t.test(rnorm(10), rnorm(10), paired=T) ttt Paired t-test data: rnorm(10) and rnorm(10) t = 1.7508, df = 9, p-value = 0.1139 alternative hypothesis: true difference in means is not equal to 0 95 percent confidence interval: -0.1750263 1.3735176 sample estimates: mean of the differences 0.5992456 str(ttt) List of 9 $ statistic : Named num 1.75 ..- attr(*, names)= chr t $ parameter : Named num 9 ..- attr(*, names)= chr df $ p.value: num 0.114 $ conf.int : atomic [1:2] -0.175 1.374 ..- attr(*, conf.level)= num 0.95 $ estimate : Named num 0.599 ..- attr(*, names)= chr mean of the differences $ null.value : Named num 0 ..- attr(*, names)= chr difference in means $ alternative: chr two.sided $ method : chr Paired t-test $ data.name : chr rnorm(10) and rnorm(10) - attr(*, class)= chr htest ttt$statistic t 1.750790 The output is list and you can call any part of it by its name or by [] braces. HTH Petr On 8 Aug 2005 at 16:26, Guido Parra Vergara wrote: Hi, how can I extract the t-values after running a pairwise.t.test? The output just list the p-values. Many thanks for your help. Cheers Guido Guido J. Parra School of Tropical Environment Studies and Geography James Cook University Townsville Queensland 4811 Phone: 61 7 47815824 Fax:61 7 47814020 Mobile: 0437630843 e-mail: [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Petr Pikal [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to test this
This is two tests: Whether the slope != 1 and whether the intercept != 0. To do this, include an offset in your model: fit - lm(y ~ x + offset(x), data=dat) HTH, Simon. At 03:44 PM 3/08/2005, [EMAIL PROTECTED] wrote: Dear there, I am wondering how to test whether a simple linear regression model (e.g. y=1.05x) is significantly different from a 1 to 1 line (i.e. y=x). Thanks. Regards, Jin [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Anova's in R
If I understand you correctly, this looks like a split-plot design. The anova is: aov(response~burning*temperature+Error(site), data=dataset) Alternatively in lme: lme(response~burning*temperature, random=~1|site, data=dataset) At 03:09 PM 29/07/2005, you wrote: Hello. I am looking for some help using anova's in RGui. My experiment ie data, has a fixed burning treatment (Factor A) 2 levels, unburnt/burnt. Nested within each level of Factor A are 2 random sites (Factor B). All sites are crossed with a fixed temperature treatment (Factor C) 2 levels, 0 degreesC/2 degreesC, with many replicates of these temperature treatments randomly located at each site. I am trying the following aov(dependent variable~burning*temperature*site+Error(replicate),data=dataset) and variations on that, however can't get it quite right the F ratios are not correct. I imagine this is a fairly common experimental design in ecology and would ask that anyone who has some advice please reply to this email? Thank-you, Frith __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Error in FUN(newX[, i], ...) : `x' must be atomic
Aah, OK. Thanks for the correction. I learn something new every day. That's why I like this list! Serves me right for relying on my Lisp knowledge. Cheers, Simon. At 07:26 PM 27/07/2005, Peter Dalgaard wrote: Simon Blomberg [EMAIL PROTECTED] writes: Actually, atoms are originally a Lisp concept. Objects are either atoms or not. Atoms are data types that cannot be taken apart, such as numbers or symbols. Lists and vectors (of length 1) are examples of non-atomic data types. Did you pass a vector to FUN? Actually, vectors ARE atomic in R, so your definition is somewhat off target. is.atomic(rnorm(5)) [1] TRUE In the extreme, the only true atom is the bit, everything else can be taken apart - doubles into (sign,exponent,mantissa) etc. So languages *define* their own atoms, as objects that are not composed of other objects in the language. (And even that is a partial lie for R, because atoms can have attributes. Atomicity is purely based on the type of an object. The help page has a list of the atomic types.) Cheers, Simon. At 12:22 PM 27/07/2005, Srinivas Iyyer wrote: Hello Group, What is the meaning of the error. is there any place to look for this. I guess 'atomic' seems to be OOP related concept. thank you srini __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Error in FUN(newX[, i], ...) : `x' must be atomic
Actually, atoms are originally a Lisp concept. Objects are either atoms or not. Atoms are data types that cannot be taken apart, such as numbers or symbols. Lists and vectors (of length 1) are examples of non-atomic data types. Did you pass a vector to FUN? Cheers, Simon. At 12:22 PM 27/07/2005, Srinivas Iyyer wrote: Hello Group, What is the meaning of the error. is there any place to look for this. I guess 'atomic' seems to be OOP related concept. thank you srini __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Nested ANOVA with a random nested factor (how to use the lme function?)
At 01:59 PM 18/07/2005, Addison, Prue wrote: Hi, I am having trouble using the lme function to perform a nested ANOVA with a random nested factor. My design is as follows: Location (n=6) (Random) Site nested within each Location (n=12) (2 Sites nested within each Location) (Random) Dependent variable: sp (species abundance) By using the aov function I can generate a nested ANOVA, however this assumes that my nested factor is fixed. summary(aov(sp~Location/Site, data=mavric)) Df Sum Sq Mean Sq F value Pr(F) Location 4 112366 28092 1.2742 0.2962 Location:Transect 5 121690 24338 1.1039 0.3736 Residuals 40 881875 22047 Yes, this is equivalent to aov(sp ~ Location + Location:Site,...) I have tried the following lme function to specify that Site is random: lme1 - lme(sp~Location, random=~1|Site, data=mavric) Here, Location is fixed, and Site is a grouping factor. There are fixed and random components to the intercept. lme2 - lme(sp~Location, random=~1|Location/Site, data=mavric) Here you have Location as a fixed effect, the intercept is random and the grouping is Site %in % Location. anova(lme1) numDF denDF F-value p-value (Intercept) 140 3.418077 0.0719 Location4 5 1.152505 0.4294 How can you have 6 sites, but only 4 df for Location (should be 5?) This gives me the correct F-value for Location from MSLocation/MSLocation:Transect, but the p-value doesn't seem to be correct (by my calculations in Microsoft Excel it should be 0.345) I don't know your data, or your calculations. Microsoft Excel does not fill me with confidence. anova(lme2) Warning in pf(q, df1, df2, lower.tail, log.p) : NaNs produced Warning: NAs introduced by coercion numDF denDF F-value p-value (Intercept) 140 0.459966 0.5015 Location4 0 0.155091 NaN ? I don't know what this output means Division by zero, somewhere? :-) anova(lme1,lme2) Model df AIC BIClogLik Test L.Ratio p-value lme1 1 7 603.7534 616.4000 -294.8767 lme2 2 8 605.7534 620.2067 -294.8767 1 vs 2 1.815674e-05 0.9966 ? I also don't know what this output means. you are testing the difference between the models, using a likelihood ratio test. The difference is not significant, so the conventional wisdom is to choose the simpler model (lme1). Note that the AIC and BIC are lower (better) for lme1, too. Can anyone tell me if there is a way to use the lme() function in order to obtain the same output as the aov() function (above), but so it correctly calculates the MS, F and p values for my main Location factor? The following code shows agreement: dat - data.frame(Location=factor(rep(1:6, each=2)), Site=factor(rep(1:2, 12)), sp=rnorm(12)) fit - aov(sp ~ Location + Error(Site), data=dat) fit2 - lme(sp ~ Location, random=~1|Site, data=dat) summary(fit) anova(fit2) HTH, Simon. Thanks, Prue This e-mail is solely for the named addressee and may be confidential. You should only read, disclose, transmit, copy, distribute, act in reliance on or commercialise the contents if you are authorised to do so. If you are not the intended recipient of this e-mail, please notify [EMAIL PROTECTED] by e-mail immediately, or notify the sender and then destroy any copy of this message. Views expressed in this e-mail are those of the individual sender, except where specifically stated to be those of an officer of Museum Victoria. Museum Victoria does not represent, warrant or guarantee that the integrity of this communication has been maintained nor that it is free from errors, virus or interference. Museum Victoria +61 3 8341 11 Nicholson St Carlton Victoria www.museum.vic.gov.au [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] help: how to get the position of a value in a matrix
See ?which Hint: arr.ind=TRUE Simon. At 09:28 AM 14/07/2005, wu sz wrote: Hello, I have a data set matrix of 1200 * 15. How can I get the position of a specific value in the matrix? I use seq(along = x)[x value] to look for the position of the value in the matrix, but seq can just find the sequence position row by row in the matrix, not a real position (like rowNumber, colNumber). Is any function for that? Thank you, Shengzhe __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] crossed random fx nlme lme4
At 09:35 AM 14/07/2005, Emilio A. Laca wrote: I need to specify a model similar to this lme.formula(fixed = sqrt(lbPerAc) ~ y + season + y:season, data = cy, random = ~y | observer/set, correlation = corARMA(q = 6)) except that observer and set are actually crossed instead of nested. Does this work for you? (following PB pp 162-3 and an R-help archive search on crossed random effects)... fit - lme(sqrt(lbPerAc) ~ y * season, random=list(pdBlocked(pdIdent(~y), pdIdent(observer-1), pdIdent(set-1))), correlation=corARMA(q = 6), data=cy) lme isn't very well set up for crossed random effects. It's easier in lmer. I don't think lmer can handle alternative correlation structures yet, though. (Prof. Bates?) HTH, Simon. observer and set are factors y and lbPerAc are numeric If you know how to do it or have suggestions for reading I will be grateful. eal ps I have already read Pinheiro Bates, the jan 05 newsletter, and several postings. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] CIs in predict?
At 08:40 AM 12/07/2005, Guy Forrester wrote:
Dear All,
I am trying to put some Confidence intervals on some regressions from a
linear model with no luck. I can extract the fitted values using
'predict', but am having difficulty in getting at the confidence
intervals, or the standard errors.
Any suggestions would be welcome
Cheers
Guy
Using Version 2.1.0 (2005-04-18) on a PC
vol.mod3 - lm(log.volume~log.area*lake,data=vol)
summary(vol.mod3)
plot(c(1.3,2.5),c(-0.7,0.45),type=n,xlab=Log area,ylab=Log volume)
areapred.a - seq(min(vol$log.area[vol$lake==a]),
max(vol$log.area[vol$lake==a]), length=100)
areapred.b - seq(min(vol$log.area[vol$lake==b]),
max(vol$log.area[vol$lake==b]), length=100)
preda - predict(vol.mod3,
data.frame(log.area=areapred.a,interval=confidence ,lake=rep(a,100)))
You have interval=confidence inside your call to data.frame, not inside
your call to predict. Hence you are creating a data frame with a variable
called interval, with one level called confidence, and predict does not see
interval=confidence at all! See ?predict.lm.
HTH,
Simon.
#This gives the fitted values as predicted, but no CIs
preda
12345
6789
-0.562577529 -0.553263576 -0.543949624 -0.534635671 -0.525321718
-0.516007765 -0.506693813 -0.497379860 -0.488065907
10 11 12 13 14
15 16 17 18
-0.478751955 -0.469438002 -0.460124049 -0.450810097 -0.441496144
-0.432182191 -0.422868239 -0.413554286 -0.404240333
19 20 21 22 23
24 25 26 27
-0.394926380 -0.385612428 -0.376298475 ETC ETC
#As does this, but with no SEs
preda - predict(vol.mod3, data.frame(log.area=areapred.a,se.fit=T
,lake=rep(a,100)))
preda
12345
6789 10
-0.562577529 -0.553263576 -0.543949624 -0.534635671 -0.525321718
-0.516007765 -0.506693813 -0.497379860 -0.488065907 -0.478751955
11 12 13 14 15
16 17 18 19 20
-0.469438002 -0.460124049 -0.450810097 ETC ETC
Guy J Forrester
Biometrician
Manaaki Whenua - Landcare Research
PO Box 69, Lincoln, New Zealand.
Tel. +64 3 325 6701 x3738
Fax +64 3 325 2418
E-mail [EMAIL PROTECTED]
www.LandcareResearch.co.nz
WARNING: This email and any attachments may be confidential ...{{dropped}}
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Finding out collinearity in regression
At 01:45 PM 30/06/2005, Young Cho wrote: Hi, I am trying to find out a collinearity in explanatory variables with alias(). I creat a dataframe: dat - ds[,sapply(ds,nlevels)=2] dat$Y - Response Explanatory variables are factor and response is continuous random variable. When I run a regression, I have the following error: fit - aov( Y ~ . , data = dat) Error in contrasts-(`*tmp*`, value = contr.treatment) : contrasts can be applied only to factors with 2 or more levels 1. Sounds like at least one of your factors has only one level. This should be easy to spot. 2. Have you considered package perturb? HTH, Simon. I think there is a dependency in explanatory variables. So, I wanted to use alias to find out a dependency in design matrix but I can't because I cannot create fit in the first place. One of examples I found is: carprice1.lm - lm(gpm100 ~ Type+Min.Price+Price+Max.Price+Range.Price,data=carprice) alias(carprice1.lm) But, what if I can create lm object ? Then is there a way to find out a dependency in design matrix? Thanks a lot for help in advance! -Young. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] sample R code for multiple imputation
There is also the mice package at http://www.multiple-imputation.com/ which uses mcmc but is different to Schafer's packages. Simon. At 09:30 AM 29/06/2005, James Reilly wrote: Schafer's functions have been ported to R in the packages norm, cat, mix and pan. Their documentation includes sample code illustrating how to use them. The aregImpute function in Hmisc provides a range of other imputation models, if you're not set on using MCMC. The mitools package might also be useful for dealing with the imputed values. Hope this helps, James On 29/06/2005 9:32 a.m., David Hwang wrote: Hi, I have a big dataset which has many missing values and want to implement Multiple imputation via Monte carlo markov chain by following J Schafer's Analysis of incomplete multivariate data. I don't know where to begin and is looking for a sample R code that implements multiple imputation with EM, MCMC, etc Any help / suggestion will be greatly appreciated. David Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] quotient and remainder
You may find the Reference Cards in the Contributed Documentation section of CRAN to be of interest. Simon. At 12:07 PM 24/06/2005, zhihua li wrote: Dear Dimitris, I've read the introduction to R thoroughly and gooogled in the internet about my question, but got no answer. I think it would be great if there's a doc grouping R functions into different functional categories. Thanks a lot for your replies! From: Dimitris Rizopoulos [EMAIL PROTECTED] To: zhihua li [EMAIL PROTECTED] CC: r-help@stat.math.ethz.ch Subject: Re: [R] quotient and remainder Date: Thu, 23 Jun 2005 09:01:08 +0200 11%/%5 [1] 2 11%%5 [1] 1 Best, Dimitris p.s., I'd suggest you to take a look at the An Introduction to R doc Dimitris Rizopoulos Ph.D. Student Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven, Belgium Tel: +32/16/336899 Fax: +32/16/337015 Web: http://www.med.kuleuven.be/biostat/ http://www.student.kuleuven.ac.be/~m0390867/dimitris.htm - Original Message - From: zhihua li [EMAIL PROTECTED] To: r-help@stat.math.ethz.ch Sent: Thursday, June 23, 2005 8:37 AM Subject: [R] quotient and remainder hi netters Is there a function in R that can compute the quotient and remainder of a division calculation? such that when 11 is given as the dividend and 5 the divider, the function returns 2(quotient) and 1(remainder). Thanks a lot! _ å è´¹ä¸è½½ MSN Explorer: http://explorer.msn.com/lccn/ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Centre for Resource and Environmental Studies The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au F: +61 2 6125 0757 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] logistic regression
predict.glm by default produces predictions on the scale of the linear predictors. If in a logistic regression, you want the predictions to be on the response scale [0,1], use x - predict(logistic.model, medians, type=response) for example. See ?predict.glm for details. Cheers, Simon. Hi I am working on corpora of automatically recognized utterances, looking for features that predict error in the hypothesis the recognizer is proposing. I am using the glm functions to do logistic regression. I do this type of thing: * logistic.model = glm(formula = similarity ~., family = binomial, data = data) and end up with a model: summary(logistic.model) Call: glm(formula = similarity ~ ., family = binomial, data = data) Deviance Residuals: Min 1Q Median 3Q Max -3.1599 0.2334 0.3307 0.4486 1.2471 Coefficients: Estimate Std. Error z value Pr(|z|) (Intercept) 11.1923783 4.6536898 2.405 0.01617 * length-0.3529775 0.2416538 -1.461 0.14410 meanPitch -0.0203590 0.0064752 -3.144 0.00167 ** minimumPitch 0.0257213 0.0053092 4.845 1.27e-06 *** maximumPitch -0.0003454 0.0030008 -0.115 0.90838 meanF1 0.0137880 0.0047035 2.931 0.00337 ** meanF2 0.0040238 0.0041684 0.965 0.33439 meanF3-0.0075497 0.0026751 -2.822 0.00477 ** meanF4-0.0005362 0.0007443 -0.720 0.47123 meanF5-0.0001560 0.0003936 -0.396 0.69187 ratioF2ToF10.2668678 2.8926149 0.092 0.92649 ratioF3ToF11.7339087 1.7655757 0.982 0.32607 jitter-5.2571384 10.8043359 -0.487 0.62656 shimmer -2.3040826 3.0581950 -0.753 0.45120 percentUnvoicedFrames 0.1959342 1.3041689 0.150 0.88058 numberOfVoiceBreaks -0.1022074 0.0823266 -1.241 0.21443 percentOfVoiceBreaks -0.0590097 1.2580202 -0.047 0.96259 meanIntensity -0.0765124 0.0612008 -1.250 0.21123 minimumIntensity 0.1037980 0.0331899 3.127 0.00176 ** maximumIntensity -0.0389995 0.0430368 -0.906 0.36484 ratioIntensity-2.0329346 1.2420286 -1.637 0.10168 noSyllsIntensity 0.1157678 0.0947699 1.222 0.22187 startSpeech0.0155578 0.1343117 0.116 0.90778 speakingRate -0.2583315 0.1648337 -1.567 0.11706 --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for binomial family taken to be 1) Null deviance: 2462.3 on 4310 degrees of freedom Residual deviance: 2209.5 on 4287 degrees of freedom AIC: 2257.5 Number of Fisher Scoring iterations: 6 I have seen models where almost all the features are showing one in a thousand significance but I accept that I could improve my model by normalizing some of the features (some are left skewed and I understand that I will get a better fir by taking their logs, for example). What really worries me is that the logistic function produces predictions that appear to fall well outside 0 to 1. If I make a dataset of the medians of the above features and use my logistic.model on it, it produces a figure of: x = predict(logistic.model, medians) x [1] 2.82959 which is well outside the range of 0 to 1. The actual distribution of all the predictions is: summary(pred) Min. 1st Qu. MedianMean 3rd Qu.Max. -1.516 2.121 2.720 2.731 3.341 6.387 I can get the model to give some sort of prediction by doing this: pred = predict(logistic.model, data) pred[pred = 1.5] = 0 pred[pred 1.5] = 1 t = table(pred, data[,24]) t pred 01 0 102 253 1 255 3701 classAgreement(t) $diag [1] 0.8821619 $kappa [1] 0.949 $rand [1] 0.7920472 $crand [1] 0.1913888 but as you can see I am using a break point well outside the range 0 to 1 and the kappa is rather low (I think). I am a bit of a novice in this, and the results worry me. Can anyone comment if the results look strange, or if they know I am doing something wrong? Stephen -- No virus found in this outgoing message. Checked by AVG Anti-Virus. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] turning labels into a vector
Try WBS.labs - as.vector(WBS) WBS.labs[1] substr(WBS.labs[1],1,1) etc.. Cheers, Simon. Hello 1/ 'priors' is a table looking like: W123 T678 S789 23 42 11 12 35 9 etc 2/ WBS - labels(priors) gives me a result of class list and length 1 looking like: W123 T678 S789 I want to read W123 into X[1] as W, T687 into X[2] as T and S789 into X[3] as S using substr(X[1],1,1) but I'm having trouble extracting each group of 4 digits from WBS Any help would be gratefully accepted. thanks Meredith __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] turning labels into a vector
Umm. Again... Try WBS.labs - as.vector(WBS[[1]]) WBS.labs[1] substr(WBS.labs[1],1,1) etc.. Cheers, Simon. Hello 1/ 'priors' is a table looking like: W123 T678 S789 23 42 11 12 35 9 etc 2/ WBS - labels(priors) gives me a result of class list and length 1 looking like: W123 T678 S789 I want to read W123 into X[1] as W, T687 into X[2] as T and S789 into X[3] as S using substr(X[1],1,1) but I'm having trouble extracting each group of 4 digits from WBS Any help would be gratefully accepted. thanks Meredith __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Double hurdle model in R
http://pscl.stanford.edu/content.html has code for poisson and negative binomial hurdle models, as well as for zip models. Package zicounts on CRAN may also be of interest. HTH, Simon. I am interested in utilizing this so called double hurdle model in my study. We can write the model in the following way: if (z'a + u 0 x'b + e 0) y = x'b + e, else y = 0 In the model, consumption y is the (left-) censored dependent variable. e and u are the normally distributed error terms. z'a is the participation equation and x'b is the expenditure equation. If we would remove the participation equation from the model we would end up with a type I tobit-model. In the tobit-model the same set of paprameters and variables determine both the discrete probability of non-zero outcome and the level of positive expenditure. In the double-hurdle-model we have separate parametrization of the participation and consumption decisions. I can estimate tobit-model using function survreg(). But what about this double hurdle thing? Has somebody written a R-function for estimating this sort of a model? Best regards, Kyösti Kurikka __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Iterative process for reading in text files
Do you mean that you have one file where the data is ordered by group? If so, using scan and a loop is an easy option. If you have the data in separate files by group, just make a character vector with your filenames as elements and use a loop to read the data using read.table at each iteration. Or if they all have similar names such as group1.txt, group2.txt, etc. you could construct the filenames on the fly using paste. HTH, Simon. Hello Instead of reading in group1.txt I want to read in groups1 for the first iteration of i, then groups2 for the second and so on. Obviously I can't use groups(i) but assume there is a way to do this. group-read.table(C:/Data/April 2005/group1.txt,header=T) thanks in advance Meredith __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] extracting correlations from nlme
Hi, I would like to know how (if) I can extract some of the information from the summary of my nlme. This is R. There is no if. Only how. at present, I get a summary looking something like this: summary(fit.nlme) [snip] I would like to extract the table of correlations, but have not been able to do it. ss - summary(fit.nlme) ss$corFixed Cheers, Simon. Any assistance, much appreciated. Cheers, Evelyn Evelyn Hall PhD Student Faculty of Veterinary Science The University of Sydney PMB 3 Camden, NSW, 2570 Phone: +61 2 9036 7736 Email: [EMAIL PROTECTED] [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R package that has (much) the same capabilities as SAS v9 PROC GENMOD
The questioner clearly wants generalized linear mixed models. lmer in package lme4 may be more appropriate. (Prof. Bates is a co-author.). glmmPQL should do the same job, though, but with less accuracy. Simon. check glm() On Apr 4, 2005 6:46 PM, William M. Grove [EMAIL PROTECTED] wrote: I need capabilities, for my data analysis, like the Pinheiro Bates S-Plus/R package nlme() but with binomial family and logit link. I need multiple crossed, possibly interacting fixed effects (age cohort of twin when entered study, sex of twin, sampling method used to acquire twin pair, and twin zygosity), a couple of random effects other than the cluster variable, and the ability to have a variable of the sort that PB call outer to the clustering variable---zygosity. Dependent variables are all parental (mom, dad separately of course) psychiatric diagnoses. In my data, twin pair ID is the clustering variable; correlations are expected to be exchangeable but substantially different between members of monozygotic twin pairs and members of dizygotic twin pairs. Hence, in my analyses, the variable that's outer to twin pair is monozygotic vs. dizygotic which of course applies to the whole pair. nlme() does all that but requires quasi-continuous responses, according to the preface/intro of PB's mixed models book and what I infer from online help (i.e., no family= or link= argument). The repeated() library by Lindsey seems to handle just one nested random effect, or so I believe I read while scanning backlogs of the R-Help list. glmmPQL() is in the ballpark of what I need, but once again seems to lack the outer variable specification that nlme() has, and which PROC GENMOD also has---and which I need. I read someplace of yags() that apparently uses GEE to estimate parameters of nonlinear models including GLIMs/mixed models, just the way PROC GENMOD (and many another program) does. But on trying to install it (either v4.0-1.zip or v4.0-2.tar.gz from Carey's site, or Ripley's Windows port) from a local, downloaded zip file (or tar.gz file converted to zip file), I always get an error saying: Error in file(file, r) : unable to open connection In addition: Warning message: cannot open file `YAGS/DESCRIPTION' with no obvious solution. So I can't really try it out to see if it does what I want. You may ask: Why not just use GENMOD and skip the R hassles? Because I want to embed the GLIM/mixed model analysis in a stratified resampling bootstrapping loop. Very easy to implement in R, moderately painful to do in SAS. Can anybody give me a lead, or some guidance, about getting this job done in R? Thanks in advance for your help. Regards, Will Grove | Iohannes Paulus PP. II, xxx Psychology Dept. | U. of Minnesota | -+ X-headers have PGP key info.; Call me at 612.625.1599 to verify key fingerprint before accepting signed mail as authentic! br x-sigsepp/x-sigsep Will Grovenbsp;nbsp;nbsp;nbsp;nbsp;nbsp; | Iohannes Paulus PP. II, xxx br Psychology Dept. |br U. of Minnesotanbsp; |br -+br br X-headers have PGP key info.; Call me at 612.625.1599 to verify key fingerprintbr before accepting signed mail as authentic!br br /body /html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- WenSui Liu, MS MA Senior Decision Support Analyst Division of Health Policy and Clinical Effectiveness Cincinnati Children Hospital Medical Center __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] LME correlation structures: user defined
Let me modify my question about user-defined covariance structures for LME models: Can somebody tell me how I can see the code for the definition of the correlation structures that come with the NLME package. Specifically I like to see the code for the functions coef, corMatrix, and intialize for any of the pre-defined correlation structures, and use this as a template to define a new correlation structure. So how do I see e.g. the code for the method initialize for the correlation structure corExp or corARMA? I'm interested in this too. Go to CRAN and download the source code for the nlme package. ungzip and untar it. Got to the R subdirectory. Inside that is a file called corStruct.R with all the method definitions for the built-in corStruct classes. Reading those should help. Cheers, Simon. thank you in advance! Michael Jerosch-Herold PS: Oh, and if somebody could still send me example code for a user defined correlation structure I would much appreciate it, as my previous requests for help have not yielded any response. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] (no subject)
I assume opened a data set means that you have attached a data frame. If you add new variables to the data frame (e.g. by transforming a variable in that data.frame), you will have to detach() and re- attach() it in order to get access to the variables without using the $ operator. I think this is in the manual somewhere. Cheers, Simon. Dear R I recently created some variables in R as in I opened a data set and then produced log base 10 transformations on some of the variables. When I ask R to do a simple x, y plot it recognises the raw data but does not recognise the log transformed variables. It says plot(logbrw, ParaSleep, type=n) Error in plot(logbrw, ParaSleep, type = n) : Object logbrw not found text(logbrw, ParaSleep, labels=row.names(sleep), cex=0.8, col=blue) Error in text(logbrw, ParaSleep, labels = row.names(sleep), cex = 0.8, : Object logbrw not found So do I have to somehow change the data for R to recognise the newly created variables? What should I do? brett __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Anova with Scheffe Tests
Hi, I don't think there are any packages on CRAN that implement Scheffe's test. If you don't mind using another multiple comparisons procedure, you could look at ?TukeyHSD and/or the multcomp package. Alternatively, you could write your own function to do Scheffe's test. At least one other person has done that. See the following post in the R-help archive http://finzi.psych.upenn.edu/R/Rhelp02a/archive/19393.html. I can't vouch for whether that person's function works properly, but it shouldn't be hard to hand-check it, and improve it. You could search R-help yourself and maybe come up with other solutions. Cheers, Simon. Hi R-people, I am wanting to run Factorial ANOVA followed by Scheffe tests on some spatial subjective data. I'm comparing X-Y independent coordinates against x-y dependent coordinates. There are only four independent spatial coordinates that form a square. I am wondering whether I am doing the right thing, because there doesn't seem to be a simple way of doing this. I have attempted to read `Practical regression and ANOVA using R' and am still confused. In good ol' Statview (now dearly departed) to complete a Scheffe test you selected the independent variables and dependent variable and it produced a table with the pairwise comparisons of the levels of the factor. I'm looking for a system that is as basic, but can be done using R and has documentation so I'm not guessing what I'm doing. I'd rather not have to do plots in R and then run over to dead software to do Scheffe's if possible. I checked on google and there seems to be code for a couple of functions out there, but I need something that has a manual. Is there a Scheffe function out there that is reasonably well documented, or should I consider some other method of dealing with this data. We have been using Scheffe for this type of analysis as I was under the impression it was very conservative. Tukey's HSD seems to be conservative as well. Should I try this? Is there a different approacch that is better and where can I read about it. Thanks for any help you can provide. Sam __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] Anova with Scheffe Tests
I stand corrected, although confidence.ellipse is a plotting function, and may not be quite what the questioner had in mind. Cheers, Simon. See confidence.ellipse() in the car() package. (Found from an R site search on Scheffe) -- Bert Gunter Genentech Non-Clinical Statistics South San Francisco, CA The business of the statistician is to catalyze the scientific learning process. - George E. P. Box -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Simon Blomberg Sent: Tuesday, March 01, 2005 4:25 PM To: [EMAIL PROTECTED]; R-help@stat.math.ethz.ch Subject: Re: [R] Anova with Scheffe Tests Hi, I don't think there are any packages on CRAN that implement Scheffe's test. If you don't mind using another multiple comparisons procedure, you could look at ?TukeyHSD and/or the multcomp package. Alternatively, you could write your own function to do Scheffe's test. At least one other person has done that. See the following post in the R-help archive http://finzi.psych.upenn.edu/R/Rhelp02a/archive/19393.html. I can't vouch for whether that person's function works properly, but it shouldn't be hard to hand-check it, and improve it. You could search R-help yourself and maybe come up with other solutions. Cheers, Simon. Hi R-people, I am wanting to run Factorial ANOVA followed by Scheffe tests on some spatial subjective data. I'm comparing X-Y independent coordinates against x-y dependent coordinates. There are only four independent spatial coordinates that form a square. I am wondering whether I am doing the right thing, because there doesn't seem to be a simple way of doing this. I have attempted to read `Practical regression and ANOVA using R' and am still confused. In good ol' Statview (now dearly departed) to complete a Scheffe test you selected the independent variables and dependent variable and it produced a table with the pairwise comparisons of the levels of the factor. I'm looking for a system that is as basic, but can be done using R and has documentation so I'm not guessing what I'm doing. I'd rather not have to do plots in R and then run over to dead software to do Scheffe's if possible. I checked on google and there seems to be code for a couple of functions out there, but I need something that has a manual. Is there a Scheffe function out there that is reasonably well documented, or should I consider some other method of dealing with this data. We have been using Scheffe for this type of analysis as I was under the impression it was very conservative. Tukey's HSD seems to be conservative as well. Should I try this? Is there a different approacch that is better and where can I read about it. Thanks for any help you can provide. Sam __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] an R script editor for Mac
I'm surprised nobody has mentioned Alpha. It has highlighting, indenting, parenthesis matching, excellent integration with R (or that other commercial version of R). There are versions for Classic (Alpha8), OS X (Alphax), as well as *NIX and Windows (AlphatTk). Alpha is shareware, based on the open source AlphTcl library. http://alphatcl.sourceforge.net/wikit/ Cheers, Simon. On Jan 21, 2005, at 5:35 AM, Jacques VESLOT wrote: Dear all, Could someone please make me know if there is a nice script editor available under Mac, similar to Crimson, that offers R syntax highlighting (and pairs of parentheses underlining) ? Thanks in advance, Jacques VESLOT Cirad __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Calculate Mean of Column Vectors?
z - apply(y, 2, mean) Cheers, Simon. Hello, I've got an array defined as y - rnorm(3000), dim(y) - c(3, 1000). I'd like to produce a 1000-element vector z that is the mean of the corresponding elements of y (like z[1,1] - mean(y[1,1], y[2,1], y[3,1])), but being new to R, I'm not sure how to do this for all elements at once (or, at least, simply). Any help is appreciated. Thanks, Tom __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. Visiting Fellow School of Botany Zoology The Australian National University Canberra ACT 0200 Australia T: +61 2 6125 8057 email: [EMAIL PROTECTED] F: +61 2 6125 5573 CRICOS Provider # 00120C [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] summary.manova and rank deficiency
Sorry for reposting. I didn't receive any of the replies to my original message (must be an email problem at my end). I will read the responses in the archive. Thanks! Simon. Simon Blomberg, PhD Depression Anxiety Consumer Research Unit Centre for Mental Health Research Australian National University http://www.anu.edu.au/cmhr/ [EMAIL PROTECTED] +61 (2) 6125 3379 -Original Message- From: Simon Blomberg Sent: Monday, 24 November 2003 9:09 AM To: [EMAIL PROTECTED] Subject: [R] summary.manova and rank deficiency Hi all, I have received the following error from summary.manova: Error in summary.manova(manova.test, test = Pillai) : residuals have rank 36 64 The data is simulated data for 64 variables. The design is a 2*2 factorial with 10 replicates per treatment. Looking at the code for summary.manova, the error involves a problem with qr(). Does anyone have a suggestion as to how to deal with this error? The analysis works fine when there are only 32 variables. Thanks, Simon. Simon Blomberg, PhD Depression Anxiety Consumer Research Unit Centre for Mental Health Research Australian National University http://www.anu.edu.au/cmhr/ [EMAIL PROTECTED] +61 (2) 6125 3379 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
