Re: [R] D(dnorm...)?
Hi Spencer, I think if you have a problem that needs a lot of symbolic manipulation you are probably better off driving it directly from something like Maple or Mathematica (I prefer maple, actually) than trying to drive it from R. It just gets too clumsy. On the other hand it is very handy having a simple differentiator available in R, like D(...) for small jobs that are not worth taking to a big system like Maple. The point I was trying to make in the previous message was that with a little thought it could be made a lot more convenient. This arose in connexion with a real problem. We needed to differentiate a pdf that had the normal density function in it, but was otherwise quite simple and we had to hack the code in another system (not unlike R, as it happens) to handle it. The hack was quite small and it became clear that with a slight change of design users would not need to do hacks like this for simple extensions such as the one we needed. As it was a hack, we only put it in for the standard density and I suspect that is the reason why even now the derivative tables in both R and the other system (not unlike R) only handle normal density and distribution funcitons in one variable. I'm sort of avoiding your question, because I don't know how hard it would be to link R with Yacas, either, but if you really wanted to go that way I see that Yacas can be driven over the net via a Java applet. Something like this might provide the simplest link, if not the most efficient. But note that Yacas uses the eccentric Mathematica notation, where the functions are Capitalised, for example, as nouns are in German. That's a small bother you could do without, too. Regards, Bill Venables. -Original Message- From: Spencer Graves [mailto:[EMAIL PROTECTED] Sent: Thursday, 26 January 2006 2:14 PM To: Venables, Bill (CMIS, Cleveland) Cc: [EMAIL PROTECTED]; r-help@stat.math.ethz.ch; [EMAIL PROTECTED] Subject: Re: [R] D(dnorm...)? Hello, Bill: I'm not qualified to make this suggestion since I'm incapable of turning it into reality, but what about creating a link between R and one of the Mathematica clones like Yacas? I can immagine that it could be substantially more difficult than linking R to other software like Excel, but ... . Spencer Graves [EMAIL PROTECTED] wrote: Yes Bert, this time you are missing something (unusually) ... As Brian Ripley pointed out 'dnorm' is in the derivative table, *but* only as a function of one variable. So if you want to find the derivative of dnorm(x, mean, sigma) you have to write it as 1/sigma * dnorm((x - mu)/sigma). Here is a little example: D(Quote(pnorm((x-mu)/sigma)), x) dnorm((x - mu)/sigma) * (1/sigma) D(D(Quote(pnorm((x-mu)/sigma)), x), mu) (x - mu)/sigma * (dnorm((x - mu)/sigma) * (1/sigma)) * (1/sigma) --- Like Brian, I recall the suggestion that we make D(...) extensible. I still think it is a good idea and worth considering. Under one scheme you would specify an object such as Fnorm - structure(quote(pnorm(x, mu, sigma)), deriv = list(x = Quote(dnorm(x, mu, sigma)/sigms), mu = Quote(-dnorm(x, mu, sigma)/sigma), sigma = Quote(-(x - mu)*dnorm(x, mu, sigma)/sigma^2), class = dfunction) ane write a generic differentiate function with a dfunction method and D as the default. I don't think it's quite that easy, but the plan is clear enough. Bill. -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Berton Gunter Sent: Thursday, 26 January 2006 8:58 AM To: 'Spencer Graves'; r-help@stat.math.ethz.ch Subject: Re: [R] D(dnorm...)? dnorm() is an internal function, so I don't see how D (or deriv) can do anything with it symbolically. Am I missing something? -- Bert -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Spencer Graves Sent: Wednesday, January 25, 2006 2:43 PM To: r-help@stat.math.ethz.ch Subject: [R] D(dnorm...)? Can someone help me understand the following: D(expression(dnorm(x, mean)), mean) [1] 0 sessionInfo() R version 2.2.1, 2005-12-20, i386-pc-mingw32 attached base packages: [1] methods stats graphics grDevices utils datasets [7] base By my computations, this should be something like ((mean-x)/sd^2)*dnorm(...). Thanks for your help. Spencer Graves __ R-help@stat.math.ethz.ch 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
[R] question about system command
Hi R Users
I am going to write a very short script for my small pilot simulation study.
I need to call a DOS program with different input data files in the middle
of for loop.
There are 100 input data files (e.g., input001.dat, input002.dat, ,
input100.dat)
for (1 in i :100)
{
system('C:\\Program Files\\DOSPROGRAM\\RUN.exe input001.dat')
}
I like to change the input data file name as for loop runs. How can I do
that ?
Thanks in advance
TM
__
R-help@stat.math.ethz.ch 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] question about system command
Something like this.
See ?formatC, ?paste
for (i in 1:100)
{
system(paste(C:\\Progra~1\\DOSPROGRAM\\RUN.exe input,
formatC(i, digits=2, flag='0'),
'.dat',
sep=''))
}
--Matt
Statistician
Amgen, Inc
-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] Behalf Of Taka Matzmoto
Sent: Wednesday, January 25, 2006 9:31 PM
To: r-help@stat.math.ethz.ch
Subject: [R] question about system command
Hi R Users
I am going to write a very short script for my small pilot simulation study.
I need to call a DOS program with different input data files in the middle
of for loop.
There are 100 input data files (e.g., input001.dat, input002.dat, ,
input100.dat)
for (1 in i :100)
{
system('C:\\Program Files\\DOSPROGRAM\\RUN.exe input001.dat')
}
I like to change the input data file name as for loop runs. How can I do
that ?
Thanks in advance
TM
__
R-help@stat.math.ethz.ch 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] question about system command
It is a little easier to use sprintf, as in
system(sprintf(%s input%03d.dat,C:\\Progra~1\\DOSPROGRAM\\RUN.exe, i))
With formatC it is not obvious why digits=2 is needed, but it would be
more obvious why width=3 would work.
On Wed, 25 Jan 2006, Austin, Matt wrote:
Something like this.
See ?formatC, ?paste
for (i in 1:100)
{
system(paste(C:\\Progra~1\\DOSPROGRAM\\RUN.exe input,
formatC(i, digits=2, flag='0'),
'.dat',
sep=''))
}
--Matt
Statistician
Amgen, Inc
-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] Behalf Of Taka Matzmoto
Sent: Wednesday, January 25, 2006 9:31 PM
To: r-help@stat.math.ethz.ch
Subject: [R] question about system command
Hi R Users
I am going to write a very short script for my small pilot simulation study.
I need to call a DOS program with different input data files in the middle
of for loop.
There are 100 input data files (e.g., input001.dat, input002.dat, ,
input100.dat)
for (1 in i :100)
{
system('C:\\Program Files\\DOSPROGRAM\\RUN.exe input001.dat')
}
I like to change the input data file name as for loop runs. How can I do
that ?
Thanks in advance
TM
__
R-help@stat.math.ethz.ch 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
--
Brian D. Ripley, [EMAIL PROTECTED]
Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel: +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax: +44 1865 272595
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] [R-SIG-Mac] Hist for different levels of a factor
The list of your interest is R-help not R-sig-mac stefano Il giorno 26/gen/06, alle ore 01:20, Sylvain Charlat ha scritto: Hi, Is there any simple way to get histogram for different levels of factor? Say you have the following data set: Island Sp.diam Moorea 1.21 Moorea 1.27 Moorea 1.28 Moorea 1.22 Moorea 1.28 Rurutu 1.5 Rurutu 1.67 Rurutu 1.75 Rurutu 1.55 Rurutu 1.7 Rurutu 1.55 Rurutu 1.59 Rurutu 1.66 Rurutu 1.7 Is there anything better than: hist(Sp.diam[factor(Island)==Moorea],main=Island=Moorea) followed by hist(Sp.diam[factor(Island)==Rurutu],main=Island=Rurutu) Sorry for this very stupid and basic question... Thanks for any help, Sylvain. --- Sylvain Charlat [EMAIL PROTECTED] / [EMAIL PROTECTED] http://www.ucl.ac.uk/~ucbtghu/Sylvain.htm - Department of Biology, University College London 4 Stephenson Way, London, NW1 2HE, UK - Fieldwork address: Gump Station, University of California Berkeley BP 244 Maharepa, 98728 Moorea, French Polynesia - Phone (Moorea): (+689) 56 43 97 or 26 90 18 or 56 52 87 - Fax (Moorea): (+689) 56 32 72 - Phone / Fax (London): (+44)20 76795072 / +(44) 20 76795052 --- __ _ Nouveau : téléphonez moins cher avec Yahoo! Messenger ! Découvez les tarifs exceptionnels pour appeler la France et l'international. Téléchargez sur http://fr.messenger.yahoo.com ___ R-SIG-Mac mailing list R-SIG-Mac@stat.math.ethz.ch https://stat.ethz.ch/mailman/listinfo/r-sig-mac __ R-help@stat.math.ethz.ch 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] construct a bundle, subdirs do not exist?
Hi, Sorry to bother, but I checked around and did not succed creating a bundle from six existing packages (which are checkable, installable, etc. individually). I carefully followed the procedure given in ch. 1.1.5 Package bundles. However, I am getting [EMAIL PROTECTED]:~/R/Sources R CMD check cwhmisc * checking for working latex ... OK * using log directory '/home/woodstock/hoffmann/R/Sources/cwhmisc.Rcheck' * using R version 2.2.1, 2005-12-20 * checking for file 'cwhmisc/DESCRIPTION' ... OK * this is package '' version '1.0.0' * checking if this is a source package ... OK /usr/local/lib/R/bin/INSTALL: /home/woodstock/hoffmann/R/Sources/cwhmisc/cwhmath,: does not exist /usr/local/lib/R/bin/INSTALL: /home/woodstock/hoffmann/R/Sources/cwhmisc/cwhstring,: does not exist /usr/local/lib/R/bin/INSTALL: /home/woodstock/hoffmann/R/Sources/cwhmisc/cwhstat,: does not exist /usr/local/lib/R/bin/INSTALL: /home/woodstock/hoffmann/R/Sources/cwhmisc/cwhplot,: does not exist /usr/local/lib/R/bin/INSTALL: /home/woodstock/hoffmann/R/Sources/cwhmisc/cwhprint,: does not exist /usr/local/lib/R/bin/INSTALL: /home/woodstock/hoffmann/R/Sources/cwhmisc/cwhool: does not exist /usr/local/lib/R/bin/INSTALL: null directory ERROR Installation failed. *But I have* [EMAIL PROTECTED]:~/R/Sources/cwhmisc ls -lFa total 144 drwxr-xr-x 8 hoffmann wsl 1024 Jan 26 07:51 ./ drwxrwxrwx 26 hoffmann wsl 3072 Jan 26 08:29 ../ -rwxr--r-- 1 hoffmann wsl 430 Jan 26 07:51 DESCRIPTION* drwxr-xr-x 4 hoffmann wsl 1024 Jan 26 07:03 cwhmath/ drwxr-xr-x 5 hoffmann wsl 1024 Jan 26 07:04 cwhplot/ drwxr-xr-x 4 hoffmann wsl 1024 Jan 26 07:04 cwhprint/ drwxr-xr-x 4 hoffmann wsl 1024 Jan 26 07:04 cwhstat/ drwxr-xr-x 4 hoffmann wsl 1024 Jan 26 07:04 cwhstring/ drwxr-xr-x 4 hoffmann wsl 1024 Jan 26 07:13 cwhtool/ i.e. the subdirectories are *present*! And the packages with their [EMAIL PROTECTED]:~/1Sources/cwhmisc/cwhmath ls -lFa total 96 drwxr-xr-x 4 hoffmann wsl 1024 Jan 26 07:03 ./ drwxr-xr-x 8 hoffmann wsl 1024 Jan 26 07:51 ../ -rwxr--r-- 1 hoffmann wsl 142 Jan 25 16:24 DESCRIPTION.in* -rwxr--r-- 1 hoffmann wsl 1121 Jan 20 14:58 INDEX* drwxr-xr-x 3 hoffmann wsl 1024 Jan 25 16:21 R/ drwxr-xr-x 2 hoffmann wsl 1024 Jan 25 16:21 man/ [EMAIL PROTECTED]:~/1Sources/cwhmisc/cwhmath cat DESCRIPTION.in Package: CWHmath Description: Miscellaneous mathematical functions for general use Title: Miscellaneous mathematical functions for general use are fully present. Is this information enough to find the flaw in my bundle construction? Thank for help Christian -- Dr. Christian W. Hoffmann, Swiss Federal Research Institute WSL Mathematics + Statistical Computing Zuercherstrasse 111 CH-8903 Birmensdorf, Switzerland Tel +41-44-7392-277 (office) -111(exchange) Fax +41-44-7392-215 (fax) [EMAIL PROTECTED] http://www.wsl.ch/staff/christian.hoffmann International Conference 5.-7.6.2006 Ekaterinburg Russia Climate changes and their impact on boreal and temperate forests http://ecoinf.uran.ru/conference/ __ R-help@stat.math.ethz.ch 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] D(dnorm...)?
While symbolic computation is handy, I actually think a more pressing addition to R is some kind of automatic differentiation facility, particularly 'reverse mode' AD, which can be spectacular. There are free tools available for it as well, though I don't know how well developed they are. See: http://www-unix.mcs.anl.gov/autodiff/AD_Tools/ I admit this is not quite the same thing, but for statistical computations this is, in my experience, the key thing you need. (Well, for frequentist estimation at any rate...) There are commercial systems that use this idea already, of course. Two that I know of are 'ADMB' (and its associated 'ADMB-RE' for random effects) estimation and of course the 'S-NUOPT' module for another system not unlike R. ADMB is, frankly, difficult to use but it performs so well and so quickly once you get it going nothing else seems to come close to it. I has become almost a de-facto standard at the higher end of the fishery stock assessment game, for example, where they are always fitting huge, highly complex and very non-linear models. Bill V. -Original Message- From: Berwin A Turlach [mailto:[EMAIL PROTECTED] On Behalf Of Berwin A Turlach Sent: Thursday, 26 January 2006 4:50 PM To: Spencer Graves Cc: Venables, Bill (CMIS, Cleveland); r-help@stat.math.ethz.ch; [EMAIL PROTECTED]; [EMAIL PROTECTED] Subject: Re: [R] D(dnorm...)? G'day Spencer, SG == Spencer Graves [EMAIL PROTECTED] writes: SG I'm not qualified to make this suggestion since I'm incapable SG of turning it into reality, [...] This statement applies to me too, nevertheless I would like to point out the following GPL library: http://www.gnu.org/software/libmatheval/ I am wondering since some time how hard it would be to incorporate that library into R and let it provide symbolic differentiation capabilities for R. Cheers, Berwin __ R-help@stat.math.ethz.ch 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] read.table problem
Dear R useRs,
I have big (23000 rows), vertical bar delimited file:
e.g.
A1|Text a,Text b, Text c|345
A2|Text bla|456
...
..
.
Try using
A - read.table('filename.txt', header=FALSE,sep='\|')
process stop at line 11975 with warning message:
number of items read is not a multiple of the number of columns
I have no problems with processing similar file, which is only 1
rows long?
Any suggestion what's the problem here. Thank's in advance.
Cheers, Andrej
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] reducing learning curves?
Hi all, I am really new to the R language. I am a long time Matlab and C++ user and I was forced to learn R because I am taking a statistics class. I am seeking to reduce the learning curve to as smooth as possible. Are there any addon/plug-in features that can reduce the learning curve, for example, the following features can be very helpful for new learners: 1. Matlab-like command line auto-completion: Matlab has huge amount of command and nobody is able to remember them off the head. So a nice feature of Matlab command line is that I just need to type the first a few letters and then I press TAB key, there will be a list of possible commands popping up so I just need to select one. This helps a lot in terms of learning for new comers. A more advanced command auto-completion is Visual C++-like, which is implemented in program editor. It helps a lot while doing programming; 2. A good IDE editor with embedded inline debugger: can be as good as VC++, but also can be as simple as Matlab's debugger, which can breakpoint and trace line-by-line... the editor can do syntax correction, syntax check, syntax highlighting, code formatting, etc. Could you please recommend some good addon/plugins that have the above features? Could you please also suggest some tips/tools/tricks that can help me reduce the learning curve? Thank you very much! Michael. [[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] read.table problem
Andrej Kastrin [EMAIL PROTECTED] writes:
Dear R useRs,
I have big (23000 rows), vertical bar delimited file:
e.g.
A1|Text a,Text b, Text c|345
A2|Text bla|456
...
..
.
Try using
A - read.table('filename.txt', header=FALSE,sep='\|')
process stop at line 11975 with warning message:
number of items read is not a multiple of the number of columns
I have no problems with processing similar file, which is only 1
rows long?
Any suggestion what's the problem here. Thank's in advance.
Well, the most obvious suspect is the file...
You might try
table(count.fields('filename.txt',sep='|'))
(There's no point in escaping the vertical bar)
Also, beware of quote symbols in the file, and possibly consider using
fill=TRUE.
--
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
Re: [R] Can R handle medium and large size data sets?
Martin == Martin Lam [EMAIL PROTECTED]
on Tue, 24 Jan 2006 12:13:07 -0800 (PST) writes:
Martin Dear Gueorgui,
Is it true that R generally cannot handle medium sized
data sets(a couple of hundreds of thousands observations)
and threrefore large date set(couple of millions of
observations)?
Martin It depends on what you want to do with the data sets.
Martin Loading the data sets shouldn't be any problem I
Martin think. But using the data sets for analysis using self
Martin written R code can get (very) slow, since R is an
Martin interpreted language (correct me if I'm wrong).
(Since you asked for it ;-) )
Yes, you are wrong to quite some extent (you are partially
right, too): Of course one *can* write ``self written R code''
that is very slow, and yes, we have seen such code more than
once. However, 98% of the problems {never trust a statistic
unless you mad it up ... :-) :-) } are relatively easily
solvable very efficiently with R.
You are right that it is easier to write slow code in an
interpreted language than in a compiled one.
E.g., not making use of vectorized operations in R is one famous
recipe to produce slow code pretty successfully ...
Martin To increase speed you will often need to experiment with
Martin the R code. For example, what I've noticed is that
Martin processing data sets as matrices works much faster
Martin than data.frame().
yes, indeed; see also the other answers to Gueorgui's question.
Martin Writing your code in C(++), compile it and include
Martin it in your R code is often the best way.
Martin HTH,
Martin 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
Re: [R] R-help Digest, Vol 35, Issue 24
I've been reluctant to step into this topic, but now feel that it may be helpful to make a certain point. On the internet, for the most part, the person behind the email is invisible and intangible. It is therefore possible, when someone puts their foot down, to stamp inadvertently on someone else's already broken toes. A friend of mine, very intelligent, very knowledgeable and creative, very articulate, nevertheless when writing uses spelling which can be a close approximation to random, and some interesting variants of grammar and vocabulary as well. The reason: dyslexia. While most of us hit the wrong keys at times (and when we read back over what we've written tend to see what we intended to write rather than what we did write), and when backed against the wall would admit that we could have got it right if we had paid better attention, there are some people who can't help getting it wrong. But, on the internet, one cannot readily recognise who they are (though in some cases, if one knows the signs, one may guess). Best wishes to all, Ted. E-Mail: (Ted Harding) [EMAIL PROTECTED] Fax-to-email: +44 (0)870 094 0861 Date: 25-Jan-06 Time: 10:06:35 -- 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
[R] documentation mistake
dear r-help. i want top report a mistake in the documentation help(splineDesign). the last sentence of value is Each B-spline is defined by a set of 'ord' successive knots so the total number of B-splines is 'length(knots)-ord'. it is not correct! one b-spline is defines by a set of 'ord+1' successive knots. take an easy example: degree 0 is a constant spline it has order 1 and you need 2 knots to define the boundaries do i have to report to somebody else? yours sincerely, michael kobl (dipl. math.) [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] How to use rfm.test ? (Package MarkedPointProcess)
I would like to compute the MC test (rfm.test) available in the package MarkedPointProcess (for the data BITOEK for example) in order to test the dependence between the marks and their locations. Why the syntax of rfm.test is false here? I have the message : ** ML WARNING! Forbidden values! -- if there are too many warnings try narrower lower and upper bounds for the variables. -9.802347e-17 1364.372 [ 0 0.01808660 , 9.990131 1364.372 ] ** The program : - # Packages # library(spatstat) library(RandomFields) library(adapt) library(MarkedPointProcess) # Program # data(BITOEK) win - ripras(steigerwald$coord) PointProcess - ppp(x=steigerwald$coord[,1],y=steigerwald$coord[,2],window=win,marks=steigerwald$diam) plot(PointProcess) rfm.test(coord=steigerwald$coord,steigerwald$diam,MCrepetitions=19) - regards. Florent Bonneu Laboratoire de Statistique et Probabilités bureau 148 bât. 1R2 Université Toulouse 3 118 route de Narbonne - 31062 Toulouse cedex 9 [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] Question about Aggregate
hello, Suppose you a monthly series you want to aggregate at a quaterly frequency with the start and the end of your series in the middle of the quarter. For example 2001M2 2001M3 2001M4 2001M5 2001M6 2001M7 12 13 12 14 16 15 how can you get something like : 2001Q1 2001Q2 200Q3 NA 14 NA or 2001Q1 2001Q2 200Q3 12.5 14 15 Thanks in advance Matthieu [[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] reducing learning curves?
Hello Michael,
you might want to utilise Emacs/ESS. ESS provides auto-completion by using
TAB for a process buffer '*R*' and C-cTAB for a source file '*.R'
(ess-mode).
As far as debugging is concerned, R offers:
?browser
?debug
?trace
for example. Additionally, there is a CRAN package named 'debug' available
and an article in RNews:
Mark Bravington. Debugging without (too many) tears. R News, 3(3):29-32,
December 2003, about it.
HTH,
Bernhard
-Ursprüngliche Nachricht-
Von: Michael [mailto:[EMAIL PROTECTED]
Gesendet: Mittwoch, 25. Januar 2006 09:10
An: R-help@stat.math.ethz.ch
Betreff: [R] reducing learning curves?
Hi all,
I am really new to the R language. I am a long time Matlab and C++ user and
I was forced to learn R because I am taking a statistics class.
I am seeking to reduce the learning curve to as smooth as possible.
Are there any addon/plug-in features that can reduce the learning curve, for
example, the following features can be very helpful for new learners:
1. Matlab-like command line auto-completion: Matlab has huge amount of
command and nobody is able to remember them off the head. So a nice feature
of Matlab command line is that I just need to type the first a few letters
and then I press TAB key, there will be a list of possible commands
popping up so I just need to select one. This helps a lot in terms of
learning for new comers. A more advanced command auto-completion is Visual
C++-like, which is implemented in program editor. It helps a lot while doing
programming;
2. A good IDE editor with embedded inline debugger: can be as good as VC++,
but also can be as simple as Matlab's debugger, which can breakpoint and
trace line-by-line... the editor can do syntax correction, syntax check,
syntax highlighting, code formatting, etc.
Could you please recommend some good addon/plugins that have the above
features?
Could you please also suggest some tips/tools/tricks that can help me reduce
the learning curve?
Thank you very much!
Michael.
[[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
*
Confidentiality Note: The information contained in this mess...{{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] Interpolating spline problems and akima
Hi everyone
I was using spline to interpolate single or two consecutive missing data points
in time series. However, when it comes to longer gaps in the data the spline
function generate new data for both my known and unknown data (see below).
Aside from not understanding why this happens, I thought thought I might try
function aspline in library (akima). However, I cannot install or find this
library/package.
Has anyone any instructions to install this package or help with my problem?
Using this works
day-c(1,2,3,4,5,6,7,8,9,10,11,12,13)
datfile-c(16,23,39,37,45,34,NA,3,7,8,15,20,21)
for (i in 1:13) {
+ dat-datfile
+ dat.ok - !is.na(dat)
+ dat2-dat[dat.ok]
+ day2-day[dat.ok]
+
+ dat.interp-spline(day2,dat2,length(day))
+ outt-rep(NA,length(day))
+ for (j in 1:length(dat.interp$x)) {
+ outt[j]-dat.interp$y[j]
+ } }
but this does not?
day-c(1,2,3,4,5,6,7,8,9,10)
datfile-c(17,6,29,22,90,72,20,4,4,NA)
for (i in 1:10) {
+ dat-datfile
+ dat.ok - !is.na(dat)
+ dat2-dat[dat.ok]
+ day2-day[dat.ok]
+
+ dat.interp-spline(day2,dat2,length(day))
+ outl-rep(NA,length(day))
+ for (j in 1:length(dat.interp$x)) {
+ outl[j]-dat.interp$y[j]
+ } }
...
Dr Tom C Cameronoffice: 0113 34 32837 (10.23 Miall)
Ecology Evolution Res. Group. lab: 0113 34 32884 (10.20 Miall)
School of Biological Sciences Mobile: 07966160266
University of Leeds email: [EMAIL PROTECTED]
Leeds LS2 9JT
LS2 9JT
__
R-help@stat.math.ethz.ch 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] persp() and character labels for axis
I want to plot dates on the y-axis of a persp() plot. persp(x=1:30,y=days,y=yields) axis(2, 1:5, LETTERS[1:5]) does not work. On the mailinglist I found old messages, that said, that text() does not apply (yet) for 3-d plots. And the same question ( http://www.r-project.org/nocvs/mail/r-help/2002/9009.html ) was not answered. Perhaps this time? ;) __ R-help@stat.math.ethz.ch 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] persp() and character labels for axis
Thomas Steiner wrote: I want to plot dates on the y-axis of a persp() plot. persp(x=1:30,y=days,y=yields) axis(2, 1:5, LETTERS[1:5]) does not work. On the mailinglist I found old messages, that said, that text() does not apply (yet) for 3-d plots. And the same question ?persp points you to ?trans3d which is useful to calculate coordinates for calls to 2D functions such as text(). Uwe Ligges ( http://www.r-project.org/nocvs/mail/r-help/2002/9009.html ) was not answered. Perhaps this time? ;) __ R-help@stat.math.ethz.ch 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 use rfm.test ? (Package MarkedPointProcess)
Dear Florent Bonneu, the optim algorithm with parameter method=L-BFGS-B (used in rfm.test) does not stay always exactly within the given bounds during the search of the optimum. This happens more frequently when the bounds are too wide. Here. rfm.test notices that -9.802347e-17 is less then the given lower bound of value 0. If the messages appear frequently it appeared to me that the algorithms runs into boundary, local minima more frequently. So, if the messages appear only a very few times, ignore them. If they appear several times, use closer bounds. Cheers, Martin Florent Bonneu wrote: I would like to compute the MC test (rfm.test) available in the package MarkedPointProcess (for the data BITOEK for example) in order to test the dependence between the marks and their locations. Why the syntax of rfm.test is false here? I have the message : ** ML WARNING! Forbidden values! -- if there are too many warnings try narrower lower and upper bounds for the variables. -9.802347e-17 1364.372 [ 0 0.01808660 , 9.990131 1364.372 ] ** The program : - # Packages # library(spatstat) library(RandomFields) library(adapt) library(MarkedPointProcess) # Program # data(BITOEK) win - ripras(steigerwald$coord) PointProcess - ppp(x=steigerwald$coord[,1],y=steigerwald$coord[,2],window=win,marks=steigerwald$diam) plot(PointProcess) rfm.test(coord=steigerwald$coord,steigerwald$diam,MCrepetitions=19) - regards. Florent Bonneu Laboratoire de Statistique et Probabilités bureau 148 bât. 1R2 Université Toulouse 3 118 route de Narbonne - 31062 Toulouse cedex 9 [EMAIL PROTECTED] -- Prof. Dr. habil. Martin Schlather Helmut-Schmidt-Universität Universität der Bundeswehr Hamburg Fachbereich WOW Holstenhofweg 85 D -- 22043 Hamburg email: [EMAIL PROTECTED] http://www2.unibw-hamburg.de/schlath phone: +49 (0)40 6541 2873 Fax : +49 (0)40 6541 2023 __ R-help@stat.math.ethz.ch 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 ANCOVA: still confused
OK, we're getting somewhere. First, it looks as though (by the error message) that you have a big dataset. My first recommendation is to use lmer instead of lme, you will see a significant benefit in terms of computional speed. For the model this would be lmer(rtot ~ sexv +(purban|box:chick) + (purban|box), bb, na.action=na.omit) Now, you have run out of memory. I don't know what operating system you are using, so go and see the appropriate FAQ for increasing memory for your OS. Second, I made a mistake in my reply. Your random statement should be random=~purban|box/chick denoting that chicks are nested in boxes, not boxes nested in chicks, sorry about that. Now, why is it that each chick within box has the same value for purban? If this is so, why are you fitting that as a random effect? It seems not to vary across individual chicks, right? It seems there is only an effect of box and not an effect for chicks. Why not just fit a random effect only for box such as: rtot.lme - lme(fixed=rtot~sexv, random=~purban2|box, na.action=na.omit,bb) or in lmer lmer(rtot ~ sexv + (purban|box), bb, na.action=na.omit) Harold -Original Message- From: Jeffrey Stratford [mailto:[EMAIL PROTECTED] Sent: Tue 1/24/2006 8:57 PM To: Doran, Harold; r-help@stat.math.ethz.ch Cc: Subject:RE: [R] nested ANCOVA: still confused R-users and Harold. First, thanks for the advice; I'm almost there. The code I'm using now is library(nlme) bb - read.csv(E:\\eabl_feather04.csv, header=TRUE) bb$sexv - factor(bb$sexv) rtot.lme - lme(fixed=rtot~sexv, random=~purban2|chick/box, na.action=na.omit, data=bb) A sample of the data looks like this box chick rtotpurban2 sexv 1 1 6333.51 0.026846f 1 2 8710.8840.026846m 2 1 5810.0070.161074f 2 2 5524.33 0.161074f 2 3 4824.4740.161074f 2 4 5617.6410.161074f 2 5 6761.7240.161074f 4 1 7569.6730.208054m 4 2 7877.0810.208054m 4 4 7455.55 0.208054f 7 1 5408.2870.436242m 10 1 6991.7270.14094 f 12 1 8590.2070.134228f 12 2 7536.7470.134228m 12 3 5145.3420.134228m 12 4 6853.6280.134228f 15 1 8048.7170.033557m 15 2 7062.1960.033557m 15 3 8165.9530.033557m 15 4 8348.58 0.033557m 16 2 6534.7750.751678m 16 3 7468.8270.751678m 16 4 5907.3380.751678f 21 1 7761.9830.221477m 21 2 6634.1150.221477m 21 3 6982.9230.221477m 21 4 7464.0750.221477m 22 1 6756.7330.281879f 23 2 8231.4960.134228m The error I'm getting is Error in logLik.lmeStructInt(lmeSt, lmePars) : Calloc could not allocate (590465568 of 8) memory In addition: Warning messages: 1: Fewer observations than random effects in all level 2 groups in: lme.formula(fixed = rtot ~ sexv, random = ~purban2 | chick/box, 2: Reached total allocation of 382Mb: see help(memory.size) There's nothing special about chick 1, 2, etc. These were simply the order of the birds measured in each box so chick 1 in box 1 has nothing to do with chick 1 in box 2. Many thanks, Jeff Jeffrey A. Stratford, Ph.D. Postdoctoral Associate 331 Funchess Hall Department of Biological Sciences Auburn University Auburn, AL 36849 334-329-9198 FAX 334-844-9234 http://www.auburn.edu/~stratja Doran, Harold [EMAIL PROTECTED] 01/24/06 2:04 PM Dear Jeff: I see the issues in your code and have provided what I think will solve your problem. It is often much easier to get help on this list when you provide a small bit of data that can be replicated and you state what the error messages are that you are receiving. OK, with that said, here is what I see. First, you do not need to use the syntax bb$sex in your model, this can be significantly simplified. Second, you do not have a random statement in your model. Here is your original model: lme(bb$rtot~bb$sex, bb$purban|bb$chick/bb$box, na.action=na.omit) Here is what it should be: lme(fixed = rtot~sex, random=~purban|chick/box, na.action=na.omit, data=bb) Notice there is a fixed and random call. You can simplify this as lme(rtot~sex, random=~purban|chick/box, na.action=na.omit, bb) Note, you can eliminate the fixed= portion but not the random statement. Last, if you want to do this in lmer, the newer function for mixed models in the Matrix package,
Re: [R] Question about Aggregate
Try this: library(zoo) z - zooreg(c(12, 13, 12, 14, 16, 15), start = c(2001, 2), freq = 12) aggregate(z, trunc(4 * time(z))/4, mean) 2001(1) 2001(2) 2001(3) 12.514.015.0 On 1/25/06, Matthieu Cornec [EMAIL PROTECTED] wrote: hello, Suppose you a monthly series you want to aggregate at a quaterly frequency with the start and the end of your series in the middle of the quarter. For example 2001M2 2001M3 2001M4 2001M5 2001M6 2001M7 12 13 12 14 16 15 how can you get something like : 2001Q1 2001Q2 200Q3 NA 14 NA or 2001Q1 2001Q2 200Q3 12.5 14 15 Thanks in advance Matthieu [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] persp() and character labels for axis
I want to plot dates on the y-axis of a persp() plot. ?persp points you to ?trans3d which is useful to calculate coordinates for calls to 2D functions such as text(). ?trans3d gives this: No documentation for 'trans3d' in specified packages and libraries: you could try 'help.search(trans3d)' Even loading grDevices and graphics libraries did not help. http://finzi.psych.upenn.edu/R/library/grDevices/html/trans3d.html seems to be what I should look for, but there it says noting about text either (I want to label axis! ,something like t=1:10 plot(x=t,y=t^2,type=l) axis(1,labels=LETTERS[1:5]) and more proper...) Thomas __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] persp() and character labels for axis
So your R is outdated. Please upgrade! Uwe Ligges Thomas Steiner wrote: I want to plot dates on the y-axis of a persp() plot. ?persp points you to ?trans3d which is useful to calculate coordinates for calls to 2D functions such as text(). ?trans3d gives this: No documentation for 'trans3d' in specified packages and libraries: you could try 'help.search(trans3d)' Even loading grDevices and graphics libraries did not help. http://finzi.psych.upenn.edu/R/library/grDevices/html/trans3d.html seems to be what I should look for, but there it says noting about text either (I want to label axis! ,something like t=1:10 plot(x=t,y=t^2,type=l) axis(1,labels=LETTERS[1:5]) and more proper...) Thomas __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] importing a VERY LARGE database from Microsoft SQL into R
That's good news that the 3GB switch is the default in R2.2.1, I will upgrade from R2.2.0 today. Its wasn't hard to modify the header, but there was at least one person I used to email the modified file to because he did not have access to the editbin software, so now we don't have to do that anymore. Thanks for mentioning the windows setting, I forgot to mention that you also have to modify the Boot.ini file as described here: http://www.microsoft.com/whdc/system/platform/server/PAE/PAEmem.mspx Does anyone have code that keeps generating random data until the memory is full and then tells you how much memory was successfully used? I could try writing it, but if someone has already done it, thats all the better! Thanks, Roger On 1/25/06, Prof Brian Ripley [EMAIL PROTECTED] wrote: On Tue, 24 Jan 2006, roger bos wrote: This question comes up a number of times what most people will tell you is that even if you get all you data into R you won't be able to do anything with it. By that I mean, you need about 3 or 4 times as much memory as the size of your data object because R will need to create copies of it. I can tell you what I do in case it helps. I also have a SQL Server database and the good thing about having the data in that format is that you probably don't need all of the data all of the time. First of all, a windows machine can handle up to 4GB of RAM, but most software cannot use all of Or even up to 32Gb of RAM. But a single process is limited to 3Gb of user address space, and to 2Gb unless you tell the OS to allow more. it by default. I have 4GB and I also use the windows binary, so that means that whenever I download a new version of R, I have to modify the header file to the it LARGEADDRESSAWARE. Using this trick, I can load up my big matrix into R to the point where task manager shows that R is using about 1.7GB of memory. Despite such large objects, I am able to do regressions and other things with the data, so I am quite happy. If you need more details just let me know. This is the default in the current R 2.2.1, but you also need to set the /3GB flag in Windows. I found the limit was about 1.7Gb without doing that, more like 2.5Gb when I did. The details are in the rw-FAQ Q2.9 (and vary by R version, so please consult that in the version of R you use). On 1/24/06, r user [EMAIL PROTECTED] wrote: I am using R 2.1.1 in a Windows Xp environment. I need to import a large database from Microsoft SQL into R. I am currently using the sqlQuery function/command. This works, but I sometimes run out of memory if my database is too big, or it take quite a long time for the data to import into R. Is there a better way to bring a large SQL database into R? IS there an efficient way to convert the data into R format prior to bringing it into R? (E.g. directly from Microsoft SQL?) __ R-help@stat.math.ethz.ch 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 -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 [[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] reducing learning curves?
In regards to text editors: If you are a Unix user, I'd recommend Emacs (although it has its own large learning curve.) On Windows I use PSpad (www.pspad.com) because it is easy to use and learn and has some of the features you request: syntax highlighting, code completion, code builder, among many other features I find useful. ~Nick -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Behalf Of Michael Sent: Wednesday, January 25, 2006 12:10 AM To: R-help@stat.math.ethz.ch Subject: [R] reducing learning curves? Hi all, I am really new to the R language. I am a long time Matlab and C++ user and I was forced to learn R because I am taking a statistics class. I am seeking to reduce the learning curve to as smooth as possible. Are there any addon/plug-in features that can reduce the learning curve, for example, the following features can be very helpful for new learners: 1. Matlab-like command line auto-completion: Matlab has huge amount of command and nobody is able to remember them off the head. So a nice feature of Matlab command line is that I just need to type the first a few letters and then I press TAB key, there will be a list of possible commands popping up so I just need to select one. This helps a lot in terms of learning for new comers. A more advanced command auto-completion is Visual C++-like, which is implemented in program editor. It helps a lot while doing programming; 2. A good IDE editor with embedded inline debugger: can be as good as VC++, but also can be as simple as Matlab's debugger, which can breakpoint and trace line-by-line... the editor can do syntax correction, syntax check, syntax highlighting, code formatting, etc. Could you please recommend some good addon/plugins that have the above features? Could you please also suggest some tips/tools/tricks that can help me reduce the learning curve? Thank you very much! Michael. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Within-Subjects ANOVA comparisons of individual means
On Tue, 24 Jan 2006, Spencer Graves wrote: 1. Did you try summary(aovRes, ...) rather than summary.aov(aovRes, ...)? From summary.aov, I got the same error message you did, but from summary(aovRes, ...), I got something that looked like what you were expecting. 2. To understand this, look at the code for summary: It consists essentially of 'UseMethod(summary)'. To trace methods dispatch here, note that 'class(aovRes)' is a 2-vector: c(aovlist, listof). When I requested 'methods(summary), I got a long list, the first two were 'summary.aov' and 'summary.aovlist'. Since aovRes was NOT of class 'aov' but instead of class 'aovlist', summary(aovRes, ...) uses 'summary.aovlist'. You got an error message, because summary.aov was expecting an argument of class 'aov', and 'aovRes' didn't have the structure it required. 3. To find an attribute contrasts in aovRes, I requested 'str(aovRes)' and searched for contrasts. I found it in two places: attr(attr(aovRes, error.qr)$qr, contrasts) attr(aovRes, contrasts) 4. I wondered why you said, From my understanding, TukeyHSD is not appropriate in this context. Then I discovered that TukeyHSD is officially a generic function with a method only for objects of class aov. Since aovRes is NOT of class aov, it can't use that function. However, it looks to me like the same algorithm could be used for what you want, e.g., by copying TukeyHSD.aov into a script file, changing the name to TukeyHSD.aovlist and then walking through the code line by line, and changing things as necessary to produce the desired results. If I wanted to do this, I might first check with the author of TukeyHSD (Douglas Bates). He might know some reason why it is inappropriate or some other special concerns of which I'm unaware. Or he might have a not-quite-fully debugged version someplace that could help you. It is a lot more complicated than that. You can probably adjust TukeyHSD to work on a single stratum of an aovlist object, but it is not so obvious that you can retrieve the information needed from the aovlist object so you may need to refit. However, my reason for not pursuing that is more philosophical: these post hoc tests are intended to allow for multiple testing, and that would only make sense if you allow for all the tests done in all the strata. But here there is only stratum, and that suggests the wrong analysis is being done. As I understand this design, it is a randomized block design with subject as block and interval as treatment. The classic analysis is a fixed-effect one (subject + interval), and TukeyHSD can be applied to that. E.g. fm - aov(dv ~ subject + interval) summary(fm) Df Sum Sq Mean Sq F valuePr(F) subject 6 2.4853 0.4142 1.30670.2845 interval 5 13.8174 2.7635 8.7178 3.417e-05 *** Residuals 30 9.5098 0.3170 (should might look familiar) TukeyHSD(fm, interval) Tukey multiple comparisons of means 95% family-wise confidence level Fit: aov(formula = dv ~ subject + interval) $interval diff lwr upr p adj 2-1 -0.7477000 -1.663060 0.16766002 0.1608612 3-1 -1.0704286 -1.985789 -0.15506855 0.0145883 4-1 -1.4761143 -2.391474 -0.56075426 0.0004035 5-1 -1.5005143 -2.415874 -0.58515426 0.0003225 6-1 -1.6827143 -2.598074 -0.76735426 0.601 3-2 -0.3227286 -1.238089 0.59263145 0.8884095 4-2 -0.7284143 -1.643774 0.18694574 0.1814407 5-2 -0.7528143 -1.668174 0.16254574 0.1557224 6-2 -0.9350143 -1.850374 -0.01965426 0.0430664 4-3 -0.4056857 -1.321046 0.50967431 0.7563642 5-3 -0.4300857 -1.345446 0.48527431 0.7095947 6-3 -0.6122857 -1.527646 0.30307431 0.3475861 5-4 -0.0244000 -0.939760 0.89096002 0.994 6-4 -0.2066000 -1.121960 0.70876002 0.9821152 6-5 -0.1822000 -1.097560 0.73316002 0.9898235 However, I find it suspicious that your treatment effects are in fact ordered, and suspect there is something more to the squeezed out here. You appear to be defining difference contrasts without using contr.sdif (in MASS). Once you have a single-stratum fit, a lot more tools become available (se.contrast is one). But for starters fm2 - lm(dv ~ subject+interval, contrasts=list(interval=contr.sdif)) summary(fm2) Coefficients: Estimate Std. Error t value Pr(|t|) ... interval2-1 -0.747700.30095 -2.484 0.0188 * interval3-2 -0.322730.30095 -1.072 0.2921 interval4-3 -0.405690.30095 -1.348 0.1877 interval5-4 -0.024400.30095 -0.081 0.9359 interval6-5 -0.182200.30095 -0.605 0.5495 I hope that is enough food for thought (it is not an invitation for further free consultancy). Steffen Katzner wrote: I am having problems with comparing individual means in a within-subjects ANOVA. From my understanding, TukeyHSD is not appropriate in this context. So I am trying to compute contrasts, as follows: seven subjects participated in each of 6 conditions (intervals).
Re: [R] Is there no definition for global variables in R?
I don't know german, but try z - whatever instead of just z - whatever. ?- On 1/24/06, Christian Hinz [EMAIL PROTECTED] wrote: Hello @all R-Help-User. I need a global variable in R. The variable ought to be known for every functions and subfunctions. It is only to comparison purposes of the numeric algorithms. Is there a possibility? please answer in german if possible. thank you in advance. Christian Hinz -- * Christian Hinz Wickrather Str. 230 41236 Mönchengladbach 02166/125369 __ R-help@stat.math.ethz.ch 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
Re: [R] exporting dates into Microsoft SQL Server
I haven't been able to figure that one out either, but I have a work around. Lets say I have a table named roger_return that has a column named datadate that is a smalldatetime. I can't get sqlSave to save to that table, so I just save it to a new table, say roger_return2. Then I alter the column and copy the data from that table into the original one. It sounds like more work, but once you have the code you just re-use it. Here are the lines that do what I mentioned: go - sqlQuery(x, alter table roger_return2 alter column datadate smalldatetime NOT NULL) go - sqlQuery(x, insert into roger_return select * from roger_return2) HTH, Roger On 1/23/06, r user [EMAIL PROTECTED] wrote: I am running R 2.1.1 in a Windows XP environment. I wish to use the sqlSave command to export a dataframe into Microsoft SQL. My dataframe is called temp and has 2 columns, monthenddate and value. Monthenddate is in 'POSIXct', format. (i.e. 'POSIXct', format: chr 1984-01-31 1984-01-31 1984-01-31 1984-01-31 ...). How can I export this dataframe into SQL and have the format in SQL by one of the standard SQL date formats? I am using the following r code: db - odbcConnect(testserver) sqlSave(db, temp) __ R-help@stat.math.ethz.ch 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
Re: [R] How to use rfm.test ? (Package MarkedPointProcess)
I don't find where it's possible to configure closer bounds in the algorithm. I have obtained some results with rfm.test(coord=steigerwald$coord,steigerwald$diam,MCrepetitions=19,n.hypo=10) but I don't know how to interpret them. Where is the result of the MC test ? Could you give me the syntax of rfm.test for this example (Steigerwald) that you have used in your paper Detecting dependence between Marks and Locations of Marked Point Processes ? regards. Martin Schlather wrote: Dear Florent Bonneu, the optim algorithm with parameter method=L-BFGS-B (used in rfm.test) does not stay always exactly within the given bounds during the search of the optimum. This happens more frequently when the bounds are too wide. Here. rfm.test notices that -9.802347e-17 is less then the given lower bound of value 0. If the messages appear frequently it appeared to me that the algorithms runs into boundary, local minima more frequently. So, if the messages appear only a very few times, ignore them. If they appear several times, use closer bounds. Cheers, Martin Florent Bonneu wrote: I would like to compute the MC test (rfm.test) available in the package MarkedPointProcess (for the data BITOEK for example) in order to test the dependence between the marks and their locations. Why the syntax of rfm.test is false here? I have the message : ** ML WARNING! Forbidden values! -- if there are too many warnings try narrower lower and upper bounds for the variables. -9.802347e-17 1364.372 [ 0 0.01808660 , 9.990131 1364.372 ] ** The program : - # Packages # library(spatstat) library(RandomFields) library(adapt) library(MarkedPointProcess) # Program # data(BITOEK) win - ripras(steigerwald$coord) PointProcess - ppp(x=steigerwald$coord[,1],y=steigerwald$coord[,2],window=win,marks=steigerwald$diam) plot(PointProcess) rfm.test(coord=steigerwald$coord,steigerwald$diam,MCrepetitions=19) - regards. Florent Bonneu Laboratoire de Statistique et Probabilités bureau 148 bât. 1R2 Université Toulouse 3 118 route de Narbonne - 31062 Toulouse cedex 9 [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] Problems running R immediately after installing
Hello, I am an employee at The George Washington University, and we have recently had trouble getting R to launch without error directly after installing. There are no errors displayed during the install, but when the program is run, a warning is displayed for ever package that is referenced in defaultPackages which says the package cannot be found. I have checked the .libPaths() value, and it is pointing to the correct location (C:\...\R-2.2.0\library). In addition, I have attempted to re-install the software with the same result consistently. I have attempted searching online for this problem, and all posts online point to the .libPaths() value being incorrect, but I have checked this value and reset it a number of times, but the problem continues to persist. Is there any known solution to this problem or another pool of answers I should look at? Thanks. Brent Roland Desktop Management Assistant The George Washington University P: (202) 994-4142 C: (202) 746-5957 [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Savitsky Golay smoothing
Hi, I am trying to process some spectral data using R, and I would like to use Savitsky Golay smoothing. Is this already implemented in R or one of the optional packages, or do I have to implement is myself? Thanks in advance, Dirk -- Dirk De Becker Work: Kasteelpark Arenberg 30 3001 Heverlee phone: ++32(0)16/32.14.44 fax: ++32(0)16/32.85.90 Home: Waversebaan 90 3001 Heverlee phone: ++32(0)16/23.36.65 [EMAIL PROTECTED] mobile phone: ++32(0)498/51.19.86 Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm __ R-help@stat.math.ethz.ch 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] Savitsky Golay smoothing
Try: RSiteSearch(Savitzky-Golay) On 1/25/06, Dirk De Becker [EMAIL PROTECTED] wrote: Hi, I am trying to process some spectral data using R, and I would like to use Savitsky Golay smoothing. Is this already implemented in R or one of the optional packages, or do I have to implement is myself? Thanks in advance, Dirk -- Dirk De Becker Work: Kasteelpark Arenberg 30 3001 Heverlee phone: ++32(0)16/32.14.44 fax: ++32(0)16/32.85.90 Home: Waversebaan 90 3001 Heverlee phone: ++32(0)16/23.36.65 [EMAIL PROTECTED] mobile phone: ++32(0)498/51.19.86 Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm __ R-help@stat.math.ethz.ch 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] Unequal sample sizes when calculating power
Power calculations two sample test for proportions is very useful. Is there
a way however, to get away from the two samples being of the same size. What
would happen if one had n=15 in the one sample and n=45 in the other sample.
Farrel Buchinsky, MD
Pediatric Otolaryngologist
Allegheny General Hospital
Pittsburgh, PA
**
This email and any files transmitted with it are confidentia...{{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] is.integer() function
Hi R users I have a simple question to ask. Why am I getting FALSE on this is.integer(10) [1] FALSE 10 is a integer number. Thanks in advance. M __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] [R-pkgs] amap
Dear R users, Version 0.7-1 of the amap package has been uploaded to CRAN. Amap package includes standard hierarchical clustering and k-means. We optimize implementation (with a parallelized hierarchical clustering) and allow the possibility of using different distances like Eulidean or Spearman (rank-based metric). We implement a principal component analysis (with robusts methods). News of version 0.7-1: a new method of clustering optimal partiton A major advantage of the method is that the number K of classes is computed from the data, and ONLY from the data. -- Antoine Lucas Centre de génétique Moléculaire, CNRS 91198 Gif sur Yvette Cedex Tel: (33)1 69 82 38 89 Fax: (33)1 69 82 38 77 ___ R-packages mailing list [EMAIL PROTECTED] https://stat.ethz.ch/mailman/listinfo/r-packages __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] is.integer() function
Becaues is.integer shows the internal representation, which is not an integer but a double (real number). Some functions create integer vectors, for example the : notation: is.integer(1:10) [1] TRUE You can create an integer vector with the integer() function: integer(10) [1] 0 0 0 0 0 0 0 0 0 0 is.integer(integer(10)) [1] TRUE Gabor On Tue, Jan 24, 2006 at 09:29:51PM -0600, Taka Matzmoto wrote: Hi R users I have a simple question to ask. Why am I getting FALSE on this is.integer(10) [1] FALSE 10 is a integer number. Thanks in advance. M __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Csardi Gabor [EMAIL PROTECTED]MTA RMKI, ELTE TTK __ R-help@stat.math.ethz.ch 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] importing a VERY LARGE database from Microsoft SQL into R
On Wed, 25 Jan 2006, roger bos wrote: Does anyone have code that keeps generating random data until the memory is full and then tells you how much memory was successfully used? I could try writing it, but if someone has already done it, thats all the better! In recent versions of R gc() reports the maximum memory used, so all you need to do is write the code to use lots of memory. -thomas __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] is.integer() function
S Poetry page 125. Patrick Burns [EMAIL PROTECTED] +44 (0)20 8525 0696 http://www.burns-stat.com (home of S Poetry and A Guide for the Unwilling S User) Taka Matzmoto wrote: Hi R users I have a simple question to ask. Why am I getting FALSE on this is.integer(10) [1] FALSE 10 is a integer number. Thanks in advance. M __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] is.integer() function
Hi Numbers are normaly stored as double so you has to declare it as an integer. str(10) num 10 str(as.integer(10)) int 10 HTH Petr On 24 Jan 2006 at 21:29, Taka Matzmoto wrote: From: Taka Matzmoto [EMAIL PROTECTED] To: r-help@stat.math.ethz.ch Date sent: Tue, 24 Jan 2006 21:29:51 -0600 Subject:[R] is.integer() function Hi R users I have a simple question to ask. Why am I getting FALSE on this is.integer(10) [1] FALSE 10 is a integer number. Thanks in advance. M __ R-help@stat.math.ethz.ch 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
Re: [R] is.integer() function
On Tue, 24 Jan 2006, Taka Matzmoto wrote: Hi R users I have a simple question to ask. Why am I getting FALSE on this is.integer(10) [1] FALSE 10 is a integer number. is.integer() asks how a number is stored, not whether the number happens to be an integer. Explicit numbers typed into R are stored as double precision even if they happen to be whole numbers. This is historical but is also useful: the double precision type can hold integers up to 2^52 exactly and the integer type can only hold integers up to 2^31. -thomas __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] is.integer() function
Gabor Csardi wrote: Becaues is.integer shows the internal representation, which is not an integer but a double (real number). Some functions create integer vectors, Some functions that you might think create integer vectors and even seem to say they create integer vectors dont create integer vectors: 'ceiling' takes a single numeric argument 'x' and returns a numeric vector containing the smallest integers not less than the corresponding elements of 'x'. ceiling(0.5) [1] 1 is.integer(ceiling(0.5)) [1] FALSE is.integer(1:3) [1] TRUE is.integer(ceiling(1:3)) [1] FALSE This could possibly be a documentation problem, since ?ceiling is using 'integer' in the sense of 'whole number', whereas ?is.integer is concerned with internal representation (aka 'storage mode') This seems to be an endless source of confusion to anyone who didn't start their programming days in Fortran, C, or assembly language (or other strongly-typed language, I guess). Barry __ R-help@stat.math.ethz.ch 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] combining variables with PCA
Christian, One of the arguments to prcomp() is retx, with a default value of TRUE. As explained in the help file, if retx is TRUE the prcomp object reurned by the function contains the projection of the original data along the principal components (which many of us call scores). Thus z - prcomp(cbind(runif(20),runif(20))) z$x PC1 PC2 [1,] -0.176564769 -0.175010202 [2,] 0.286746995 -0.464200150 [3,] -0.010385993 0.415020704 [4,] 0.002108131 -0.150929991 [5,] -0.030974478 0.005799164 [6,] 0.516400555 0.051731728 [7,] -0.273545829 0.082611939 [8,] 0.028956640 -0.095674198 [9,] -0.343703128 0.148523268 [10,] -0.267409591 0.236258325 [11,] -0.362410088 0.278893388 [12,] 0.198170683 -0.282954098 [13,] -0.195901187 0.030789448 [14,] 0.544997090 -0.089702393 [15,] -0.361799979 -0.055526763 [16,] 0.117021002 0.362033198 [17,] -0.385818070 -0.459613063 [18,] -0.182413187 -0.218456938 [19,] 0.386595491 0.337704022 [20,] 0.509929714 0.042702611 You can use the columns of the x component of your prcomp object as the independent variable in your regression analysis or whatever. David W. Roberts office 406-994-4548 Professor and Head FAX 406-994-3190 Department of Ecology email [EMAIL PROTECTED] Montana State University Bozeman, MT 59717-3460 Christian Jones wrote: hello R_team having perfomed a PCA on my fitted model with the function: data- na.omit(dataset) data.pca-prcomp(data,scale =TRUE), I´ve decided to aggregate two variables that are highly correlated. My first question is: How can I combine the two variables into one new predictor? and secondly: How can I predict with the newly created variable in a new dataset? Guess I need the predict and new data command but I´m having problems with the syntax and the help function was not sufficient on this issue. many thanks in advance Christian [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] reducing learning curves?
Hello, If you work under Windows, you can find a lot of useful tools in SciViews-R (http://www.sciviews.org/SciViews-R) and Tinn-R (http://www.sciviews.org/Tinn-R). For instance, you have: - syntax coloring, - code completion, - calltips (tips displaying the syntax of a function as you type it), - object explorer with lots of useful shortcuts in the object's context menu, - electronic reference cards, - viewing and reproting features, - etc... For a nice, graphical, debugger, look at debug package (you have to install it from CRAN and load it using: library(debug) Then, try: ?mtrace Best, Philippe Grosjean Michael wrote: Hi all, I am really new to the R language. I am a long time Matlab and C++ user and I was forced to learn R because I am taking a statistics class. I am seeking to reduce the learning curve to as smooth as possible. Are there any addon/plug-in features that can reduce the learning curve, for example, the following features can be very helpful for new learners: 1. Matlab-like command line auto-completion: Matlab has huge amount of command and nobody is able to remember them off the head. So a nice feature of Matlab command line is that I just need to type the first a few letters and then I press TAB key, there will be a list of possible commands popping up so I just need to select one. This helps a lot in terms of learning for new comers. A more advanced command auto-completion is Visual C++-like, which is implemented in program editor. It helps a lot while doing programming; 2. A good IDE editor with embedded inline debugger: can be as good as VC++, but also can be as simple as Matlab's debugger, which can breakpoint and trace line-by-line... the editor can do syntax correction, syntax check, syntax highlighting, code formatting, etc. Could you please recommend some good addon/plugins that have the above features? Could you please also suggest some tips/tools/tricks that can help me reduce the learning curve? Thank you very much! Michael. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Unequal sample sizes when calculating power
Farrel Buchinsky [EMAIL PROTECTED] writes: Power calculations two sample test for proportions is very useful. Is there a way however, to get away from the two samples being of the same size. What would happen if one had n=15 in the one sample and n=45 in the other sample. Take a look at https://www.stat.math.ethz.ch/pipermail/r-devel/2003-September/027458.html That's for t-tests, but similar modifications might be used for power.prop.test -- 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
Re: [R] read.table problem
Hello,
Well... the error message is explicit enough: number of items read is
not a multiple of the number of columns means that you do not have the
right number of items around line 11975 (not the same number as in the
11974 previous lines)! This is an error in you file.
Best,
Philippe Grosjean
Andrej Kastrin wrote:
Dear R useRs,
I have big (23000 rows), vertical bar delimited file:
e.g.
A1|Text a,Text b, Text c|345
A2|Text bla|456
...
..
.
Try using
A - read.table('filename.txt', header=FALSE,sep='\|')
process stop at line 11975 with warning message:
number of items read is not a multiple of the number of columns
I have no problems with processing similar file, which is only 1
rows long?
Any suggestion what's the problem here. Thank's in advance.
Cheers, Andrej
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R vs. Excel (R-squared)
I was quite interested in this thread (discussion),
once that I am chemistry student and I work with Mixtures Designs that are
models without intercept.
I thought quite attention the follow afirmation:
' Thus SST, the corrected total
sum of squares, should be used when you have a model with an intercept
term but the uncorrected total sum of squares should be used when you
do not have an intercept term. ' (Douglas Bates)
I have as reference a book called:
Experiments with Mixtures: Designs, Models, and the Analysis of Mixture
Data
second edition
John A. Cornell
(Professor of Statistics in University Of Florida)
In this book, pg 42: item 2.7 - THE ANALYSIS OF VARIANCE TABLE,
I have the model below:
y(x) = 11.7x1 + 9.4x2 + 16.4x3 + 19.0x1x2 + 11.4x1x3 - 9.6x2x3
with the follow ANOVA Table:
source of variationD.F.SSMS
Regressionp-1SSR=\sum( y_{pred} - y_{mean} )^2ssR/(p-1)
ResidualN-pSSE=\sum( y_{exp} - y_{pred} )^2ssE/(N-p)
TotalN-1SSE=\sum( y_{exp} - y_{mean} )^2
pred = predicted
exp = experimental
and in many others books.
I always see the ANOVA Table of Mixtures systems with SST, the
corrected total
sum of squares ( N-1 degrees freedom ).
I would like to ask:
1) What is approach ( point view ) more adequate ?
2) Could someone indicate some reference about this subject ?
Thanks a lot.
Regards
Cleber N. Borges
Dados
x1 x2 x3 y
10011
10012.4
0.50.5015
0.50.5014.8
0.50.5016.1
0108.8
01010
00.50.510
00.50.59.7
00.50.511.8
00116.8
00116
0.500.517.7
0.500.516.4
0.500.516.6
## Model
d.lm - lm( y ~ -1 + x1*x2*x3 - x1:x2:x3, data = Dados )
### Anova like in the book
d.aov - aov( y ~ x1*x2*x3 - x1:x2:x3, data = Dados )
SSR (fitted Model) = 128.296
Douglas Bates wrote:
On 1/24/06, Lance Westerhoff [EMAIL PROTECTED] wrote:
Hi-
On Jan 24, 2006, at 12:08 PM, Peter Dalgaard wrote:
Lance Westerhoff [EMAIL PROTECTED] writes:
Hello All-
I found an inconsistency between the R-squared reported in Excel vs.
that in R, and I am wondering which (if any) may be correct and if
this is a known issue. While it certainly wouldn't surprise me if
Excel is just flat out wrong, I just want to make sure since the R-
squared reported in R seems surprisingly high. Please let me know if
this is the wrong list. Thanks!
Excel is flat out wrong. As the name implies, R-squared values cannot
be less than zero (adjusted R-squared can, but I wouldn't think
that is what Excel does).
I had thought the same thing, but then I came across the following
site which states: Note that it is possible to get a negative R-
square for equations that do not contain a constant term. If R-square
is defined as the proportion of variance explained by the fit, and if
the fit is actually worse than just fitting a horizontal line, then R-
square is negative. In this case, R-square cannot be interpreted as
the square of a correlation. Since
R^2 = 1 - (SSE/SST)
I guess you can have SSE SST which would result in a R^2 of less
then 1.0. However, it still seems very strange which made me wonder
what is going on in Excel needless to say!
http://www.mathworks.com/access/helpdesk/help/toolbox/curvefit/
ch_fitt9.html
This seems to be a case of using the wrong formula. R^2 should
measure the amount of variation for which the given model accounts
relative to the amount of variation for which the *appropriate* null
model does not account. If you have a constant or intercept term in a
linear model then the null model for comparison is one with the
intercept only. If you have a linear model without an intercept term
then the appropriate null model for comparison is the model that
predicts all the responses as zero. Thus SST, the corrected total
sum of squares, should be used when you have a model with an intercept
term but the uncorrected total sum of squares should be used when you
do not have an intercept term.
It is disappointing to see the MathWorks propagating such an
elementary misconception.
__
R-help@stat.math.ethz.ch 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] (no subject)
hello, i'm doing a master thesis in biology, and i'm studying genetic differenciation between population. So i used genetic distances, called FST, defining as a proportion of variance between two population. Then i'd like to use this pairwise FST as a dependent variable in glm. But i'm not sure of the error's distribution to use... Gaussian or binomial? do you have an idea? thanks a lot! guillaume emaresi __ R-help@stat.math.ethz.ch 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] Unequal sample sizes when calculating power
On Wed, 2006-01-25 at 09:50 -0500, Farrel Buchinsky wrote: Power calculations two sample test for proportions is very useful. Is there a way however, to get away from the two samples being of the same size. What would happen if one had n=15 in the one sample and n=45 in the other sample. See ?bpower and ?bpower.sim (along with the other functions on the same help page) in Frank Harrell's 'Hmisc' package on CRAN. HTH, Marc Schwartz P.S. Go Steelers! :-) __ R-help@stat.math.ethz.ch 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] reducing learning curves?
Keeping this reference card handy might reduce it somewhat: http://www.rpad.org/Rpad/Rpad-refcard.pdf On 1/25/06, Michael [EMAIL PROTECTED] wrote: Hi all, I am really new to the R language. I am a long time Matlab and C++ user and I was forced to learn R because I am taking a statistics class. I am seeking to reduce the learning curve to as smooth as possible. Are there any addon/plug-in features that can reduce the learning curve, for example, the following features can be very helpful for new learners: 1. Matlab-like command line auto-completion: Matlab has huge amount of command and nobody is able to remember them off the head. So a nice feature of Matlab command line is that I just need to type the first a few letters and then I press TAB key, there will be a list of possible commands popping up so I just need to select one. This helps a lot in terms of learning for new comers. A more advanced command auto-completion is Visual C++-like, which is implemented in program editor. It helps a lot while doing programming; 2. A good IDE editor with embedded inline debugger: can be as good as VC++, but also can be as simple as Matlab's debugger, which can breakpoint and trace line-by-line... the editor can do syntax correction, syntax check, syntax highlighting, code formatting, etc. Could you please recommend some good addon/plugins that have the above features? Could you please also suggest some tips/tools/tricks that can help me reduce the learning curve? Thank you very much! Michael. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R-help Digest, Vol 35, Issue 24
[Gabor Grothendieck] [...] this list is inhabited by some rather rude participants but everyone puts up with them in the hope that they do have some useful remarks. I've been witnessing this list for about one year, and also read *lots* of archived messages. While it is true that a few members do not use white gloves, are rather fond on concise replies, and do express strong opinions at times, they never went overboard insulting people and always kept a reasonable measure, at least so far that I could see (yet who knows, outliers might happen! :-). (*) Our whole society is a bit shy and shivers easily when opinions are expressed nowadays, I often observed than people quickly get insecure, feel attacked, and overreact (by running away or starting a fight). there is even a group of thought that feels it is a justifiable way to keep the list volume under control. This may work because of the starred paragraph above, that is, for wrong reasons. Best is, and this often occurs on the R list, when everything (facts, opinions) is being shared efficiently, without useless arguing. Then, threads quickly fade out. -- François Pinard http://pinard.progiciels-bpi.ca __ R-help@stat.math.ethz.ch 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] About lmer output
Dear R users: I am using lmer fo fit binomial data with a probit link function: fer_lmer_PQL-lmer(fer ~ gae + ctipo + (1|perm) -1, +family = binomial(link=probit), +method = 'PQL', +data = FERTILIDAD, +msVerbose= True) The output look like this: fer_lmer_PQL Generalized linear mixed model fit using PQL Formula: fer ~ gae + ctipo + (1 | perm) - 1 Data: FERTILIDAD Family: binomial(probit link) AIC BIClogLik deviance 2728.086 2918.104 -1332.043 2664.086 Random effects: GroupsNameVarianceStd.Dev. perm (Intercept) 0.28256 0.53156 # of obs: 2802, groups: perm, 529 Estimated scale (compare to 1) 0.8958656 My question is about the meaning of Estimated scale (compare to 1) 0.8958656 I think that the scale would be 0.28256+1.0, Isn´t it? Thanks in advance, Juan Pablo Sánchez __ R-help@stat.math.ethz.ch 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] Log-Likelihood 3d-plot and contourplot / optim() starting values
Hello,
i have coded the following loglikelihood-function
# Log-Likelihood-Funktion
loglik_jm-function(N,phi,t) {
n-length(t)
i-seq(along=t)
s1-sum(log(N-(i-1)))
s2-phi*sum((N-(i-1))*t[i])
n*log(phi)+s1-s2
}
# the data
t-c(7,11,8,10,15,22,20,25,28,35)
# now i want to do a 3d-plot and a contourplot in order to see at which
values of N and phi the loglikelihood function becomes zero.
# i do this in order to get an idea where the starting values for a
optimization for the mle of N could be
# 3dplot and contourplot
phi-seq(0,1,length=50)
N-seq(101,110,length=50)
z-outer(N,phi,loglik_jm(N,phi,t))
persp(phi,N,z, theta = 30, phi = 30, expand = 0.5, col = lightblue)
contourplot(z~N*phi)
# but i get some error messages, i don't know why?
# if you are interested, the mle function for N is
ll2 - function(N,t) {
i-seq(along=t)
n-length(t)
s1-sum(1/(N-(i-1)))
s2-n/(N-(1/sum(t[i]))*(sum((i-1)*t[i])))
(s1-s2)
}
# you get this function as usual, if you set the loglikelihood equal
zero and differentiate for N
# i take the squares of ll2 in order to get the minimum easier
ll3-function(N,t) {
ll2(N,t)^2
}
# then i do an iteration to get the mle of N
xmin-optim(10,ll3,method=BFGS,control=list(reltol=.Machine$double.eps),t=t)
# the problem is that this method only works well, if the starting
values for optim() are very close to the real value!
# so i got the idea with contourplot() and persp() to see where good
starting values could be.
i would be very thankful if anyone could give me some advice!
regards
Andreas
__
R-help@stat.math.ethz.ch 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 vs. Excel (R-squared)
Cleber N. Borges wrote:
I was quite interested in this thread (discussion),
once that I am chemistry student and I work with Mixtures Designs that are
models without intercept.
I thought quite attention the follow afirmation:
' Thus SST, the corrected total
sum of squares, should be used when you have a model with an intercept
term but the uncorrected total sum of squares should be used when you
do not have an intercept term. ' (Douglas Bates)
I have as reference a book called:
Experiments with Mixtures: Designs, Models, and the Analysis of Mixture
Data
second edition
John A. Cornell
(Professor of Statistics in University Of Florida)
In this book, pg 42: item 2.7 - THE ANALYSIS OF VARIANCE TABLE,
I have the model below:
y(x) = 11.7x1 + 9.4x2 + 16.4x3 + 19.0x1x2 + 11.4x1x3 - 9.6x2x3
with the follow ANOVA Table:
source of variationD.F.SSMS
Regressionp-1SSR=\sum( y_{pred} - y_{mean} )^2ssR/(p-1)
ResidualN-pSSE=\sum( y_{exp} - y_{pred} )^2ssE/(N-p)
TotalN-1SSE=\sum( y_{exp} - y_{mean} )^2
pred = predicted
exp = experimental
and in many others books.
I always see the ANOVA Table of Mixtures systems with SST, the
corrected total
sum of squares ( N-1 degrees freedom ).
I would like to ask:
1) What is approach ( point view ) more adequate ?
With a mixture model, although you do not have a intercept term
directly in the model, it is there, occulted,as the sum of the design
variables representing the mixture is 1! So it is correct to use the
corrected sum of squares.
Kjetil
2) Could someone indicate some reference about this subject ?
Thanks a lot.
Regards
Cleber N. Borges
Dados
x1 x2 x3 y
10011
10012.4
0.50.5015
0.50.5014.8
0.50.5016.1
0108.8
01010
00.50.510
00.50.59.7
00.50.511.8
00116.8
00116
0.500.517.7
0.500.516.4
0.500.516.6
## Model
d.lm - lm( y ~ -1 + x1*x2*x3 - x1:x2:x3, data = Dados )
### Anova like in the book
d.aov - aov( y ~ x1*x2*x3 - x1:x2:x3, data = Dados )
SSR (fitted Model) = 128.296
Douglas Bates wrote:
On 1/24/06, Lance Westerhoff [EMAIL PROTECTED] wrote:
Hi-
On Jan 24, 2006, at 12:08 PM, Peter Dalgaard wrote:
Lance Westerhoff [EMAIL PROTECTED] writes:
Hello All-
I found an inconsistency between the R-squared reported in Excel vs.
that in R, and I am wondering which (if any) may be correct and if
this is a known issue. While it certainly wouldn't surprise me if
Excel is just flat out wrong, I just want to make sure since the R-
squared reported in R seems surprisingly high. Please let me know if
this is the wrong list. Thanks!
Excel is flat out wrong. As the name implies, R-squared values cannot
be less than zero (adjusted R-squared can, but I wouldn't think
that is what Excel does).
I had thought the same thing, but then I came across the following
site which states: Note that it is possible to get a negative R-
square for equations that do not contain a constant term. If R-square
is defined as the proportion of variance explained by the fit, and if
the fit is actually worse than just fitting a horizontal line, then R-
square is negative. In this case, R-square cannot be interpreted as
the square of a correlation. Since
R^2 = 1 - (SSE/SST)
I guess you can have SSE SST which would result in a R^2 of less
then 1.0. However, it still seems very strange which made me wonder
what is going on in Excel needless to say!
http://www.mathworks.com/access/helpdesk/help/toolbox/curvefit/
ch_fitt9.html
This seems to be a case of using the wrong formula. R^2 should
measure the amount of variation for which the given model accounts
relative to the amount of variation for which the *appropriate* null
model does not account. If you have a constant or intercept term in a
linear model then the null model for comparison is one with the
intercept only. If you have a linear model without an intercept term
then the appropriate null model for comparison is the model that
predicts all the responses as zero. Thus SST, the corrected total
sum of squares, should be used when you have a model with an intercept
term but the uncorrected total sum of squares should be used when you
do not have an intercept term.
It is disappointing to see the MathWorks propagating such an
elementary misconception.
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
.
Re: [R] is.integer() function
On 1/25/2006 10:57 AM, Barry Rowlingson wrote: Gabor Csardi wrote: Becaues is.integer shows the internal representation, which is not an integer but a double (real number). Some functions create integer vectors, Some functions that you might think create integer vectors and even seem to say they create integer vectors dont create integer vectors: 'ceiling' takes a single numeric argument 'x' and returns a numeric vector containing the smallest integers not less than the corresponding elements of 'x'. ceiling(0.5) [1] 1 is.integer(ceiling(0.5)) [1] FALSE is.integer(1:3) [1] TRUE is.integer(ceiling(1:3)) [1] FALSE This could possibly be a documentation problem, since ?ceiling is using 'integer' in the sense of 'whole number', whereas ?is.integer is concerned with internal representation (aka 'storage mode') Here numeric vector is being used in the R-specific technical sense as a vector of double precision values, so the documentor was trying hard to be precise. The problem is that English also admits the interpretation in a non-technical sense as a vector of numbers. I believe your country is to blame for the language. :-) Duncan Murdoch This seems to be an endless source of confusion to anyone who didn't start their programming days in Fortran, C, or assembly language (or other strongly-typed language, I guess). Barry __ R-help@stat.math.ethz.ch 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] read.table problem
By the way, you might find this sed one-liner useful:
sed -n '11981q;11970,11980p' filename.txt
It will print the offending line and its neighbors. If you're on
Windows you need to install Windows Services For Unix or Cygwin.
-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Andrej Kastrin
Sent: Wednesday, January 25, 2006 3:08 AM
To: r-help
Subject: [R] read.table problem
Dear R useRs,
I have big (23000 rows), vertical bar delimited file:
e.g.
A1|Text a,Text b, Text c|345
A2|Text bla|456
...
..
.
Try using
A - read.table('filename.txt', header=FALSE,sep='\|')
process stop at line 11975 with warning message:
number of items read is not a multiple of the number of columns
I have no problems with processing similar file, which is only 1
rows long?
Any suggestion what's the problem here. Thank's in advance.
Cheers, Andrej
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] is.integer() function
Duncan Murdoch wrote:
Here numeric vector is being used in the R-specific technical sense as
a vector of double precision values, so the documentor was trying hard
to be precise. The problem is that English also admits the
interpretation in a non-technical sense as a vector of numbers. I
believe your country is to blame for the language. :-)
But can numeric vector in the R-specific technical sense also mean
a vector of integer (representation) values? It passes is.numeric() and
is.vector():
x=1:3
is.numeric(x)
[1] TRUE
is.vector(x)
[1] TRUE
is.integer(x)
[1] TRUE
Unless by the R-specific technical sense of 'numeric vector' you
dont mean something for which is.numeric() and is.vector() are both
true. Which is perverse. But then large chunks of R are. Anyway, is this
right:
is.RSpecificTechnicalSenseNumericVector=function(v){is.numeric(v)
is.double(v)}
Barry
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] is.integer() function
My understanding is that a variable of _mode_ numeric can be
of _type_ double or _type_ integer.Just knowing that a variable
is numeric does not tell you which of these two subcategories it
falls in.
x - 1:3
mode(x) # numeric
typeof(x) # integer
y - 1.2
mode(y) # numeric
typeof(y) # double
On 1/25/06, Barry Rowlingson [EMAIL PROTECTED] wrote:
Duncan Murdoch wrote:
Here numeric vector is being used in the R-specific technical sense as
a vector of double precision values, so the documentor was trying hard
to be precise. The problem is that English also admits the
interpretation in a non-technical sense as a vector of numbers. I
believe your country is to blame for the language. :-)
But can numeric vector in the R-specific technical sense also mean
a vector of integer (representation) values? It passes is.numeric() and
is.vector():
x=1:3
is.numeric(x)
[1] TRUE
is.vector(x)
[1] TRUE
is.integer(x)
[1] TRUE
Unless by the R-specific technical sense of 'numeric vector' you
dont mean something for which is.numeric() and is.vector() are both
true. Which is perverse. But then large chunks of R are. Anyway, is this
right:
is.RSpecificTechnicalSenseNumericVector=function(v){is.numeric(v)
is.double(v)}
Barry
__
R-help@stat.math.ethz.ch 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] Question about fitting power
Hi R users I'm trying to fit a model y=ax^b. I know if I made ln(y)=ln(a)+bln(x) this is a linear regression. But I obtain differente results with nls() and lm() My commands are: nls(CV ~a*Est^b, data=limiares, start =list(a=100,b=0), trace = TRUE) for nonlinear regression and : lm(ln_CV~ln_Est, data=limiares) for linear regression Nonlinear regression model: a=738.2238151 andb=-0.3951013 Linear regression: Coefficients: EstimateStd. Errort value Pr(|t|) (Intercept) 7.8570224 0.0103680 757.8 2e-16 *** ln_Est -0.5279412 0.0008658 -609.8 2e-16 *** I think it should be a=exp((Intercept) ) = exp(7.8570224) = 2583.815 and b=ln_Est Probably I'm wrong, but why?? Thanks in advance. Ana Quiterio [[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] Question about fitting power
The two methods are fitting different models. With lm(), the model is y = a * x^b * error or, equivalently, ln(y) = ln(a) + b * ln(x) + ln(error) With nls(), the model is y = a * x^b + error Thus you will get two different estimates. Andy From: Ana Quitério Hi R users I'm trying to fit a model y=ax^b. I know if I made ln(y)=ln(a)+bln(x) this is a linear regression. But I obtain differente results with nls() and lm() My commands are: nls(CV ~a*Est^b, data=limiares, start =list(a=100,b=0), trace = TRUE) for nonlinear regression and : lm(ln_CV~ln_Est, data=limiares) for linear regression Nonlinear regression model: a=738.2238151 andb=-0.3951013 Linear regression: Coefficients: EstimateStd. Errort value Pr(|t|) (Intercept) 7.8570224 0.0103680 757.8 2e-16 *** ln_Est -0.5279412 0.0008658 -609.8 2e-16 *** I think it should be a=exp((Intercept) ) = exp(7.8570224) = 2583.815 and b=ln_Est Probably I'm wrong, but why?? Thanks in advance. Ana Quiterio [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] paste - eliminate spaces?
I am trying to combine the value of a variable and text. e.g. I want test1, with no spaces. I try: h=1 paste(test,1) But get: [1] test 1 (i.e. there is a space between test and 1) Is there a way to eliminate the space? __ R-help@stat.math.ethz.ch 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] paste - eliminate spaces?
I found the answer:
add sep= to the paste command
paste('test',1,sep=)
--- r user [EMAIL PROTECTED] wrote:
I am trying to combine the value of a variable and
text.
e.g.
I want test1, with no spaces.
I try:
h=1
paste(test,1)
But get:
[1] test 1
(i.e. there is a space between test and 1)
Is there a way to eliminate the space?
__
Do You Yahoo!?
protection around
http://mail.yahoo.com
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] paste - eliminate spaces?
Try
paste('test',1,sep='')
-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of r user
Sent: Wednesday, January 25, 2006 1:58 PM
To: rhelp
Subject: [R] paste - eliminate spaces?
I am trying to combine the value of a variable and text.
e.g.
I want test1, with no spaces.
I try:
h=1
paste('test',1)
But get:
[1] test 1
(i.e. there is a space between test' and 1)
Is there a way to eliminate the space?
__
R-help@stat.math.ethz.ch 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] paste - eliminate spaces?
I don't think the help page for paste is all that hard to read or
understand, is it? Please read about the `sep' option (and note its
default).
Andy
From: r user
I am trying to combine the value of a variable and
text.
e.g.
I want test1, with no spaces.
I try:
h=1
paste('test',1)
But get:
[1] test 1
(i.e. there is a space between test' and 1)
Is there a way to eliminate the space?
__
R-help@stat.math.ethz.ch 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] paste - eliminate spaces?
Just read more carefully the online help for paste (?paste). You have: sep: a character string to separate the terms. and the default value for sep is (a space). So, just use: paste(test, 1, sep = ) Best, Philippe Grosjean r user wrote: I am trying to combine the value of a variable and text. e.g. I want “test1”, with no spaces. I try: h=1 paste(‘test’,1) But get: [1] test 1 (i.e. there is a space between “test’“ and “1”) Is there a way to eliminate the space? __ R-help@stat.math.ethz.ch 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] About lmer output
On 1/25/06, Juan Pablo Sánchez [EMAIL PROTECTED] wrote: Dear R users: I am using lmer fo fit binomial data with a probit link function: fer_lmer_PQL-lmer(fer ~ gae + ctipo + (1|perm) -1, +family = binomial(link=probit), +method = 'PQL', +data = FERTILIDAD, +msVerbose= True) The output look like this: fer_lmer_PQL Generalized linear mixed model fit using PQL Formula: fer ~ gae + ctipo + (1 | perm) - 1 Data: FERTILIDAD Family: binomial(probit link) AIC BIClogLik deviance 2728.086 2918.104 -1332.043 2664.086 Random effects: GroupsNameVarianceStd.Dev. perm (Intercept) 0.28256 0.53156 # of obs: 2802, groups: perm, 529 Estimated scale (compare to 1) 0.8958656 My question is about the meaning of Estimated scale (compare to 1) 0.8958656 I think that the scale would be 0.28256+1.0, Isn´t it? The estimated scale is what would be the estimate of the scale parameter in the GLM family if there was a scale parameter. For the binomial and Poisson families there is no scale parameter but the Iteratively Reweighted Least Squares (IRLS) algorithm still produces an estimate of one. If the data are neither overdispersed nor underdispersed then that estimate should be close to 1. It can provide a diagnostic for the model. A value that is substantially different from 1 indicates model failure or over-modeling the data. I would say that the value of 0.896 is close to an indication of underdispersion. Frequently this is caused be including random effects associated with groups that have few observations in each group. __ R-help@stat.math.ethz.ch 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] Question about fitting power
-Original Message-
From: [EMAIL PROTECTED] [SMTP:[EMAIL PROTECTED] On Behalf Of Ana Quitério
Sent: Wednesday, January 25, 2006 3:33 PM
To: r-help@stat.math.ethz.ch
Subject: [R] Question about fitting power
From: Ana Quitério
Hi R users
I'm trying to fit a model y=ax^b.
I know if I made ln(y)=ln(a)+bln(x) this is a linear regression.
But I obtain differente results with nls() and lm()
My commands are:
nls(CV ~a*Est^b, data=limiares, start =list(a=100,b=0), trace = TRUE)
for nonlinear regression, and :
lm(ln_CV~ln_Est, data=limiares)
for linear regression
Nonlinear regression model: a=738.2238151 andb=-0.3951013
Linear regression: Coefficients:
EstimateStd. Errort value Pr(|t|)
(Intercept) 7.8570224 0.0103680 757.8 2e-16 ***
ln_Est -0.5279412 0.0008658 -609.8 2e-16 ***
I think it should be a=exp((Intercept) ) =
exp(7.8570224) = 2583.815 and b=ln_Est
Probably I'm wrong, but why??
Thanks in advance.
Ana Quiterio
In addition to what Andy said, maybe your nls() model is not natural because a
power
statistical model should produce multiplicative outcomes. An alternative is to
try the
nonlinear power model with multiplicative deviates assuming lognormality (which
is
what the linear model does but without the need for back-transformation), with
nlm(),
for instance:
fn-function(p){
+ CV_mod=p[1]*Est^p[2];
+ squdiff=(log(limiares$CV)-log(CV_mod))^2;
+ lik=(length(limiares$CV)/2)*log(sum(squdiff)/length(limiares$CV))
+ }
CV.lik-nlm(fn,p=c(100,0),hessian=TRUE)
fec.lik
covmat-solve(CV.lik$hessian)
covmat
Ruben
__
R-help@stat.math.ethz.ch 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] cox.zph
Dear R-users, I am sorry if this is obvious. I am testing the proportional hazard assumptions using cox.zph. If i am not wrong, a g(t) function must be assumed. Four possibilities available in R are km,identity and rank. may i know what functions of time are these transformation assuming? Thanks a lot in advance for your wisdom. kind regards, sing yee ling [[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] panel function with barchart (lattice)
On 1/25/06, Drew [EMAIL PROTECTED] wrote:
Folks at R help,
I can't quite get the panel function to work the way I
want within barchart.
I guess I'm still not understanding how to piece
together multiple panel
arguments, especially when groups is specified.
Example: I want to be able to add the value of yield
to each section of
each bar in this graph:
barchart(yield ~ variety | site, data = barley,
groups = year,
layout = c(1,6),
stack=TRUE,
ylab = Barley Yield (bushels/acre)
)
To do this, I add my panel function:
barchart(yield ~ variety | site, data = barley,
groups = year,
layout = c(1,6),
stack=TRUE,
ylab = Barley Yield (bushels/acre),
panel = function(x,y,subscripts,groups,...){
panel.barchart(x,y,...)
Well, panel.barchart needs the subscripts and groups arguments to draw
stacked bar charts, and you are calling it without them.
ltext(x = x, y = y, label =
round(barley$yield[subscripts],1), cex=.8)
The y values will need to be accumulated. Have you looked at what
panel.barchart does?
}
)
Then I get the values to print on each bar (which is
what I want) but the
bars no longer stack to appropriate height, and I
cannot get the subsections
of each bar to be a different color. I've tried
numerous variations of
panel.barchart, panel.superpose, etc. using examples
from ?xyplot, but
nothing quite works or I get an error message.
Any help would be appreciated.
~Nick
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
--
http://www.stat.wisc.edu/~deepayan/
__
R-help@stat.math.ethz.ch 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] cox.zph
On Wed, 25 Jan 2006, singyee ling wrote: Dear R-users, I am sorry if this is obvious. I am testing the proportional hazard assumptions using cox.zph. If i am not wrong, a g(t) function must be assumed. Four possibilities available in R are km,identity and rank. may i know what functions of time are these transformation assuming? It should be fairly clear what identity, log and rank are: g(t)=t g(t)=log(t) and the rank of the observations time, respectively. km uses the Kaplan-Meier estimator of survival at time t. -thomas __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] reducing learning curves?
In addition to the other suggestions you may want to look at JGR (http://www.rosuda.org/JGR/). It does most of what you asked (except debugging, use the debug package for that). -Original Message- From: [EMAIL PROTECTED] on behalf of Michael Sent: Wed 1/25/2006 1:09 AM To: R-help@stat.math.ethz.ch Subject: [R] reducing learning curves? Hi all, I am really new to the R language. I am a long time Matlab and C++ user and I was forced to learn R because I am taking a statistics class. I am seeking to reduce the learning curve to as smooth as possible. Are there any addon/plug-in features that can reduce the learning curve, for example, the following features can be very helpful for new learners: 1. Matlab-like command line auto-completion: Matlab has huge amount of command and nobody is able to remember them off the head. So a nice feature of Matlab command line is that I just need to type the first a few letters and then I press TAB key, there will be a list of possible commands popping up so I just need to select one. This helps a lot in terms of learning for new comers. A more advanced command auto-completion is Visual C++-like, which is implemented in program editor. It helps a lot while doing programming; 2. A good IDE editor with embedded inline debugger: can be as good as VC++, but also can be as simple as Matlab's debugger, which can breakpoint and trace line-by-line... the editor can do syntax correction, syntax check, syntax highlighting, code formatting, etc. Could you please recommend some good addon/plugins that have the above features? Could you please also suggest some tips/tools/tricks that can help me reduce the learning curve? Thank you very much! Michael. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] how to test robustness of correlation
Hi, there: As you all know, correlation is not a very robust procedure. Sometimes correlation could be driven by a few outliers. There are a few ways to improve the robustness of correlation (pearson correlation), either by outlier removal procedure, or resampling technique. I am wondering if there is any R package or R code that have incorporated outlier removal or resampling procedure in calculating correlation coefficient. Your help is greatly appreciated. Thanks. Yang Yang Qiu Integrated Data Analysis [EMAIL PROTECTED] GlaxoSmithKline [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] a nice tool
Have you checked out the Zelig package for Windows?It simplifies using R for those who are less than enthusiastic about the using a command line environment. It is a framework for using the underlying R packages already existing. The manual has a nice appendix on getting up to speed in R. Check it out at http://gking.harvard.edu/zelig/ or download from CRAN. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to test robustness of correlation
check out cov.rob() in MASS (among others, I'm sure). The procedure is far more sophisticated than outlier removal or resampling (??). References are given in the docs. -- 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 [EMAIL PROTECTED] Sent: Wednesday, January 25, 2006 12:37 PM To: r-help@stat.math.ethz.ch Subject: [R] how to test robustness of correlation Hi, there: As you all know, correlation is not a very robust procedure. Sometimes correlation could be driven by a few outliers. There are a few ways to improve the robustness of correlation (pearson correlation), either by outlier removal procedure, or resampling technique. I am wondering if there is any R package or R code that have incorporated outlier removal or resampling procedure in calculating correlation coefficient. Your help is greatly appreciated. Thanks. Yang Yang Qiu Integrated Data Analysis [EMAIL PROTECTED] GlaxoSmithKline [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] reducing learning curves?
Michael comtech.usa at gmail.com writes: Hi all, I am really new to the R language. I am a long time Matlab and C++ user and I was forced to learn R because I am taking a statistics class. I am seeking to reduce the learning curve to as smooth as possible. This cheatsheet might be helpful for a Matlab user: http://cran.us.r-project.org/doc/contrib/R-and-octave.txt and help.search(keyword) and RSiteSearch(terms of intererst) are very useful in finding tips. [ rest deleted ] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] D(dnorm...)?
Can someone help me understand the following: D(expression(dnorm(x, mean)), mean) [1] 0 sessionInfo() R version 2.2.1, 2005-12-20, i386-pc-mingw32 attached base packages: [1] methods stats graphics grDevices utils datasets [7] base By my computations, this should be something like ((mean-x)/sd^2)*dnorm(...). Thanks for your help. Spencer Graves __ R-help@stat.math.ethz.ch 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] D(dnorm...)?
dnorm() is an internal function, so I don't see how D (or deriv) can do anything with it symbolically. Am I missing something? -- Bert -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Spencer Graves Sent: Wednesday, January 25, 2006 2:43 PM To: r-help@stat.math.ethz.ch Subject: [R] D(dnorm...)? Can someone help me understand the following: D(expression(dnorm(x, mean)), mean) [1] 0 sessionInfo() R version 2.2.1, 2005-12-20, i386-pc-mingw32 attached base packages: [1] methods stats graphics grDevices utils datasets [7] base By my computations, this should be something like ((mean-x)/sd^2)*dnorm(...). Thanks for your help. Spencer Graves __ R-help@stat.math.ethz.ch 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] D(dnorm...)?
Hi, Bert: I think I was too terse: Why didn't I get an error message? When I tried the same thing with dpois, I got an error message: D(expression(dpois(x, prob)), mean) Error in D(expression(dpois(x, prob)), mean) : Function 'dpois' is not in the derivatives table With dnorm, I got 0!. Whe didn't I get an error message? Thanks for the reply. Spencer Graves Berton Gunter wrote: dnorm() is an internal function, so I don't see how D (or deriv) can do anything with it symbolically. Am I missing something? -- Bert -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Spencer Graves Sent: Wednesday, January 25, 2006 2:43 PM To: r-help@stat.math.ethz.ch Subject: [R] D(dnorm...)? Can someone help me understand the following: D(expression(dnorm(x, mean)), mean) [1] 0 sessionInfo() R version 2.2.1, 2005-12-20, i386-pc-mingw32 attached base packages: [1] methods stats graphics grDevices utils datasets [7] base By my computations, this should be something like ((mean-x)/sd^2)*dnorm(...). Thanks for your help. Spencer Graves __ R-help@stat.math.ethz.ch 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] D(dnorm...)?
On Wed, 25 Jan 2006, Spencer Graves wrote: Hi, Bert: I think I was too terse: Why didn't I get an error message? When I tried the same thing with dpois, I got an error message: D(expression(dpois(x, prob)), mean) Error in D(expression(dpois(x, prob)), mean) : Function 'dpois' is not in the derivatives table With dnorm, I got 0!. Whe didn't I get an error message? Because dnorm _is_ in the derivatives table, as a function of x. (So is pnorm.) These are functions we write as \phi and \Phi in Statistics 101. This is not particularly well documented in R, but it is in MASS4 pp. 437-8. It would be worth enhancing the help page with the list of known functions. Long ago we talked about making it extensible, but nothing happened. Thanks for the reply. Spencer Graves Berton Gunter wrote: dnorm() is an internal function, so I don't see how D (or deriv) can do anything with it symbolically. Am I missing something? -- Bert -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Spencer Graves Sent: Wednesday, January 25, 2006 2:43 PM To: r-help@stat.math.ethz.ch Subject: [R] D(dnorm...)? Can someone help me understand the following: D(expression(dnorm(x, mean)), mean) [1] 0 sessionInfo() R version 2.2.1, 2005-12-20, i386-pc-mingw32 attached base packages: [1] methods stats graphics grDevices utils datasets [7] base By my computations, this should be something like ((mean-x)/sd^2)*dnorm(...). Thanks for your help. Spencer Graves __ R-help@stat.math.ethz.ch 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 -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] D(dnorm...)?
Hi, Prof. Ripley: Thanks for the explanation. If I had read your book more carefully, I would not have needed this email exchange. Thanks again, Spencer Graves Prof Brian Ripley wrote: On Wed, 25 Jan 2006, Spencer Graves wrote: Hi, Bert: I think I was too terse: Why didn't I get an error message? When I tried the same thing with dpois, I got an error message: D(expression(dpois(x, prob)), mean) Error in D(expression(dpois(x, prob)), mean) : Function 'dpois' is not in the derivatives table With dnorm, I got 0!. Whe didn't I get an error message? Because dnorm _is_ in the derivatives table, as a function of x. (So is pnorm.) These are functions we write as \phi and \Phi in Statistics 101. This is not particularly well documented in R, but it is in MASS4 pp. 437-8. It would be worth enhancing the help page with the list of known functions. Long ago we talked about making it extensible, but nothing happened. Thanks for the reply. Spencer Graves Berton Gunter wrote: dnorm() is an internal function, so I don't see how D (or deriv) can do anything with it symbolically. Am I missing something? -- Bert -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Spencer Graves Sent: Wednesday, January 25, 2006 2:43 PM To: r-help@stat.math.ethz.ch Subject: [R] D(dnorm...)? Can someone help me understand the following: D(expression(dnorm(x, mean)), mean) [1] 0 sessionInfo() R version 2.2.1, 2005-12-20, i386-pc-mingw32 attached base packages: [1] methods stats graphics grDevices utils datasets [7] base By my computations, this should be something like ((mean-x)/sd^2)*dnorm(...). Thanks for your help. Spencer Graves __ R-help@stat.math.ethz.ch 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
[R] imbalanced classes
Hi Andy, I know this topic has been discussed before on the R-help, but I was wondering if you could offer some advice specific to my application. I'm using the R random forest package to compare two classes of data, the number of cases in each class relatively low, 28 in class 1 and 9 in class 2. I'd really like to use R environment to analyze this data, however I'm finding it difficult to put much trust in the results of my analysis. As you've stated, the classwt variables do not do much, and I've tried working with the cuttoff and sampsize variables as well, with limited success in balancing error rates between the two classes. It was unclear to me how to use the cuttoff parameter correctly. If you have any recommendations here, it would be appreciated. Additionally with the sampsize variable, I have tried a few values, for example setting sampsize = c(2, 6) and c(9, 3), etc. It wasn't clear to me if I should be sampling more from the larger class or the other way around. Lastly, I'm wondering if you are currently working or have plans to release in the near future an R version of randomForest that is equivalent to the FORTRAN rf5 package. It works wonderfully for my application, but getting data in and out of it, changing parameters, compiling is just a pain, as I'm sure you agree. Your thoughts would be greatly appreciated. Kind regards, Mark D'Ascenzo Biomedical Engineering Cornell University Ithaca, NY 14853 __ R-help@stat.math.ethz.ch 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 ANCOVA: still confused
Harold, Kingsford, and R-users, I settled on using the lmer function. I think the memory issue was more a function of my poor coding than an actual memory problem. I also switched the label from box to clutch to avoid any potential confusion with other functions. This coding seems to have worked: eabl - lmer(rtot ~ sexv + (purban2|clutch), maxIter=1000, data=bb, na.action=na.omit) However, I have two remaining questions: (1)how concerned should I be with the warning message below and (2) is there a way to invoke output to get an estimate of the effect of purban2 (the proportion of urban cover 200 m around a box) on feather color (rtot) and if there is a difference between the sexes? I used the summary function and it doesn't tell me much (see output below). I'll read up mixed models when Pinheiro arrives but any suggestions for diagnostics? I'm going to repeat this study and expand it by doubling or tripling the number of birds. Warning message: nlminb returned message false convergence (8) in: LMEoptimize-(`*tmp*`, value = list(maxIter = 200, tolerance = 1.49011611938477e-08, summary(eabl) Linear mixed-effects model fit by REML Formula: rtot ~ sexv + (purban2 | clutch) Data: bb AIC BIClogLik MLdeviance REMLdeviance 5164.284 6997.864 -2052.142 4128.792 4104.284 Random effects: Groups Name Variance Std.Dev. Corr clutch (Intercept) 502829 709.10 purban20 1341990 1158.44 -0.477 purban20.006711409 5683957 2384.10 -0.226 0.082 purban20.013422821772922 1331.51 -0.386 0.176 0.067 . . . . . # of obs: 235, groups: clutch, 74 Fixed effects: Estimate Std. Error t value (Intercept) 5950.01 241.59 24.628 sexvm1509.07 145.73 10.355 Correlation of Fixed Effects: (Intr) sexvm -0.304 Thanks many time over, Jeff Jeffrey A. Stratford, Ph.D. Postdoctoral Associate 331 Funchess Hall Department of Biological Sciences Auburn University Auburn, AL 36849 334-329-9198 FAX 334-844-9234 http://www.auburn.edu/~stratja Doran, Harold [EMAIL PROTECTED] 01/25/06 6:37 AM OK, we're getting somewhere. First, it looks as though (by the error message) that you have a big dataset. My first recommendation is to use lmer instead of lme, you will see a significant benefit in terms of computional speed. For the model this would be lmer(rtot ~ sexv +(purban|box:chick) + (purban|box), bb, na.action=na.omit) Now, you have run out of memory. I don't know what operating system you are using, so go and see the appropriate FAQ for increasing memory for your OS. Second, I made a mistake in my reply. Your random statement should be random=~purban|box/chick denoting that chicks are nested in boxes, not boxes nested in chicks, sorry about that. Now, why is it that each chick within box has the same value for purban? If this is so, why are you fitting that as a random effect? It seems not to vary across individual chicks, right? It seems there is only an effect of box and not an effect for chicks. Why not just fit a random effect only for box such as: rtot.lme - lme(fixed=rtot~sexv, random=~purban2|box, na.action=na.omit,bb) or in lmer lmer(rtot ~ sexv + (purban|box), bb, na.action=na.omit) Harold -Original Message- From: Jeffrey Stratford [mailto:[EMAIL PROTECTED] Sent: Tue 1/24/2006 8:57 PM To: Doran, Harold; r-help@stat.math.ethz.ch Cc: Subject:RE: [R] nested ANCOVA: still confused R-users and Harold. First, thanks for the advice; I'm almost there. The code I'm using now is library(nlme) bb - read.csv(E:\\eabl_feather04.csv, header=TRUE) bb$sexv - factor(bb$sexv) rtot.lme - lme(fixed=rtot~sexv, random=~purban2|chick/box, na.action=na.omit, data=bb) A sample of the data looks like this box chick rtotpurban2 sexv 1 1 6333.51 0.026846f 1 2 8710.8840.026846m 2 1 5810.0070.161074f 2 2 5524.33 0.161074f 2 3 4824.4740.161074f 2 4 5617.6410.161074f 2 5 6761.7240.161074f 4 1 7569.6730.208054m 4 2 7877.0810.208054m 4 4 7455.55 0.208054f 7 1 5408.2870.436242m 10 1 6991.7270.14094 f 12 1 8590.2070.134228f 12 2 7536.747
Re: [R] nested ANCOVA: still confused
Jeff: However, I have two remaining questions: (1)how concerned should I be with the warning message below and There was a definitive comment on this just a few days ago on the list (search the archives), the gist of it was: **very concerned** . False convergence means that you're not truly converged, the details for which I've forgotten (sigh...). Anyway, this means that your parameter estimates could be far from the correct minimized values = you could be in trouble. I can't help you any more than that, but hopefully you'll get responses from those with suitable expertise who can. Cheers, Bert (2) is there a way to invoke output to get an estimate of the effect of purban2 (the proportion of urban cover 200 m around a box) on feather color (rtot) and if there is a difference between the sexes? I used the summary function and it doesn't tell me much (see output below). I'll read up mixed models when Pinheiro arrives but any suggestions for diagnostics? I'm going to repeat this study and expand it by doubling or tripling the number of birds. Warning message: nlminb returned message false convergence (8) in: LMEoptimize-(`*tmp*`, value = list(maxIter = 200, tolerance = 1.49011611938477e-08, summary(eabl) Linear mixed-effects model fit by REML Formula: rtot ~ sexv + (purban2 | clutch) Data: bb AIC BIClogLik MLdeviance REMLdeviance 5164.284 6997.864 -2052.142 4128.792 4104.284 Random effects: Groups Name Variance Std.Dev. Corr clutch (Intercept) 502829 709.10 purban20 1341990 1158.44 -0.477 purban20.006711409 5683957 2384.10 -0.226 0.082 purban20.013422821772922 1331.51 -0.386 0.176 0.067 . . . . . # of obs: 235, groups: clutch, 74 Fixed effects: Estimate Std. Error t value (Intercept) 5950.01 241.59 24.628 sexvm1509.07 145.73 10.355 Correlation of Fixed Effects: (Intr) sexvm -0.304 Thanks many time over, Jeff Jeffrey A. Stratford, Ph.D. Postdoctoral Associate 331 Funchess Hall Department of Biological Sciences Auburn University Auburn, AL 36849 334-329-9198 FAX 334-844-9234 http://www.auburn.edu/~stratja Doran, Harold [EMAIL PROTECTED] 01/25/06 6:37 AM OK, we're getting somewhere. First, it looks as though (by the error message) that you have a big dataset. My first recommendation is to use lmer instead of lme, you will see a significant benefit in terms of computional speed. For the model this would be lmer(rtot ~ sexv +(purban|box:chick) + (purban|box), bb, na.action=na.omit) Now, you have run out of memory. I don't know what operating system you are using, so go and see the appropriate FAQ for increasing memory for your OS. Second, I made a mistake in my reply. Your random statement should be random=~purban|box/chick denoting that chicks are nested in boxes, not boxes nested in chicks, sorry about that. Now, why is it that each chick within box has the same value for purban? If this is so, why are you fitting that as a random effect? It seems not to vary across individual chicks, right? It seems there is only an effect of box and not an effect for chicks. Why not just fit a random effect only for box such as: rtot.lme - lme(fixed=rtot~sexv, random=~purban2|box, na.action=na.omit,bb) or in lmer lmer(rtot ~ sexv + (purban|box), bb, na.action=na.omit) Harold -Original Message- From: Jeffrey Stratford [mailto:[EMAIL PROTECTED] Sent: Tue 1/24/2006 8:57 PM To: Doran, Harold; r-help@stat.math.ethz.ch Cc: Subject: RE: [R] nested ANCOVA: still confused R-users and Harold. First, thanks for the advice; I'm almost there. The code I'm using now is library(nlme) bb - read.csv(E:\\eabl_feather04.csv, header=TRUE) bb$sexv - factor(bb$sexv) rtot.lme - lme(fixed=rtot~sexv, random=~purban2|chick/box, na.action=na.omit, data=bb) A sample of the data looks like this box chick rtotpurban2 sexv 1 1 6333.51 0.026846f 1 2 8710.8840.026846m 2 1 5810.0070.161074f 2 2 5524.33 0.161074f 2 3 4824.4740.161074f 2 4 5617.6410.161074f 2 5 6761.7240.161074f 4 1 7569.6730.208054m 4 2 7877.081
[R] panel.xyplot : incorrectly connecting points
R 2.2, WinXP. I am having problems getting the right kind of
xyplot( ) to be generated. The first of these works fine, but
doesn't overlay a reference grid (which I need):
xyplot(Y ~ X | Factor1, type = 'b', groups = GROUP,
col = c(1,13), pch = c(16,6), lty = 1, lwd = 2,
cex = 1.2, data = FOO.Frame,
between = list(x = .5, y = .5),
scales = list(alternating = TRUE))
The second of these displays the grid as I need, but incorrectly
connects the points from _different_ GROUP values within each
panel. I have made sure that GROUP is an ordered factor:
xyplot(PROB ~ MEAN | SD, data = SimProb,
groups = GROUP,
between = list(x = .5, y = .5),
scales = list(alternating = TRUE),
panel = function(x, y, ...)
{
panel.grid(h=-1, v=-1, lwd = 1)
panel.xyplot(x, y, type = 'b', col = 1,
lwd = 2, cex = 1.2, ...)
}
)
I am sure that this sort of question has been asked before, and
I apologize for any redundancy. I would appreciate any help.
Kind regards,
Greg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] panel.xyplot : incorrectly connecting points
On 1/25/06, Greg Tarpinian [EMAIL PROTECTED] wrote:
R 2.2, WinXP. I am having problems getting the right kind of
xyplot( ) to be generated. The first of these works fine, but
doesn't overlay a reference grid (which I need):
xyplot(Y ~ X | Factor1, type = 'b', groups = GROUP,
col = c(1,13), pch = c(16,6), lty = 1, lwd = 2,
cex = 1.2, data = FOO.Frame,
between = list(x = .5, y = .5),
scales = list(alternating = TRUE))
The simplest solution is to use
type = c('b', 'g')
The second of these displays the grid as I need, but incorrectly
connects the points from _different_ GROUP values within each
panel. I have made sure that GROUP is an ordered factor:
xyplot(PROB ~ MEAN | SD, data = SimProb,
groups = GROUP,
between = list(x = .5, y = .5),
scales = list(alternating = TRUE),
panel = function(x, y, ...)
{
panel.grid(h=-1, v=-1, lwd = 1)
panel.xyplot(x, y, type = 'b', col = 1,
lwd = 2, cex = 1.2, ...)
}
)
This should work in the devel version of R (2.3.0 to be). The reason
it doesn't work in 2.2.x is that the default panel function when
'groups' in non-null is not panel.xyplot, but rather panel.superpose.
(This is admittedly confusing, hence the change for 2.3.0.)
Deepayan
--
http://www.stat.wisc.edu/~deepayan/
__
R-help@stat.math.ethz.ch 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] regression with nestedness
Have you considered lme in library(nlme)? The companion book Pinheiro and Bates (2000) Mixed-Effects Models in S and S-Plus (Springer) is my favorite reference for this kind of thing. From what I understand of your question, you should be able to find excellent answers in this book. [You could also try lmer associated with the lme4 package. However, many of the helper functions needed to do what you want are available for lme but not yet for lmer; lmer can do things that lme can not, but it doesn't sound like you need any of those extra features. Even if you did, I would suggest you first ignore those issues and see what you can get from lme.] If this is not adequate and you would like more help from this list, please submit another post. Before you do, however, PLEASE do read the posting guide! www.R-project.org/posting-guide.html, especially the part about providing a brief but self-contained toy example illustrating what you tried to do and why it was not adequate. Jeffrey Stratford wrote: Dear R-users, I set up an experiment where I put up bluebird boxes across an urbanization gradient. I monitored these boxes and at some point I pulled a feather from a chick and a friend used spectral properties (rtot, a continuous var) to index chick health. There is an effect of sex that I would like to include but how would I set up a regression and look at the effect of urbanization (purban, a continuous var)) on feather properties of chicks within boxes. So the model should look something like rtot = sex + purban + (chick)clutch Also, when I plot purban against rtot using the plot function I get boxplots but I would like to ignore the clutch and just plot each point. I've tried type = p but this has no effect. Thanks, Jeff Jeffrey A. Stratford, Ph.D. Postdoctoral Associate 331 Funchess Hall Department of Biological Sciences Auburn University Auburn, AL 36849 334-329-9198 FAX 334-844-9234 http://www.auburn.edu/~stratja __ R-help@stat.math.ethz.ch 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] Partial Canonical Correlation
Hello, I am interested in doing a partial canonical correlation (identical to the SAS function, Proc Cancorr with the Partial statement). By this I mean, I have 3 sets of data, a vegetation matrix (columns of abundances of species in rows of plots), an environment matrix (columns of environmental variables in same rows of plots), and a space matrix (x and y locations of each of the plots). I would like to look at the relationships between these multivariate groups, with the space partialed out. This would be somewhat analogous in my mind to a partial mantel test (i.e., 3 data matrices), but would be somewhat more sophisticated than the correlations between distance matrices of the partial mantel test. In the r-help archives I have found similar things, partial correlations, canonical correlations, but not partial canonical correlations. Does somebody know how to do this in R? Thank you in advance, Eliot -- Eliot McIntire Post Doctoral Fellow Department of Ecosystems and Conservation Science College of Forestry and Conservation University of Montana, Missoula, MT 59812 406-243-5239 fax: 406-243-4557 [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] A function slightly different from diff()
Hi,
I wonder if the following function has already been implemented in (some) R
(package):
summ - function(x, lag=1) # x is a vector
{
n - length(x)
x[(1+lag):n] + x[1:(n-lag)]
}
which (I think) the only difference from diff() is + instead of -, for
differences=1.
Thank you,
Dimitri Szerman
__
R-help@stat.math.ethz.ch 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 function slightly different from diff()
This isn't a single function but its a simple expression:
x - ts(1:10) # test data
x + lag(x)
On 1/25/06, [EMAIL PROTECTED] [EMAIL PROTECTED] wrote:
Hi,
I wonder if the following function has already been implemented in (some) R
(package):
summ - function(x, lag=1) # x is a vector
{
n - length(x)
x[(1+lag):n] + x[1:(n-lag)]
}
which (I think) the only difference from diff() is + instead of -, for
differences=1.
Thank you,
Dimitri Szerman
__
R-help@stat.math.ethz.ch 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] imbalanced classes
Mark, I guess the message is meant for me (yet you sent it to R-help). If you have 10 class A and 100 class B, not setting sampsize would cause a random sample (with replacement) of 110 from the whole lot, which, of course, would give you on the average 10 times more Bs than As in the sample. If you grow a tree on such a sample, it's not going to do so well in predicting the As. However, if you set sampsize=c(10, 10), then each tree is grown on 10 randomly sampled As and 10 randomly sampled Bs, giving the tree a much better chance of giving roughly similar error rates for predicting As and Bs. If setting the sampsize to be equal doesn't quite do it, you can try setting it to the more extreme direction. As to cutoff, in a two-class problem, it's the same as setting the classification threshold to something other than 0.5. E.g., if cutoff=c(0.9, 0.1), then a case with 80% of the votes for class A would still be classified as B, because .8/.9 .2/.1. Hope that's clear. I do have to wonder, though, if you only have a total of 37 cases in the data, how can you be sure the estimates of class error rates you get will pan out on a larger test set? I would think the variability on the estimate of the class error rates is so high that it doesn't make too much sense to try to balance them too much... Just my $0.02. I do plan on implementing the weighted RF (see the To Do part of rfNews()), but don't hold your breath... Cheers, Andy From: Mark D'Ascenzo Hi Andy, I know this topic has been discussed before on the R-help, but I was wondering if you could offer some advice specific to my application. I'm using the R random forest package to compare two classes of data, the number of cases in each class relatively low, 28 in class 1 and 9 in class 2. I'd really like to use R environment to analyze this data, however I'm finding it difficult to put much trust in the results of my analysis. As you've stated, the classwt variables do not do much, and I've tried working with the cuttoff and sampsize variables as well, with limited success in balancing error rates between the two classes. It was unclear to me how to use the cuttoff parameter correctly. If you have any recommendations here, it would be appreciated. Additionally with the sampsize variable, I have tried a few values, for example setting sampsize = c(2, 6) and c(9, 3), etc. It wasn't clear to me if I should be sampling more from the larger class or the other way around. Lastly, I'm wondering if you are currently working or have plans to release in the near future an R version of randomForest that is equivalent to the FORTRAN rf5 package. It works wonderfully for my application, but getting data in and out of it, changing parameters, compiling is just a pain, as I'm sure you agree. Your thoughts would be greatly appreciated. Kind regards, Mark D'Ascenzo Biomedical Engineering Cornell University Ithaca, NY 14853 __ R-help@stat.math.ethz.ch 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] D(dnorm...)?
Yes Bert, this time you are missing something (unusually) ... As Brian Ripley pointed out 'dnorm' is in the derivative table, *but* only as a function of one variable. So if you want to find the derivative of dnorm(x, mean, sigma) you have to write it as 1/sigma * dnorm((x - mu)/sigma). Here is a little example: D(Quote(pnorm((x-mu)/sigma)), x) dnorm((x - mu)/sigma) * (1/sigma) D(D(Quote(pnorm((x-mu)/sigma)), x), mu) (x - mu)/sigma * (dnorm((x - mu)/sigma) * (1/sigma)) * (1/sigma) --- Like Brian, I recall the suggestion that we make D(...) extensible. I still think it is a good idea and worth considering. Under one scheme you would specify an object such as Fnorm - structure(quote(pnorm(x, mu, sigma)), deriv = list(x = Quote(dnorm(x, mu, sigma)/sigms), mu = Quote(-dnorm(x, mu, sigma)/sigma), sigma = Quote(-(x - mu)*dnorm(x, mu, sigma)/sigma^2), class = dfunction) ane write a generic differentiate function with a dfunction method and D as the default. I don't think it's quite that easy, but the plan is clear enough. Bill. -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Berton Gunter Sent: Thursday, 26 January 2006 8:58 AM To: 'Spencer Graves'; r-help@stat.math.ethz.ch Subject: Re: [R] D(dnorm...)? dnorm() is an internal function, so I don't see how D (or deriv) can do anything with it symbolically. Am I missing something? -- Bert -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Spencer Graves Sent: Wednesday, January 25, 2006 2:43 PM To: r-help@stat.math.ethz.ch Subject: [R] D(dnorm...)? Can someone help me understand the following: D(expression(dnorm(x, mean)), mean) [1] 0 sessionInfo() R version 2.2.1, 2005-12-20, i386-pc-mingw32 attached base packages: [1] methods stats graphics grDevices utils datasets [7] base By my computations, this should be something like ((mean-x)/sd^2)*dnorm(...). Thanks for your help. Spencer Graves __ R-help@stat.math.ethz.ch 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] D(dnorm...)?
Hello, Bill: I'm not qualified to make this suggestion since I'm incapable of turning it into reality, but what about creating a link between R and one of the Mathematica clones like Yacas? I can immagine that it could be substantially more difficult than linking R to other software like Excel, but ... . Spencer Graves [EMAIL PROTECTED] wrote: Yes Bert, this time you are missing something (unusually) ... As Brian Ripley pointed out 'dnorm' is in the derivative table, *but* only as a function of one variable. So if you want to find the derivative of dnorm(x, mean, sigma) you have to write it as 1/sigma * dnorm((x - mu)/sigma). Here is a little example: D(Quote(pnorm((x-mu)/sigma)), x) dnorm((x - mu)/sigma) * (1/sigma) D(D(Quote(pnorm((x-mu)/sigma)), x), mu) (x - mu)/sigma * (dnorm((x - mu)/sigma) * (1/sigma)) * (1/sigma) --- Like Brian, I recall the suggestion that we make D(...) extensible. I still think it is a good idea and worth considering. Under one scheme you would specify an object such as Fnorm - structure(quote(pnorm(x, mu, sigma)), deriv = list(x = Quote(dnorm(x, mu, sigma)/sigms), mu = Quote(-dnorm(x, mu, sigma)/sigma), sigma = Quote(-(x - mu)*dnorm(x, mu, sigma)/sigma^2), class = dfunction) ane write a generic differentiate function with a dfunction method and D as the default. I don't think it's quite that easy, but the plan is clear enough. Bill. -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Berton Gunter Sent: Thursday, 26 January 2006 8:58 AM To: 'Spencer Graves'; r-help@stat.math.ethz.ch Subject: Re: [R] D(dnorm...)? dnorm() is an internal function, so I don't see how D (or deriv) can do anything with it symbolically. Am I missing something? -- Bert -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Spencer Graves Sent: Wednesday, January 25, 2006 2:43 PM To: r-help@stat.math.ethz.ch Subject: [R] D(dnorm...)? Can someone help me understand the following: D(expression(dnorm(x, mean)), mean) [1] 0 sessionInfo() R version 2.2.1, 2005-12-20, i386-pc-mingw32 attached base packages: [1] methods stats graphics grDevices utils datasets [7] base By my computations, this should be something like ((mean-x)/sd^2)*dnorm(...). Thanks for your help. Spencer Graves __ R-help@stat.math.ethz.ch 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
