[R] Re: estimate the parameter of exponential distribution, etc.
Hi David, You can estimate parameters using fitdistr() in MASS package. I suggest also to read my contribute Fitting distributions with R available on CRAN: http://cran.r-project.org/doc/contrib/Ricci-distributions-en.pdf Best regards, Vito you wrote: Given a numeric vector of observations, does R have any generic way to estimate the parameters of commonly used distributions (exponential, gamma, etc.) without numerically optimizing the likelihood function? Thanks, David ___ David R. Bickel http://davidbickel.com Research Scientist Pioneer Hi-Bred International Bioinformatics Exploratory Research 7250 NW 62nd Ave., PO Box 552 Johnston, Iowa 50131-0552 515-334-4739 Tel 515-334-6634 Fax david.bickel at pioneer.com, bickel at prueba.info = Diventare costruttori di soluzioni Became solutions' constructors The business of the statistician is to catalyze the scientific learning process. George E. P. Box Top 10 reasons to become a Statistician 1. Deviation is considered normal 2. We feel complete and sufficient 3. We are 'mean' lovers 4. Statisticians do it discretely and continuously 5. We are right 95% of the time 6. We can legally comment on someone's posterior distribution 7. We may not be normal, but we are transformable 8. We never have to say we are certain 9. We are honestly significantly different 10. No one wants our jobs Visitate il portale http://www.modugno.it/ e in particolare la sezione su Palese http://www.modugno.it/archivio/palese/ __ R-help@stat.math.ethz.ch 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] Re: nonlinear least square fit of an unknown function
Hi Eric, if I understand you question, are you trying to fit an unknown distribution? Are looking for computing theorical frequencies (area under the curve)? In this case you could see: http://cran.r-project.org/doc/contrib/Ricci-distributions-en.pdf Use hist() and density() to identify your distribution and then fit parameters. Another way could be using splines (?spline). Hoping I helped you. Cordially Vito you wrote: Hi, I have a set of twelve points and wonder how I can get a function that can then be used to calculate the area under the curve (most important). Thanks. Eric = Diventare costruttori di soluzioni Became solutions' constructors The business of the statistician is to catalyze the scientific learning process. George E. P. Box Top 10 reasons to become a Statistician 1. Deviation is considered normal 2. We feel complete and sufficient 3. We are 'mean' lovers 4. Statisticians do it discretely and continuously 5. We are right 95% of the time 6. We can legally comment on someone's posterior distribution 7. We may not be normal, but we are transformable 8. We never have to say we are certain 9. We are honestly significantly different 10. No one wants our jobs Visitate il portale http://www.modugno.it/ e in particolare la sezione su Palese http://www.modugno.it/archivio/palese/ __ R-help@stat.math.ethz.ch 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] Memory error in Mac OS X Aqua GUI v1.01 with cluster package functions
The problem is probably that your 32-bit application has become fragmented in its address space. (I have been told that 64-bit builds of R on MacOS X will need 10.4, and the binary on CRAN is definitely 32-bit.) The distance matrix for your object alone takes 800Mb, so you will struggle on any 32-bit OS, let alone on a machine with `only' 1Gb RAM. However, using R's dist() rather than cluster's daisy() will need less memory for this task. This is R-help not the Bioconductor list, and outside bioinformatics the prevailing opinion is that clustering 10,000 cases is not sensible. But there are methods of clustering that do not use the distance matrix, e.g. clara() in the (unmentioned) package cluster than provides daisy(). On Tue, 22 Feb 2005, Betty Gilbert wrote: I'm sorry if the answer to my problem is buried in the archives. I have limited experience with R and I couldn't find a solution to my particular problem. I am running Mac OS X Aqua GUI v1.01 on a new G5 running os 10.3.8 with a 1.8Ghz processor and 1GB of sdram. I just downloaded bioconducter a week ago and I'm trying to cluster a matrix I created with a simulation with dimensions dim(nca35) [1] 1048112 with size object.size(nca352) [1] 1426204 I checked my ulimits variable on the unix terminal and it says it's unlimited as does But your machine and OS impose limits that ulimits is not telling you. mem.limits() nsize vsize NANA But I'm still getting errors like the following with funtions in the cluster package daisy(nca352, metric= euclidean, stand=FALSE)-dnca35 Error: cannot allocate vector of size 858213 Kb *** malloc: vm_allocate(size=878813184) failed (error code=3) *** malloc[599]: error: Can't allocate region if it helps i also checked gc() used (Mb) gc trigger (Mb) Ncells 448662 12.0 741108 19.8 Vcells 847630 6.5 135357901 1032.7 I tried the suggested unix command in the memory help doc but that doesn't work in the Aqua GUI. Can someone tell me how to change the Vcells? Although to the best of my understanding (which is limited) I shouldn't have to do that. Any suggestions would be greatly appreciated. thanks, betty -- Betty Gilbert [EMAIL PROTECTED] Taylor Lab Plant and Microbial Biology 321 Koshland Hall U.C. Berkeley Berkeley, Ca 94720 -- 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] Memory problems
R does use virtual memory, and memory.size() (Windows only) is documented to report the usage, not change the limit. Please read that help page more carefully. If you are not already doing so, try arima0 rather than arima. And do see the posting guide! On Tue, 22 Feb 2005, Konstantinos Kleisouris wrote: I use R to do some ARIMA forecasting and R runs out of memory. The problem is that I have 20160 samples(which are quite alot) and when I try to fit the model it runs out of memory. I tried with memory.size() to change the limit, but it wouldn't work. Is there anything you can suggest? Is it possible R can use virtual memory? 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
[R] Re: Does R has the function for garch-t, gjr-garch, qgarch and egarch
Hi, see: - HestonNandiGarchFit(fOptions) Heston-Nandi Garch(1,1) Modelling - HestonNandiOptions(fOptions) Option Price for the Heston-Nandi Garch Option Model - GarchModelling(fSeries) Univariate GARCH Time Series Modelling - GarchDistributionFits(fSeries) Parameter Fit of a Distribution - garch(tseries) Fit GARCH Models to Time Series Best regards Vito You wrote: Dear all, I would like to know that R has the function for garch-t,gjr- garch,qgarch and egarch. Best Regards, Luck = Diventare costruttori di soluzioni Became solutions' constructors The business of the statistician is to catalyze the scientific learning process. George E. P. Box Top 10 reasons to become a Statistician 1. Deviation is considered normal 2. We feel complete and sufficient 3. We are 'mean' lovers 4. Statisticians do it discretely and continuously 5. We are right 95% of the time 6. We can legally comment on someone's posterior distribution 7. We may not be normal, but we are transformable 8. We never have to say we are certain 9. We are honestly significantly different 10. No one wants our jobs Visitate il portale http://www.modugno.it/ e in particolare la sezione su Palese http://www.modugno.it/archivio/palese/ __ R-help@stat.math.ethz.ch 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] Re: Run Sweave and LaTeX directly from command line
On Tue, 22 Feb 2005 20:21:46 -0600 Dirk Eddelbuettel wrote:
On 22 February 2005 at 14:58, [EMAIL PROTECTED] wrote:
| side note 2: For make afficionados the follwing 2 rules in
| combination with the Sweave script from the FAQ do almost the same
| (that's what I use :-)
Sure, but it requires a (arguably small) Makefile in every sweave
project directory. Would you consider integrating Gregor's script (or
a suitable modification) instead?
As discussed in various bits and pieces in this thread: the main
building blocks are all there in R so that you could do something like
Rnw2dvi - function(x, ...) {
Sweave(paste(x, .Rnw, sep = ))
texi2dvi(paste(x, .tex, sep = ), ...)
}
But as texi2dvi() requires some attention under Windows, people might
want to use something like Gabor suggested. Others will prefer a
Makefile like Fritz, others a shell script like Gregor or...
Personally, I've got a Rnw2pdf script with
---
Sweave $1
texi2dvi --pdf --clean $(basename $1 .Rnw).tex
rm -f $(basename $1 .Rnw)-*.pdf
rm -f Rplots.ps
---
in my bin/ directory.
So I think it is probably easiest if every user sets up something he/she
is comfortable with on his/her platform.
Z
Dirk
--
Better to have an approximate answer to the right question than a
precise answer to the wrong question. -- John Tukey as quoted by
John Chambers
__
R-help@stat.math.ethz.ch 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] Easy cut paste from Excel to R?
On Fri, 18 Feb 2005, Peter Dalgaard wrote: Don MacQueen [EMAIL PROTECTED] writes: I tried Ken's suggestion read.table(pipe(pbpaste),header=TRUE) on my Mac OS X system and it worked *without* generating any warning message. If my experience represents the norm, and Ken's is the exception, it is so simple that no further contribution to R is needed, I would say. Thank you, Ken. My conjecture is that it only happens when there are fewer than 5 data lines. We still need to sort out X11. Too bad that the xclip program isn't ubiquitous. The read side in X11 is not too hard, and R-devel now has read from the primary selection via file(clipboard). (I may change that name and allow reading from other selections later: I just made small changes to the Windows code.) Xlib doesn't it seems really have a clipboard, and so it is much harder to act as the provider of the primary selection (you need to respond to X11 events) -- xclip forks to do so. BTW, xclip is at http://people.debian.org/~kims/xclip/ and seems no longer under active development (last change 18 months ago). -- 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] using 'nice' with R
Roger D. Peng writes: On a Unix like system you can do `nice +19 R' or perhaps `nice +19 R CMD BATCH commands.R'. At least on Suse (9.1) and Debian (3.0) Linux, the syntax is `nice -19 R' (i.e. with `-', not `+'.) -- Bjørn-Helge Mevik __ R-help@stat.math.ethz.ch 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] Hosting a R Graph Gallery?
Earlier I suggested as a *start* for the R-gallery to process the examples in R base and in the packages on CRAN. However, this was only as a start It seems to me that there are three types of images that would be useful in such a gallery; 1: The processed examples of base and CRAN 2: Figures based upon code AND data sent in by users 3: Figures with code but WITHOUT data sent in by users (we could have figures without code or data [pure decorations] I would not support that myself) I think a package to automate production of gallery items would certainly be useful. I would suggest the basic function in the package would be one that produces the a Gallery Item with whatever details and in whatever format is suitable. Other functions may allow automatic submission of the Gallery Item and another would allow bulk processing of packages etc by looping over the function etc (Paul Murrell's graphicsQC package already processes all the example images in a package). I agree with Paul Murrell that storage in a database would be the best option though probably more difficult to set up at the start. I agree with Sander Om that encouraging SVG would be nice. I looked at the Japanese R-Wiki mentioned by Shigeru Mase and the RGallery Eric Lecoutre mentioned and certainly these are the *sort* of things we are interested in though perhaps with different formats and language ;-) I also agree with Gabor Grothendieck that it would be useful to make use of the existing R Wiki, though this depends on how much Gallery Item processing we would expect the host to do and whether this could be arranged on the Wiki. I guess how much the server would be expected to do is a critical question, from there we could proceed to the details of a Gallery Item format and then functions and a package. Cheers, Robert Cunningham __ R-help@stat.math.ethz.ch 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] Package pixmap breaks try() under circumstances
On Wed, 23 Feb 2005, Christian Lederer wrote: Dear R users, in some circumstances, try() shows a strange behaviour, when the pixmap package is loaded. The following piece of code works as expected, if it is either sourced in an interactive session or invoked via R CMD BATCH (the try-error is printed). However, if i invoke R using ``R --vanilla source.R'', the execution halts (without printing the try-error). # source.R library(pixmap) x - numeric() y - numeric() result - try(plot(lm(x~y))) print(result) If i don't use library(pixmap), than also R --vanilla source.R works as expected. This happens with R-2.0.1 and pixmap-0.4.2 under SuSE 9.2. But this does work as expected in your setting: library(pixmap) x - numeric() y - numeric() result - try(lm(x~y)) #result - try(plot(lm(x~y))) print(result) res1 - try(plot(result)) print(res1) The failed lm() in your version seems to be muddling plot(), and the try() is on the plot(), not the lm(). Without pixmap present: getAnywhere(plot) A single object matching 'plot' was found It was found in the following places package:graphics namespace:graphics with value ... With pixmap: getAnywhere(plot) 2 differing objects matching 'plot' were found in the following places package:pixmap package:graphics namespace:graphics Use [] to view one of them This also works: library(pixmap) x - numeric() y - numeric() result - try(plot(try(lm(x~y print(result) Hope this helps. Christian :-( P.S. The reason for preferring R --vanilla script.R over R CMD BATCH is, that i have to produce png images for a server and want to avoid the bitmap device for performance reasons. With R --vanilla script.R, this is possible, using the xvnc virtual X server. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Roger Bivand Economic Geography Section, Department of Economics, Norwegian School of Economics and Business Administration, Breiviksveien 40, N-5045 Bergen, Norway. voice: +47 55 95 93 55; fax +47 55 95 93 93 e-mail: [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] Run Sweave and LaTeX directly from command line
On Tue, 22 Feb 2005, Gorjanc Gregor wrote: Friedrich, This texi2dvi is also very nice. Dou you know of any pros and cons in comparison to rubber? texi2dvi is pretty standard, part of the texinfo distribution. rubber is not in the standard Linux distros, and uses python, also not very standard, especially on Windows. Is tex2dvi also shipped with R for windows? I am not able to find a binnary in my windows installation, however I found help page. Did I missed anything? I have 2.0.1 I believe Fritz was referring to the script src/script/texi2dvi. That is far too slow on Windows even if you have the tools installed as it runs ca 2000 processes per invocation. R has a function texi2dvi in package tools that avoids this on Windows. Let me repost (from a thread on R-devel about automating Sweave scripts on Windows). % cat Rnw.bat rterm --no-save --args %1 /tmp/MakeSweave.R %1.log % cat /tmp/MakeSweave.R library(tools) args - commandArgs() inp - args[length(args)] Sweave(inp) base - sub(\.(Rnw|Rtex)$, , inp) texi2dvi(paste(base, .tex, sep=), pdf=TRUE) shell.exec(paste(base, .pdf, sep=)) # display PDF file on Windows. which shows how easy this is to script in R in a cross-platform way On Tue, 22 Feb 2005 14:29:00 +0100, Gregor GORJANC (GG) wrote: Hello! Those of you, who use Sweave a lot, will probably find my shell script usable. You can get it at: http://www.bfro.uni-lj.si/MR/ggorjan/programs/shell/Sweave.sh No warranty, however don't hesitate to contact me if you find an error or have a patch! Very nice! Side note 1: R ships a version of texi2dvi, hence you might use that one in case rubber is not found. side note 2: For make afficionados the follwing 2 rules in combination with the Sweave script from the FAQ do almost the same (that's what I use :-) %.tex: %.Rnw Sweave $ %.pdf : %.tex texi2dvi --clean --pdf $ Best, Fritz Leisch -- --- Friedrich Leisch Institut für Statistik Tel: (+43 1) 58801 10715 Technische Universität WienFax: (+43 1) 58801 10798 Wiedner Hauptstraße 8-10/1071 A-1040 Wien, Austria http://www.ci.tuwien.ac.at/~leisch __ R-help@stat.math.ethz.ch 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] using 'nice' with R
On Wed, 23 Feb 2005, Bjørn-Helge Mevik wrote: Roger D. Peng writes: On a Unix like system you can do `nice +19 R' or perhaps `nice +19 R CMD BATCH commands.R'. At least on Suse (9.1) and Debian (3.0) Linux, the syntax is `nice -19 R' (i.e. with `-', not `+'.) The syntax depends on the shell you use: tcsh for example works as Roger says, and the `nice' executable as you say. -- 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] Re: Run Sweave and LaTeX directly from command line
On Tue, 22 Feb 2005 20:21:46 -0600,
Dirk Eddelbuettel (DE) wrote:
On 22 February 2005 at 14:58, [EMAIL PROTECTED] wrote:
| side note 2: For make afficionados the follwing 2 rules in combination
| with the Sweave script from the FAQ do almost the same (that's what I
| use :-)
Sure, but it requires a (arguably small) Makefile in every sweave project
directory.
Not really, the way I do it is that I have a $HOME/etc/Makeconf with a
lot of personal rules like those for .Snw files, and in my .bashrc I
define
fmake ()
{
make -f ~/etc/Makeconf $*
}
such that calling fmake instead of make gives me my private rules in
every place I need them :-)
Would you consider integrating Gregor's script (or a suitable
modification) instead?
Yes, sure, that's why I made him aware of texi2dvi: At least on Unix
systems we make sure that that is always available, while rubber might
not be. On the other hand, the R 2-liner I mentioned in an earlier
mail might be the better way to go. I always use the make rules I
described, hence had no need for either one of them - none is in the
R sources (yet).
Best,
Fritz
__
R-help@stat.math.ethz.ch 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] Run Sweave and LaTeX directly from command line
On Tue, 22 Feb 2005 20:46:01 +0100, Gorjanc Gregor (GG) wrote: Gabor, I definitely agree about use of R, however I am not so much in R as I am in Bash. It took me less than two hours to write this script. At home I use R under windows and I solve problem of portability with Cygwin. I know that this is not optimal solution for everyone but ... Can you provide me some insights/thoughts how this script might be written in R? I am open for discussion and cooperation on this script. What do you think Friedrich? I must look at this texi2dvi. I didn't hav any experience with it jet. I will look at it. Well, I didn't know about rubber before I read your script ... seems like the 2 do the same job. Some windows installations have a tool called texify which again does the same: run latex, bibtex, etc. as many times as needed. Ad shell vs. R: In R we already have almost all that is needed: example(Sweave) library(tools) texi2dvi(Sweave-test-1.tex) should do what you want ... and a Snw2dvi function is the corresponding 2-liner. If there is need for it I'm happy to include it in package utils. HTH, Fritz __ R-help@stat.math.ethz.ch 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] data.frame error message
At 09:47 PM 2/22/2005 -0500, Liaw, Andy wrote: Following the suggestions in the Posting Guide would help us to help you much better. What command(s) did you use to get the data into R? read.table(C:/.../persist.dat,col.names=c(a,b,etc.)) What does the `.dat' file look like? Itn was created in Word, as a plain text file. 6000 24 female 0.0014 1.4 1 3.47 0.25 2478.57 6000 168 female 0.0014 1.4 1 0.73 0.05 521.43 6000 96 female 0.0014 1.4 1 0.96 0.07 685.71 6000 168 female 0.0014 1.4 1 1.36 0.10 971.43 6000 24 female 0.0014 1.4 1 2.69 0.19 1921.43 6000 96 female 0.0014 1.4 1 0.76 0.05 542.86... What are `x' and `y'? 'x' is trt and 'y' is ppmtrans. I did subset commands to find all the males with ppmsamp above a certain value. persist[1:5,] dose trtsex massg massmg vol ppm perc ppmsamp 1 6000 24 female 0.00141.4 1 3.47 0.25 2478.57 2 6000 168 female 0.00141.4 1 0.73 0.05 521.43 3 6000 96female 0.00141.4 1 0.96 0.07 685.71 4 6000 168 female 0.00141.4 1 1.36 0.10 971.43 5 6000 24female 0.00141.4 12.69 0.19 1921.43 ...subsetting commands... ppmtrans-sqrt(sqrt(ppm.male.mark$ppmsamp)) persist.male-cbind(ppm.male.mark,ppmtrans) persist.male[1:5,] dose trt sex massg massmg vol ppm perc ppmsamp ppmtrans 27 6000 96 male 0.00121.2 1 1.82 0.15 1516.67 6.240549 55 6000 168 male 0.00121.2 1 0.90 0.08 750.00 5.233176 70 6000 168 male 0.00121.2 1 1.97 0.16 1641.67 6.365338 76 6000 96 male 0.00121.2 1 4.02 0.34 3350.00 7.607837 83 6000 168 male 0.00121.2 1 1.26 0.11 1050.00 5.692425 plot(trt,ppmtrans) Error in model.frame(formula, rownames, variables, varnames, extras, extranames, : variable lengths differ After my initial post, I found that if I specified the data set in the plot arguments, I could generate a boxplot. I then performed ANOVA on the response after creating factor levels for trt, and generated a stripchart showing the data in each of the three trts. But for future reference, if I need to do linear regression with a similar data set, I'd like to know what is happening. Thank you-- NH Ellis [[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] new package
Dear R-Users, We would like to announce a new package called frailtypack. This package fits a shared gamma frailty model and Cox proportional hazards model using a Penalized Likelihood on the hazard function. Left truncated, censored data and strata (max=2) are also allowed. Future versions will include functions for estimating nested frailty models using Penalized approach. Regards, Juan R Gonzalez Cancer Prevention and Control Unit Catalan Institute of Oncology, Barcelona (Spain) and Virginie RONDEAU, Ph.D. - CR2 INSERM Equipe Mixte INSERM E03 38 (Biostatistique), Université Victor Segalen Bordeaux2, 146 rue Léo Saignat, 33076 Bordeaux Cedex, France. email : [EMAIL PROTECTED] web site : http://www.isped.u-bordeaux2.fr/RECHERCHE/BIOSTATS/FR-BIOSTATS-Accueil.htm [[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] bias of a boot statistic
Question: How can I get access to the bias value of a boot statistic? Details: Boot function: boot(data, statistic, R, sim=ordinary, stype=i, strata=rep(1,n), L=NULL, m=0, weights=NULL, ran.gen=function(d, p) d, mle=NULL, ...) When I create an object, containing the bootstrap statistic (object - boot ())I can call it and will get an output with t, bias and standarderror as follows: Bootstrap Statistics : original biasstd. error t1* 5.650.01 0.9134185 My question is now, where is the value of the bias stored? How can I get access to this value to do further caluculations? Thanks, K. St. __ R-help@stat.math.ethz.ch 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] H-F corr.: covariance matrix for interaction effect
Hi,
I'm still not quite there with my H-F (G-G) correction code. I have it
working for the main effects, but I just can't figure out how to do it
for the effect interactions. The thing I really don't know (and can't
find anything about) is how to calculate the covariance matrix for the
interaction between the two (or even n) main factors.
I've looked through some books here and I've tried everything that came
to my mind, but I can't seem to be able to figure out an algorithm that
does it for me.
Could anyone give me a hint about how I could do this?
(I'll append my code at the end, in case that helps in any way...)
Thanks
Bela
# parameters for this function are:
# S - variance matrix (created by var() )
# k - number of factor levels (i.e. dim of S)
# n - number of measurements (i.e. number of rows in original matrix)
epsi.GG.HF - function (S,k,n) {
D - (k^2 * (mean(S) - mean(diag(S)))^2)
N1 - sum(S^2)
N2 - 2 * k * sum(apply(S, 1, mean)^2)
N3 - k^2 * mean(S)^2
epsiGG - D / ((k - 1) * (N1 - N2 + N3))
epsiHF - (n * (k-1) * epsiGG - 2) / ((k-1) * ((n-1) - (k-1)*epsiGG))
c(epsiGG,epsiHF)
}
# three factors, facROI,facCond,facSubj
# facROI,facCond are main effects, facSubj is repeatedness
# G-G and H-F corrections for a main effect
# we do the gghf stuff for the ROI, which means ROIs in columns,
# subjects in rows
mtx - NULL
for (iROI in 1:length(unique( facROI ))) {
for (iSubj in 1:length(unique( facSubj ))) {
mtx - c(mtx,
mean(vecData[facROI==unique(facROI)[iROI]
facSubj==unique(facSubj)[iSubj]])
)
}
}
mtx - matrix(mtx,ncol=length(unique( facROI )),byrow=F)
GgHfROI - epsi.GG.HF(var(mtx),length(mtx[1,]),length(mtx[,1]))
print(GgHfROI)
# now for the facROI:facCond interaction...how to go about this?
__
R-help@stat.math.ethz.ch 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] bias of a boot statistic
On Wed, 23 Feb 2005, K. Steinmann wrote: Question: How can I get access to the bias value of a boot statistic? Details: Boot function: boot(data, statistic, R, sim=ordinary, stype=i, strata=rep(1,n), L=NULL, m=0, weights=NULL, ran.gen=function(d, p) d, mle=NULL, ...) When I create an object, containing the bootstrap statistic (object - boot ())I can call it and will get an output with t, bias and standarderror as follows: Bootstrap Statistics : original biasstd. error t1* 5.650.01 0.9134185 My question is now, where is the value of the bias stored? How can I get access to this value to do further caluculations? apply(object, 2, mean, na.rm = TRUE) - boot.out$t0[1] in your case. I read that from getS3method(print, boot) -- 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] data.frame error message
From: Nic Ellis At 09:47 PM 2/22/2005 -0500, Liaw, Andy wrote: Following the suggestions in the Posting Guide would help us to help you much better. What command(s) did you use to get the data into R? read.table(C:/.../persist.dat,col.names=c(a,b,etc.)) What does the `.dat' file look like? Itn was created in Word, as a plain text file. 6000 24 female 0.0014 1.4 1 3.47 0.25 2478.57 6000 168 female 0.0014 1.4 1 0.73 0.05 521.43 6000 96 female 0.0014 1.4 1 0.96 0.07 685.71 6000 168 female 0.0014 1.4 1 1.36 0.10 971.43 6000 24 female 0.0014 1.4 1 2.69 0.19 1921.43 6000 96 female 0.0014 1.4 1 0.76 0.05 542.86... What are `x' and `y'? 'x' is trt and 'y' is ppmtrans. I did subset commands to find all the males with ppmsamp above a certain value. persist[1:5,] dose trtsex massg massmg vol ppm perc ppmsamp 1 6000 24 female 0.00141.4 1 3.47 0.25 2478.57 2 6000 168 female 0.00141.4 1 0.73 0.05 521.43 3 6000 96female 0.00141.4 1 0.96 0.07 685.71 4 6000 168 female 0.00141.4 1 1.36 0.10 971.43 5 6000 24female 0.00141.4 12.69 0.19 1921.43 ...subsetting commands... ppmtrans-sqrt(sqrt(ppm.male.mark$ppmsamp)) persist.male-cbind(ppm.male.mark,ppmtrans) persist.male[1:5,] dose trt sex massg massmg vol ppm perc ppmsamp ppmtrans 27 6000 96 male 0.00121.2 1 1.82 0.15 1516.67 6.240549 55 6000 168 male 0.00121.2 1 0.90 0.08 750.00 5.233176 70 6000 168 male 0.00121.2 1 1.97 0.16 1641.67 6.365338 76 6000 96 male 0.00121.2 1 4.02 0.34 3350.00 7.607837 83 6000 168 male 0.00121.2 1 1.26 0.11 1050.00 5.692425 plot(trt,ppmtrans) Error in model.frame(formula, rownames, variables, varnames, extras, extranames, : variable lengths differ After my initial post, I found that if I specified the data set in the plot arguments, I could generate a boxplot. I then performed ANOVA on the response after creating factor levels for trt, and generated a stripchart showing the data in each of the three trts. But for future reference, if I need to do linear regression with a similar data set, I'd like to know what is happening. Seems like you still haven't told all that you did. What is ppm.male.mark? Is it an exact copy of `persist', or some subset? Did you do any subsetting after you create persiste.male? If you get that error with plot(trt, ppmtrans), but not with plot(trt, ppmtrans, data=whatever), then the first version is probably getting `trt' and `ppmtrans' that you have lying around the global environment, which are probably results from different subsetting operations (and thus having different lengths). Without more info, we can only guess. One thing to keep in mind: When using the formula interface, it will save you some hair if all the variables used in the formula are in the data frame supplied (and try to always supply the data frame). Andy Thank you-- NH Ellis [[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] bias of a boot statistic
On Wed, 23 Feb 2005 11:42:05 +0100 K. Steinmann wrote: Question: How can I get access to the bias value of a boot statistic? Details: Boot function: boot(data, statistic, R, sim=ordinary, stype=i, strata=rep(1,n), L=NULL, m=0, weights=NULL, ran.gen=function(d, p) d, mle=NULL, ...) When I create an object, containing the bootstrap statistic (object - boot())I can call it and will get an output with t, bias and standarderror as follows: Bootstrap Statistics : original biasstd. error t1* 5.650.01 0.9134185 My question is now, where is the value of the bias stored? How can I get access to this value to do further caluculations? In the univariate case, it's simply mean(object$t - object$t0) hth, Z Thanks, K. St. __ R-help@stat.math.ethz.ch 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] H-F corr.: covariance matrix for interaction effect
Hi,
I'm still not quite there with my H-F (G-G) correction code. I have it
working for the main effects, but I just can't figure out how to do it
for the effect interactions. The thing I really don't know (and can't
find anything about) is how to calculate the covariance matrix for the
interaction between the two (or even n) main factors.
I've looked through some books here and I've tried everything that came
to my mind, but I can't seem to be able to figure out an algorithm that
does it for me.
Could anyone give me a hint about how I could do this?
(I'll append my code at the end, in case that helps in any way...)
Thanks
Bela
# parameters for this function are:
# S - variance matrix (created by var() )
# k - number of factor levels (i.e. dim of S)
# n - number of measurements (i.e. number of rows in original matrix)
epsi.GG.HF - function (S,k,n) {
D - (k^2 * (mean(S) - mean(diag(S)))^2)
N1 - sum(S^2)
N2 - 2 * k * sum(apply(S, 1, mean)^2)
N3 - k^2 * mean(S)^2
epsiGG - D / ((k - 1) * (N1 - N2 + N3))
epsiHF - (n * (k-1) * epsiGG - 2) / ((k-1) * ((n-1) - (k-1)*epsiGG))
c(epsiGG,epsiHF)
}
# three factors, facROI,facCond,facSubj
# facROI,facCond are main effects, facSubj is repeatedness
# G-G and H-F corrections for a main effect
# we do the gghf stuff for the ROI, which means ROIs in columns,
# subjects in rows
mtx - NULL
for (iROI in 1:length(unique( facROI ))) {
for (iSubj in 1:length(unique( facSubj ))) {
mtx - c(mtx,
mean(vecData[facROI==unique(facROI)[iROI]
facSubj==unique(facSubj)[iSubj]])
)
}
}
mtx - matrix(mtx,ncol=length(unique( facROI )),byrow=F)
GgHfROI - epsi.GG.HF(var(mtx),length(mtx[1,]),length(mtx[,1]))
print(GgHfROI)
# now for the facROI:facCond interaction...how to go about this?
__
R-help@stat.math.ethz.ch 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] Slightly off topic but concerning R#DSC-2005
On Wed, 23 Feb 2005 11:52:26 +0100, Ritter, Christian C GSMCIL-GSTMS/2 [EMAIL PROTECTED] wrote : Has anyone seen an official announcement to DSC-2005. I saw and email exchange on the R-news list in mid January in which a date was announced unofficially. (to DSC-organizers: we need and official announcement to request funding and travel permission and flights from Europe to Seattle are starting to fill up, so this is urgent). For the benefit of others, here is the posting giving the time and place. I don't think the official call for papers has come out yet. Duncan Murdoch From: Thomas Lumley tlumley_at_u.washington.edu Date: Thu 13 Jan 2005 - 02:23:50 EST On Wed, 12 Jan 2005, Marc Schwartz wrote to r-help: I have not seen anything posted yet for DSC 2005, unless I missed it someplace. DSC 2005 will be held in Seattle, at the University of Washington, August 13-15. This date is immediately after the Joint Statistical Meetings, and was chosen for the convenience of our European colleagues who might also be attending JSM. A call for papers and more information will be posted Real Soon Now. -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] Run Sweave and LaTeX directly from command line
Hi! I really started a nice debate. I have now removed rubber from my script. I found tex2dvi equally good and as Brian pointed it is more frequent than rubber. I also added some parts from Achim (rm ...). http://www.bfro.uni-lj.si/MR/ggorjan/programs/shell/Sweave.sh Shell vs R: It is nice to have such a script/function in R so it is easier to port it to Windows. However I find it crucial that user has ability to run this script from command line directly not only within R, which should not be the problem from what I have read from this mails. Maybe someone can put all this ideas in one pot. Friedrich? -- Lep pozdrav / With regards, Gregor GORJANC --- University of Ljubljana Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan Zootechnical Department mail: gregor.gorjanc at bfro.uni-lj.si Groblje 3 tel: +386 (0)1 72 17 861 SI-1230 Domzale fax: +386 (0)1 72 17 888 Slovenia, Europe __ R-help@stat.math.ethz.ch 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] data.frame error message
Hi Nic Be careful with variables with same names in your environment as variables in data.frame. If you have a variable with the same name e.g. ppmtrans in data frame and in environment, your commands use variable from environment. So probably trt is from data frame and ppmtrans from environment. See what ls() will show you. In models you can specify data=your.data.frame argument and all variables in formula are preferably used from this data frame, if I am not mistaken. Cheers Petr On 23 Feb 2005 at 5:42, Nic Ellis wrote: At 09:47 PM 2/22/2005 -0500, Liaw, Andy wrote: Following the suggestions in the Posting Guide would help us to help you much better. What command(s) did you use to get the data into R? read.table(C:/.../persist.dat,col.names=c(a,b,etc.)) What does the `.dat' file look like? Itn was created in Word, as a plain text file. 6000 24 female 0.0014 1.4 1 3.47 0.25 2478.57 6000 168 female 0.0014 1.4 1 0.73 0.05 521.43 6000 96 female 0.0014 1.4 1 0.96 0.07 685.71 6000 168 female 0.0014 1.4 1 1.36 0.10 971.43 6000 24 female 0.0014 1.4 1 2.69 0.19 1921.43 6000 96 female 0.0014 1.4 1 0.76 0.05 542.86... What are `x' and `y'? 'x' is trt and 'y' is ppmtrans. I did subset commands to find all the males with ppmsamp above a certain value. persist[1:5,] dose trtsex massg massmg vol ppm perc ppmsamp 1 6000 24 female 0.00141.4 1 3.47 0.25 2478.57 2 6000 168 female 0.00141.4 1 0.73 0.05 521.43 3 6000 96female 0.00141.4 1 0.96 0.07 685.71 4 6000 168 female 0.00141.4 1 1.36 0.10 971.43 5 6000 24female 0.00141.4 12.69 0.19 1921.43 ...subsetting commands... ppmtrans-sqrt(sqrt(ppm.male.mark$ppmsamp)) persist.male-cbind(ppm.male.mark,ppmtrans) persist.male[1:5,] dose trt sex massg massmg vol ppm perc ppmsamp ppmtrans 27 6000 96 male 0.00121.2 1 1.82 0.15 1516.67 6.240549 55 6000 168 male 0.00121.2 1 0.90 0.08 750.00 5.233176 70 6000 168 male 0.00121.2 1 1.97 0.16 1641.67 6.365338 76 6000 96 male 0.00121.2 1 4.02 0.34 3350.00 7.607837 83 6000 168 male 0.00121.2 1 1.26 0.11 1050.00 5.692425 plot(trt,ppmtrans) Error in model.frame(formula, rownames, variables, varnames, extras, extranames, : variable lengths differ After my initial post, I found that if I specified the data set in the plot arguments, I could generate a boxplot. I then performed ANOVA on the response after creating factor levels for trt, and generated a stripchart showing the data in each of the three trts. But for future reference, if I need to do linear regression with a similar data set, I'd like to know what is happening. Thank you-- NH Ellis [[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 Petr Pikal [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] large data set, and RDBMS
I have this somewhat large data set that is given to me in a fixed width format. The file itself is already 100MB (Maybe R can actually handle this but I am trying to gain some experience in postgres and RODBC). I am using postgres to preprocess the file and connect to the database through the RODBC package. My question is much before the processing of the database in R (Iknow this is somewhat off topic). I have more than 6000 variables in the data and postgres is outputing and error that it cannot handle more than 1600 columns in one table. Is there any Linux database managers that actually do do this.; Of course I can break it into 4 or 5 tables but I am wondering if there is a better solution. System is a Linux 2.6.8 psql (PostgreSQL) 7.4.7 Thank you, Jean __ R-help@stat.math.ethz.ch 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] H-F corr.: covariance matrix for interaction effect
Bela Bauer [EMAIL PROTECTED] writes: Hi, I'm still not quite there with my H-F (G-G) correction code. I have it working for the main effects, but I just can't figure out how to do it for the effect interactions. The thing I really don't know (and can't find anything about) is how to calculate the covariance matrix for the interaction between the two (or even n) main factors. I've looked through some books here and I've tried everything that came to my mind, but I can't seem to be able to figure out an algorithm that does it for me. Could anyone give me a hint about how I could do this? (I'll append my code at the end, in case that helps in any way...) I have given it to you before: My plan is to drop the explicit formula involving on/off diagonal elements of S and go directly at Box (1954), theorems 3.1 and 6.1, involving eigenvalues of TST', where T is the relevant residual operator. In the case where one of the factors have only two levels, I believe you just take differences and use the usual formula, but more than two levels is tricky. -- O__ Peter Dalgaard Blegdamsvej 3 c/ /'_ --- Dept. of Biostatistics 2200 Cph. N (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] r: ridge regression
hello all some help required once again! does anyone recall the equations for the following ridge constants? 1. hoerl and kennard (1970) 2. hoerl, kennard and baldwin (1975) 3. lawless and wang could you also specify whether or not one has to transform the X and Y variables. if so , how and in which cases. a worked example with a data set would be most helpful. thanking you in advance *** Allan__ R-help@stat.math.ethz.ch 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] filling columns in frame according to another column frame
R-help, I have a frame which I want to fill up conditioning to another data frame column. The one I want to fill up is as follows (basically an empty one): test2 cm 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 The other looks like : test1 cm 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 38 0 0 1 0 0 0 6 0 0 0 0 0 0 0 0 0 39 0 0 1 0 0 0 0 0 0 0 0 0 6 0 0 0 40 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 41 0 0 2 0 0 0 0 0 0 0 6 0 0 0 0 0 43 0 0 1 0 0 0 4 0 0 0 0 0 0 0 0 0 44 0 0 4 0 0 0 5 0 0 0 0 0 0 0 0 0 45 0 0 2 0 0 0 0 0 0 0 6 0 0 0 0 0 47 0 0 3 0 0 0 0 0 0 0 0 0 0 0 0 0 48 0 0 2 0 0 0 0 0 0 6 0 0 0 0 0 0 49 0 0 2 0 0 0 0 0 0 6 0 0 0 0 0 0 50 0 0 3 0 0 0 0 0 0 3 0 0 0 0 0 0 51 0 0 2 0 0 0 0 0 0 3 0 0 0 0 0 0 Length of both frames are different ( test2 = 150 and test1 = 70 ) The key column is 'cm' I have tried someting (fill column '3' in test2): test2 [, '3' ]- ifelse ( test2$'cm' %in% test1$'cm' , test1$'3' , 0) but the result is wrong. Any suggestions? Thanks in advance version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major2 minor0.1 year 2004 month11 day 15 language R __ R-help@stat.math.ethz.ch 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] to print dataframe
Dear all, Is it possible to print a dataframe without the row numbers? For example if I have a dataframe like that: df - data.frame(name1=sample(LETTERS,10),name2=sample(c(0,1),10,replace=TRUE)) after printing name1 name2 1 O 1 2 H 0 3 R 0 4 T 0 5 V 1 6 E 0 7 W0 8 P 1 9 G 0 10 J 1 But I would like the dataframe printed like that name1 name2 O 1 H 0 R 0 T 0 V 1 E 0 W 0 P 1 G 0 J 1 I look at ?print.dataframe ?data.frame but I can't find anything. Somebody could help me? Thanks in advance. __ Tiscali Adsl 3 Mega Flat, 3 MESI GRATIS! Con Tiscali Adsl 3 Mega Flat navighi in Rete alla supervelocita' a soli 29.95 euro al mese senza limiti di tempo. Attivati entro il 28 Febbraio 2005, 3 MESI sono GRATIS Scopri come http://abbonati.tiscali.it/adsl/sa/2flat_tc/ __ R-help@stat.math.ethz.ch 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] model.matrix for a factor effect with no intercept
I was surprised by this (in R 2.0.1): a - ordered(-1:1) a [1] -1 0 1 Levels: -1 0 1 model.matrix(~ a) (Intercept) a.La.Q 1 1 -7.071068e-01 0.4082483 2 1 -9.073800e-17 -0.8164966 3 1 7.071068e-01 0.4082483 attr(,assign) [1] 0 1 1 attr(,contrasts) attr(,contrasts)$a [1] contr.poly model.matrix(~ -1 + a) a-1 a0 a1 1 1 0 0 2 0 1 0 3 0 0 1 attr(,assign) [1] 1 1 1 attr(,contrasts) attr(,contrasts)$a [1] contr.poly Without the intercept, treatment contrasts seem to have been used (this despite the contr.poly in the contrasts attribute). It's not restricted to ordered factors. For example, if Helmert contrasts are used for nominal factors, the same sort of thing happens. I suppose it is a deliberate feature (perhaps to protect the user from accidentally fitting models that make no sense? or maybe some better reason?) -- is it explained somewhere? David __ R-help@stat.math.ethz.ch 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] filling columns in frame according to another column frame
I am confused. Are you saying that your two data frames are of different dimensions ? In any case what I think what you are looking for is which. # generate the conditioning matrix a - matrix( sample(0:1, 9, replace=TRUE), nc=3 ) a [,1] [,2] [,3] [1,]111 [2,]100 [3,]001 # find the index where zero is present ( w - which( a == 0, arr.ind=T ) ) row col [1,] 3 1 [2,] 2 2 [3,] 3 2 [4,] 2 3 # generate the matrix of interest ( b - matrix(1:9, nc=3, byrow=T) ) [,1] [,2] [,3] [1,]123 [2,]456 [3,]789 # values that will be used to impute the zero's in a b[w] [1] 7 5 8 6 # impute the values of a with b where a is zero a[w] - b[w] # the result a [,1] [,2] [,3] [1,]111 [2,]156 [3,]781 Regards, Adai On Wed, 2005-02-23 at 14:14 +, Luis Ridao Cruz wrote: R-help, I have a frame which I want to fill up conditioning to another data frame column. The one I want to fill up is as follows (basically an empty one): test2 cm 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 The other looks like : test1 cm 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 38 0 0 1 0 0 0 6 0 0 0 0 0 0 0 0 0 39 0 0 1 0 0 0 0 0 0 0 0 0 6 0 0 0 40 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 41 0 0 2 0 0 0 0 0 0 0 6 0 0 0 0 0 43 0 0 1 0 0 0 4 0 0 0 0 0 0 0 0 0 44 0 0 4 0 0 0 5 0 0 0 0 0 0 0 0 0 45 0 0 2 0 0 0 0 0 0 0 6 0 0 0 0 0 47 0 0 3 0 0 0 0 0 0 0 0 0 0 0 0 0 48 0 0 2 0 0 0 0 0 0 6 0 0 0 0 0 0 49 0 0 2 0 0 0 0 0 0 6 0 0 0 0 0 0 50 0 0 3 0 0 0 0 0 0 3 0 0 0 0 0 0 51 0 0 2 0 0 0 0 0 0 3 0 0 0 0 0 0 Length of both frames are different ( test2 = 150 and test1 = 70 ) The key column is 'cm' I have tried someting (fill column '3' in test2): test2 [, '3' ]- ifelse ( test2$'cm' %in% test1$'cm' , test1$'3' , 0) but the result is wrong. Any suggestions? Thanks in advance version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major2 minor0.1 year 2004 month11 day 15 language R __ R-help@stat.math.ethz.ch 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] Memory error in Mac OS X Aqua GUI v1.01 with cluster pack age functions
It's trying to allocate about 850 MB. And that's just the object that broke the camel's back. You probably really are out of memory. You could increase swap space and cross your fingers, but probably daisy creates the 10481 x 10481 distance matrix, which would be about 800 MB since each entry is 8 bytes. It may even create multiple copies. You might try increasing RAM to 4 GB or a clustering routine that doesn't need a distance matrix, like k-means (which is based on Euclinean distance as well). Reid Huntsinger -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Betty Gilbert Sent: Tuesday, February 22, 2005 5:57 PM To: r-help@stat.math.ethz.ch Subject: [R] Memory error in Mac OS X Aqua GUI v1.01 with cluster package functions I'm sorry if the answer to my problem is buried in the archives. I have limited experience with R and I couldn't find a solution to my particular problem. I am running Mac OS X Aqua GUI v1.01 on a new G5 running os 10.3.8 with a 1.8Ghz processor and 1GB of sdram. I just downloaded bioconducter a week ago and I'm trying to cluster a matrix I created with a simulation with dimensions dim(nca35) [1] 1048112 with size object.size(nca352) [1] 1426204 I checked my ulimits variable on the unix terminal and it says it's unlimited as does mem.limits() nsize vsize NANA But I'm still getting errors like the following with funtions in the cluster package daisy(nca352, metric= euclidean, stand=FALSE)-dnca35 Error: cannot allocate vector of size 858213 Kb *** malloc: vm_allocate(size=878813184) failed (error code=3) *** malloc[599]: error: Can't allocate region if it helps i also checked gc() used (Mb) gc trigger (Mb) Ncells 448662 12.0 741108 19.8 Vcells 847630 6.5 135357901 1032.7 I tried the suggested unix command in the memory help doc but that doesn't work in the Aqua GUI. Can someone tell me how to change the Vcells? Although to the best of my understanding (which is limited) I shouldn't have to do that. Any suggestions would be greatly appreciated. thanks, betty -- Betty Gilbert [EMAIL PROTECTED] Taylor Lab Plant and Microbial Biology 321 Koshland Hall U.C. Berkeley Berkeley, Ca 94720 __ R-help@stat.math.ethz.ch 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] Need your help in calculating the p-value
Hi, I am using R to perform wilcox.test and wondering if you know how the p-value in wilcox.test is calculated? Thank you for your help:-) Regards, Latha __ R-help@stat.math.ethz.ch 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] filling columns in frame according to another column frame
I think I did not explain very well what my problem is, Both frames have the same number of columns but different number of rows. The point is to compare 'cm' column in test1 with 'cm' in test2 as follows first element in test2$'cm' versus first element in test1$'cm' first element in test2$'cm' versus second element in test1$'cm' first element in test2$'cm' versus element element in test1$'cm' . if any of the above matches returns the value in column , lets say 10, of element in test1$'cm' , if not then 0 Luis Adaikalavan Ramasamy [EMAIL PROTECTED] 23/02/2005 15:07:55 I am confused. Are you saying that your two data frames are of different dimensions ? In any case what I think what you are looking for is which. # generate the conditioning matrix a - matrix( sample(0:1, 9, replace=TRUE), nc=3 ) a [,1] [,2] [,3] [1,]111 [2,]100 [3,]001 # find the index where zero is present ( w - which( a == 0, arr.ind=T ) ) row col [1,] 3 1 [2,] 2 2 [3,] 3 2 [4,] 2 3 # generate the matrix of interest ( b - matrix(1:9, nc=3, byrow=T) ) [,1] [,2] [,3] [1,]123 [2,]456 [3,]789 # values that will be used to impute the zero's in a b[w] [1] 7 5 8 6 # impute the values of a with b where a is zero a[w] - b[w] # the result a [,1] [,2] [,3] [1,]111 [2,]156 [3,]781 Regards, Adai On Wed, 2005-02-23 at 14:14 +, Luis Ridao Cruz wrote: R-help, I have a frame which I want to fill up conditioning to another data frame column. The one I want to fill up is as follows (basically an empty one): test2 cm 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 The other looks like : test1 cm 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 38 0 0 1 0 0 0 6 0 0 0 0 0 0 0 0 0 39 0 0 1 0 0 0 0 0 0 0 0 0 6 0 0 0 40 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 41 0 0 2 0 0 0 0 0 0 0 6 0 0 0 0 0 43 0 0 1 0 0 0 4 0 0 0 0 0 0 0 0 0 44 0 0 4 0 0 0 5 0 0 0 0 0 0 0 0 0 45 0 0 2 0 0 0 0 0 0 0 6 0 0 0 0 0 47 0 0 3 0 0 0 0 0 0 0 0 0 0 0 0 0 48 0 0 2 0 0 0 0 0 0 6 0 0 0 0 0 0 49 0 0 2 0 0 0 0 0 0 6 0 0 0 0 0 0 50 0 0 3 0 0 0 0 0 0 3 0 0 0 0 0 0 51 0 0 2 0 0 0 0 0 0 3 0 0 0 0 0 0 Length of both frames are different ( test2 = 150 and test1 = 70 ) The key column is 'cm' I have tried someting (fill column '3' in test2): test2 [, '3' ]- ifelse ( test2$'cm' %in% test1$'cm' , test1$'3' , 0) but the result is wrong. Any suggestions? Thanks in advance version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major2 minor0.1 year 2004 month11 day 15 language R __ R-help@stat.math.ethz.ch 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] RE: R-help Digest, Vol 24, Issue 22
I used SAS to analyze the data initially, since the data set was made up of several files when I received it, and I'm still not very good at manipulating data in R. I have posted the data set from one location at the following address: http://uwstudentweb.uwyo.edu/A/AKNISS/sxherb.txt var=cultivar trt=herbicide treatment yield=response variable of interest All plot# from 101 to 104 are rep 1, 201-204 rep 2, and 301 to 304 rep 3. It was the only file that was in an easy format for R to read at the moment, and was probably the most reliable trial of the two locations. I would like to use power.anova.test() with this data set to plan next years study (to get a sample size for each herb*var combination), but I'm not quite sure how that is done for an interaction effect. Do I just use the MS for herb*var as the between group variance and the MSE as the within group variance? Or do I need to somehow include other variance parameters in the model? The model for this location (split-block design): yield = rep + herb + var + herb*var ## all are fixed effects rep*herb = error term for herb rep*var = error term for cultivar residual = error term for herb*var I hope this attempt at my question was a little more clear. I appreciate any help that is offered. Andrew Kniss Assistant Research Scientist University of Wyoming Department of Plant Sciences 1000 E. Univesity Ave Laramie, WY 82071 USA [EMAIL PROTECTED] -Original Message- From: John Maindonald [mailto:[EMAIL PROTECTED] Sent: Tuesday, February 22, 2005 3:37 PM To: r-help@stat.math.ethz.ch Cc: [EMAIL PROTECTED] Subject: Re: R-help Digest, Vol 24, Issue 22 You need to give the model formula that gave your output. There are two sources of variation (at least), within and between locations; though it looks as though your analysis may have tried to account for this (but if so, the terms are not laid out in a way that makes for ready interpretation. The design is such (two locations) that you do not have much of a check that effects are consistent over locations. You need to check whether results really are similar for all cultivars and for all herbicides, so that it is legitimate to pool as happens in the overall analysis. If a herbicide:cultivar combination has little effect the variability may be large, while if it has a dramatic effect (kills everything!), there may be no variability to speak of. John Maindonald. On 22 Feb 2005, at 10:06 PM, [EMAIL PROTECTED] wrote: To: 'Bob Wheeler' [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Subject: RE: [R] power.anova.test for interaction effects Reply-To: [EMAIL PROTECTED] It's a rather complex model. A 37*4 factorial (37 cultivars[var]; 4 herbicide treatments[trt]) with three replications[rep] was carried out at two locations[loc], with different randomizations within each rep at each location. Source DF Error Term MS Loc 1 Trt*rep(loc)12314 Rep(loc) 4 Trt*rep(loc)1230.5 Trt 3 Trt*rep(loc)64.72 Trt*loc 3 Trt*rep(loc)33.42 Trt*rep(loc) 12 Residual76.78 Var 36 Var*trt*loc 93.91 Var*trt 108 Var*trt*loc 12.06 Var*trt*loc 144 Residual43.09 Residual575 NA 21.23 -Original Message- From: Bob Wheeler [mailto:[EMAIL PROTECTED] Sent: Monday, February 21, 2005 4:33 PM To: [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Subject: Re: [R] power.anova.test for interaction effects Your F value is so low as to make me suspect your model. Where did the 144 denominator degrees of freedom come from? John Maindonald email: [EMAIL PROTECTED] phone : +61 2 (6125)3473fax : +61 2(6125)5549 Centre for Bioinformation Science, Room 1194, John Dedman Mathematical Sciences Building (Building 27) Australian National University, Canberra ACT 0200. __ R-help@stat.math.ethz.ch 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 Building R on AIX 5.2.0.0 (Re-post)
Hello, I am trying to build R 2.0.1 on an AIX 5.2.0.0 machine using gcc 3.3.2: $ oslevel 5.2.0.0 $ gcc -v Reading specs from /usr/local/lib/gcc-lib/powerpc-ibm-aix5.2.0.0/3.3.2/specs Configured with: ../gcc-3.3.2/configure : (reconfigured) ../gcc-3.3.2/configure --disable-nls : (reconfigured) ../gcc-3.3.2/configure --disable-nls Thread model: aix gcc version 3.3.2 Configure goes okay, but I get an error that I don't quite know how to interpret during make. I've included the summary output from the end of configure as well as the error that I get during make below. Any suggestions/recommendations are very much appreciate: I'm stuck on ideas for what could be going wrong. Paul $ ./configure --prefix=/db2blaste/R snip R is now configured for powerpc-ibm-aix5.2.0.0 Source directory: . Installation directory:/db2blast/R C compiler:gcc -mno-fp-in-toc -g -O2 C++ compiler: g++ -g -O2 Fortran compiler: g77 -g -O2 Interfaces supported: X11 External libraries: Additional capabilities: PNG, JPEG Options enabled: R profiling Recommended packages: yes configure: WARNING: you cannot build DVI versions of the R manuals configure: WARNING: you cannot build info or html versions of the R manuals configure: WARNING: you cannot build PDF versions of the R manuals configure: WARNING: I could not determine a browser configure: WARNING: I could not determine a PDF viewer $ make snip gcc -Wl,-bM:SRE -Wl,-H512 -Wl,-T512 -Wl,-bnoentry -Wl,-bexpall -Wl,- bI:. ./../../etc/R.exp -L/usr/local/lib -o lapack.so -Wl,-bI:../../../etc/Rlapack.exp Lapack.lo rgeev.lo rsyev.lo -L../../../lib -lRlapack -L/usr/local/lib -L/usr/ local/lib/gcc-lib/powerpc-ibm-aix5.2.0.0/3.3.2 -L/usr/local/lib/gcc-lib/powe rpc- ibm-aix5.2.0.0/3.3.2/../../.. -lfrtbegin -lg2c -lm -lgcc_s /usr/local/lib/gcc- lib/powerpc-ibm-aix5.2.0.0/3.3.2/libgcc.a -lg -ldl -ltermcap -lm -lc ld: 0706-006 Cannot find or open library file: -l Rlapack ld:open(): A file or directory in the path name does not exist. collect2: ld returned 255 exit status make: 1254-004 The error code from the last command is 1. Stop. make: 1254-004 The error code from the last command is 2. Stop. make: 1254-004 The error code from the last command is 1. Stop. make: 1254-004 The error code from the last command is 1. Stop. make: 1254-004 The error code from the last command is 1. Stop. __ R-help@stat.math.ethz.ch 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] model.matrix for a factor effect with no intercept
MASS4 p.150 White Book p.38 Those are the only two reasonably comprehensive accounts that I am aware of (and they have only partial overlap). The underlying motivation is to span the _additional_ vector space covered by the term, the complement to what has gone before. Put another way, as each term is added, only enough columns are added to the model matrix to span the same space as if dummy coding had been used for that term and its predecessors. So think of this as a way to produce a parsimonious (usually full-rank) basis for the model space. On Wed, 23 Feb 2005, David Firth wrote: I was surprised by this (in R 2.0.1): a - ordered(-1:1) a [1] -1 0 1 Levels: -1 0 1 model.matrix(~ a) (Intercept) a.La.Q 1 1 -7.071068e-01 0.4082483 2 1 -9.073800e-17 -0.8164966 3 1 7.071068e-01 0.4082483 attr(,assign) [1] 0 1 1 attr(,contrasts) attr(,contrasts)$a [1] contr.poly model.matrix(~ -1 + a) a-1 a0 a1 1 1 0 0 2 0 1 0 3 0 0 1 attr(,assign) [1] 1 1 1 attr(,contrasts) attr(,contrasts)$a [1] contr.poly Without the intercept, treatment contrasts seem to have been used (this despite the contr.poly in the contrasts attribute). It's not restricted to ordered factors. For example, if Helmert contrasts are used for nominal factors, the same sort of thing happens. I suppose it is a deliberate feature (perhaps to protect the user from accidentally fitting models that make no sense? or maybe some better reason?) -- is it explained somewhere? -- 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
[R] Sweave and \input or \include LaTeX commands
Hello! I was just wondering if Sweave can work with \input or \include LaTeX commands. So, is it aware of such a possible hierarchy in documents. I would test that, but I don't have such a report available at the moment. I thought of that when I was writting shell script for Sweave from command line and I have solved that part there. -- Lep pozdrav / With regards, Gregor GORJANC --- University of Ljubljana Biotechnical Faculty URI: http://www.bfro.uni-lj.si Zootechnical Departmentemail: gregor.gorjanc at bfro.uni-lj.si Groblje 3 tel: +386 (0)1 72 17 861 SI-1230 Domzalefax: +386 (0)1 72 17 888 Slovenia __ R-help@stat.math.ethz.ch 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] Need your help in calculating the p-value
Latha Raja wrote: Hi, I am using R to perform wilcox.test and wondering if you know how the p-value in wilcox.test is calculated? You can get the source code of the wilcox.test function: methods(wilcox.test) [1] wilcox.test.default* wilcox.test.formula* Non-visible functions are asterisked getAnywhere(wilcox.test.default) ... You will see that it uses 'psignrank' (for the one sample case). Christophe Pallier www.pallier.org __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] filling columns in frame according to another column frame
So you want to grep for each pattern as indicated by columns of test2 in the columns of test1. Something like your initial approach. p - c(1,3,5) # pattern like test2$'cm' x - sample(1:10) # datalike test1$'cm' x [1] 3 4 9 8 7 5 10 2 1 6 x %in% p [1] TRUE FALSE FALSE FALSE FALSE TRUE FALSE FALSE TRUE FALSE (w - which( x %in% p )) [1] 1 6 9 If no matches are found, you get numeric(0) for 'w' and I usually test this using length(w)==0. On second thought, I tried your original solution and it appears to work but not tested thoroughly. Can you provide a _simple_ example where this does not work. Thank you. z - 1:length(x) # like test1$'3' z [1] 1 2 3 4 5 6 7 8 9 10 ifelse( x %in% p, z, 0 ) [1] 1 0 0 0 0 6 0 0 9 0 BTW, why are you quoting your names as in 'cm' and '3' and are test1, test2 really data.frames ? What does class(test1) and class(test1$'cm') say ? Regards, Adai On Wed, 2005-02-23 at 15:38 +, Luis Ridao Cruz wrote: I think I did not explain very well what my problem is, Both frames have the same number of columns but different number of rows. The point is to compare 'cm' column in test1 with 'cm' in test2 as follows first element in test2$'cm' versus first element in test1$'cm' first element in test2$'cm' versus second element in test1$'cm' first element in test2$'cm' versus element element in test1$'cm' . if any of the above matches returns the value in column , lets say 10, of element in test1$'cm' , if not then 0 Luis Adaikalavan Ramasamy [EMAIL PROTECTED] 23/02/2005 15:07:55 I am confused. Are you saying that your two data frames are of different dimensions ? In any case what I think what you are looking for is which. # generate the conditioning matrix a - matrix( sample(0:1, 9, replace=TRUE), nc=3 ) a [,1] [,2] [,3] [1,]111 [2,]100 [3,]001 # find the index where zero is present ( w - which( a == 0, arr.ind=T ) ) row col [1,] 3 1 [2,] 2 2 [3,] 3 2 [4,] 2 3 # generate the matrix of interest ( b - matrix(1:9, nc=3, byrow=T) ) [,1] [,2] [,3] [1,]123 [2,]456 [3,]789 # values that will be used to impute the zero's in a b[w] [1] 7 5 8 6 # impute the values of a with b where a is zero a[w] - b[w] # the result a [,1] [,2] [,3] [1,]111 [2,]156 [3,]781 Regards, Adai On Wed,2005-02-23 at 14:14 +, Luis Ridao Cruz wrote: R-help, I have a frame which I want to fill up conditioning to another data frame column. The one I want to fill up is as follows (basically an empty one): test2 cm 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 The other looks like : test1 cm 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 38 0 0 1 0 0 0 6 0 0 0 0 0 0 0 0 0 39 0 0 1 0 0 0 0 0 0 0 0 0 6 0 0 0 40 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 41 0 0 2 0 0 0 0 0 0 0 6 0 0 0 0 0 43 0 0 1 0 0 0 4 0 0 0 0 0 0 0 0 0 44 0 0 4 0 0 0 5 0 0 0 0 0 0 0 0 0 45 0 0 2 0 0 0 0 0 0 0 6 0 0 0 0 0 47 0 0 3 0 0 0 0 0 0 0 0 0 0 0 0 0 48 0 0 2 0 0 0 0 0 0 6 0 0 0 0 0 0 49 0 0 2 0 0 0 0 0 0 6 0 0 0 0 0 0 50 0 0 3 0 0 0 0 0 0 3 0 0 0 0 0 0 51 0 0 2 0 0 0 0 0 0 3 0 0 0 0 0 0 Length of both frames are different ( test2 = 150 and test1 = 70 ) The key column is 'cm' I have tried someting (fill column '3' in test2): test2 [, '3' ]- ifelse ( test2$'cm' %in% test1$'cm' , test1$'3' , 0) but the result is wrong. Any suggestions? Thanks in advance version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major2 minor0.1 year 2004 month11 day 15 language R __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide!
Re: [R] Need your help in calculating the p-value
See the details and references section of help(wilcox.test). A small section from the details might be relevant : By default (if 'exact' is not specified), an exact p-value is computed if the samples contain less than 50 finite values and there are no ties. Otherwise, a normal approximation is used. Also type in stats:::wilcox.test.default or getAnywhere (wilcox.test.default) to see the actual codes. Regards, Adai On Wed, 2005-02-23 at 10:28 -0500, Latha Raja wrote: Hi, I am using R to perform wilcox.test and wondering if you know how the p-value in wilcox.test is calculated? Thank you for your help:-) Regards, Latha __ R-help@stat.math.ethz.ch 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] RODBC type conversion bug
I run R 2.0.1 on Debian and connect to Informix database via RODBC. In the table below the column month is of type char(1). RODBC seems to be converting this column to boolean if the value is F or T. This is the data in my table: id month year 25 F 2005 26 Z 2005 When I select * for id 25 I get 25 FALSE 2005 When I select * for id 16 I get 26Z 2005 Is there a fix for this issue? I invoke odbc connect with: db - odbcConnect(mydb, believeNRows=FALSE) __ R-help@stat.math.ethz.ch 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] Graphics (crashes under Windows)
In message [EMAIL PROTECTED], r-help- [EMAIL PROTECTED] writes The R platform that I installed on my Windows XP crashes everytime that I try to run some sophisticated graphics (e.g. Demo Graphics). Is that to do with the configuration? Shall I reinstall it? You may have a buggy video driver. If you go to Control Panel, Display, Settings, Advanced, Troubleshoot, and reduce the hardware acceleration, it may fix the problem. (Maybe it is worth adding this trick to the R for Windows FAQ?) -- Graham Jones, author of SharpEye Music Reader http://www.visiv.co.uk 21e Balnakeil, Durness, Lairg, Sutherland, IV27 4PT, Scotland, UK __ R-help@stat.math.ethz.ch 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] stopping a function
I've looked for this information in all the R help sources I could find and found nothing. Is it possible to use some function key to stop the execution of some R command without ending the R session (Windows, R 1.91)? I've several times started functions that for various reasons are not executing properly and it would be nice to stop them without killing the R session. I've been using taskmgr to end the R session but I then lose all objects created during the session. Brian Brian S. Cade U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: [EMAIL PROTECTED] tel: 970 226-9326 __ R-help@stat.math.ethz.ch 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] problems with nonlinear fits using nls
Corey Bradshaw wrote:
Hello colleagues,
I am attempting to determine the nonlinear least-squares estimates of
the nonlinear model parameters using nls. I have come across a common
problem that R users have reported when I attempt to fit a particular
3-parameter nonlinear function to my dataset:
Error in nls(r ~ tlm(a, N.fix, k, theta), data = tlm.data, start =
list(a = a.st, :
step factor 0.000488281 reduced below `minFactor' of 0.000976563
Despite modifying minFactor using nls.control, I am unable to counter
the apparent singularity in the model fit. I have also tried changing
the tolerance and start parameter values to no avail. If anyone can
provide a relatively simple solution (perhaps adjusting the gradient,
but I'm not sure how to do this), I would be most appreciative. My
dataset is:
tlm.data
r N.fix
1 -0.5240708576
2 0.1053605245
3 -0.1743533950
4 0.1941560142
5 0.4870149851
6 -0.5068176083
7 -0.1743533950
8 0.5527898242
9 0.1521918273
10 0.4989911785
11 0.10821358 140
12 -0.83034830 156
13 -0.3074847068
14 -0.2231435550
15 0.0487901640
16 -0.0487901642
17 0.7537718040
18 -0.1251631485
19 -0.3662443975
My function is:
tlm - function(a,N,k,theta) (a*(1-((N/k)^theta)))
The nls fit I've coded is:
tlm.fit - try(nls(r~tlm(a,N.fix,k,theta), data=tlm.data,
start=list(a=a.st,k=k.st,theta=1),
trace=TRUE,
control=nls.control(maxiter=6000,tol=1e-05,minFactor=1/1024)))
I'm using start values parsed in from another (previous, but not shown)
model fit. In this case,
a.st
[1] 0.3812922
k.st
[1] 64.66529
I happen to know the true values for the optimised parameters (from
another application), but I can't get nls to reproduce them. They are:
a = 2.0466
k = 60.8275
theta = 0.2277
Any ideas?
Regards,
Corey Bradshaw
It is possible to fit this model to these data using nls as shown in the
enclosed transcript. You have one conditionally linear parameter ('a')
in the model so I used the plinear algorithm and I also generated
analytic derivatives for the tls function using the deriv function.
There are several things to note:
- Your data are very noisy. It is not surprising that it is difficult
to fit a 3-parameter nonlinear model to such data.
- The fitted model has negative values for both 'a' and 'theta'.
- The estimates are highly imprecise.
summary(fm1)
Formula: r ~ tlm(N.fix, k, theta)
Parameters:
Estimate Std. Error t value Pr(|t|)
k 49.0724 6.9153 7.096 2.53e-06
theta -4.6676 5.1805 -0.9010.381
.lin -0.2333 0.1685 -1.3850.185
Residual standard error: 0.3662 on 16 degrees of freedom
Correlation of Parameter Estimates:
k theta
theta 0.8353
.lin -0.3458 -0.7104
R : Copyright 2004, The R Foundation for Statistical Computing
Version 2.0.1 (2004-11-15), ISBN 3-900051-07-0
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for a HTML browser interface to help.
Type 'q()' to quit R.
dat - read.table(/tmp/tlm.dat, header = TRUE)
dat
r N.fix
1 -0.5240708576
2 0.1053605245
3 -0.1743533950
4 0.1941560142
5 0.4870149851
6 -0.5068176083
7 -0.1743533950
8 0.5527898242
9 0.1521918273
10 0.4989911785
11 0.10821358 140
12 -0.83034830 156
13 -0.3074847068
14 -0.2231435550
15 0.0487901640
16 -0.0487901642
17 0.7537718040
18 -0.1251631485
19 -0.3662443975
plot(r ~ N.fix, dat)
tlm - deriv( ~ 1-((N/k)^theta), c(k, theta), function(N, k, theta){})
tlm
function (N, k, theta)
{
.expr1 - N/k
.expr2 - .expr1^theta
.value - 1 - .expr2
.grad - array(0, c(length(.value), 2), list(NULL, c(k,
theta)))
.grad[, k] - .expr1^(theta - 1) * (theta * (N/k^2))
.grad[, theta] - -(.expr2 * log(.expr1))
attr(.value, gradient) - .grad
.value
}
fm1 - nls(r ~ tlm(N.fix, k, theta), dat, start=c(k = 65, theta = 1), alg =
plinear, trace = TRUE)
2.347334 : 65.000 1.000 0.3836488
2.233973 : 56.633165 -0.431863 -1.282247
2.187513 : 54.8608041 -1.3747423 -0.4908351
2.167521 : 52.8247351 -2.1635728 -0.3592099
2.155828 : 51.4475043 -2.8951236 -0.3003772
2.149315 : 50.4174840 -3.5536098 -0.2682023
2.146427 : 49.7548459 -4.0636580 -0.257
2.145514 : 49.3861622 -4.3828441 -0.2405601
2.145301 : 49.2043941 -4.5476404 -0.2362051
2.145260 : 49.1221764 -4.6225967 -0.2343549
2.145252 : 49.0870152 -4.6544862 -0.2335967
2.145251 : 49.0724421 -4.6676430 -0.2332897
newdat - list(N.fix = seq(40, 160))
[R] corCompSymm in nlme package
We are trying to use the corCompSymm function in nlme The example from the help pages for the corAR1 function gives the following corAR1(0.2, form = ~ 1 | Mare) Correlation structure of class corAR1 representing Phi 0.2 We are expecting a somewhat similar correlation specification with the help page example for corCompSymm, but just get an error instead corCompSymm(0.5, form = ~ 1 | Subject) Error in names-.default(`*tmp*`, value = Rho) : names attribute [1] must be the same length as the vector [0] this error is concistent in R1.8 and R2.0 on windows XP and also R2.0.1 on Redhat Thanks in advance for any help Dave and Jerome Concerns about content should be sent to [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: RODBC works as documented (was [R] RODBC type conversion bug)
On Wed, 23 Feb 2005, Omar Lakkis wrote: I run R 2.0.1 on Debian and connect to Informix database via RODBC. In the table below the column month is of type char(1). RODBC seems to be converting this column to boolean if the value is F or T. Sounds reasonable. So would read.table, and that equally is not a bug but as documented. This is the data in my table: id month year 25 F 2005 26 Z 2005 When I select * for id 25 I get 25 FALSE 2005 When I select * for id 16 I get 26Z 2005 Is there a fix for this issue? Yes, for you to read the help page (as the posting guide asks you to do before posting). See ?sqlGetResults and its description of 'as.is'. -- 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
[R] Problem saving logic regression result equation to disk file
I want to get some simple logic regression examples to work before exploring a hard problem. I can get results, but I'm having some problems using cat to save the logic regression equation to a disk file. Consider this: # Simple Logic Regression Example # efg, 23 Feb 2005 library(LogicReg) # Create simulated data with known logic equation: # noise logic matrix X - matrix(as.numeric(runif(160) 0.5), 20,8) colnames(X) - paste(X, 1:ncol(X), sep=) rownames(X) - paste(case, 1:nrow(X), sep=) # Define expected result: Y = (NOT X2) AND X6 Y - as.numeric(!X[,2] X[,6]) # set seed for reproducible test set.seed(19937) # 100 interations too few: some results in single node with |Parameter| 1 Annealing - logreg.anneal.control(start = -1, end = -4, iter = 500, update = 50) logicfit - logreg(resp=Y, bin=X, type = REGRESSION.TYPE-2, select = FIT.SINGLE.MODEL-1, ntrees=1, nleaves=2, # force shape of final tree anneal.control=Annealing) # I don't always want to see the plot plot(logicfit) # I'd like to write my regression equation to a file and # then run many times to test my parameter selection # with a known case before exploring unknown cases logicfit # In this case I want either of these equivalent answers # (equivalent via DeMorgan's Theorem), and no others, # such as single node results. I want to run this say 100s (later 1000s) of times and look at the variation in the results. I want to figure out what parameters I should use so I only see these results: score 0 +1 * (X6 and (not X2)) -1 * ((not X6) or X2) # I can't use cat to write this model to a file: cat(logicfit) Error in cat(list(...), file, sep, fill, labels, append) : argument 1 not yet handled by cat summary(logicfit) Length Class Mode nsample 1-none- numeric nbinary 1-none- numeric nseparate1-none- numeric type 1-none- character select 1-none- character anneal.control 5-none- list tree.control 4-none- list seed 1-none- numeric choice 1-none- numeric nleaves 1-none- numeric ntrees 1-none- numeric penalty 1-none- numeric response20-none- numeric binary 160-none- numeric separate 1-none- numeric censor 20-none- numeric weight 20-none- numeric model5logregmodel list call 8-none- call # Just the logicfit$model would be good enough but I can't cat that either: logicfit$model +1 * (X6 and (not X2)) cat(logicfit$model) Error in cat(list(...), file, sep, fill, labels, append) : argument 1 not yet handled by cat Using sink to get this result seems to be a huge kludge: sink(saveresults.txt) logicfit$model sink() results - readLines(saveresults.txt) results [1] +1 * (X6 and (not X2)) # FINALLY something I could write this result to a file:. cat(results, \n) +1 * (X6 and (not X2)) What is a simple way to get my logic regression equation as a string that I can cat without dealing with the internal data structures that are present here? Thanks for any help with this. efg -- Earl F. Glynn Scientific Programmer Bioinformatics Department Stowers Institute for Medical Research __ R-help@stat.math.ethz.ch 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] stopping a function
On Wed, 23 Feb 2005, Brian S Cade wrote: I've looked for this information in all the R help sources I could find and found nothing. Is it possible to use some function key to stop the execution of some R command without ending the R session (Windows, R 1.91)? R 1.91 is not expected until 2050, if ever. I've several times started functions that for various reasons are not executing properly and it would be nice to stop them without killing the R session. I've been using taskmgr to end the R session but I then lose all objects created during the session. It's in the README, and for RGui, on the menus! Esc in RGui, Ctrl-c in Rterm. Not all computations are interruptible (most compiled code is not), but if they are, these will do it. -- 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
[R] Course***R/S System: Advanced Programming, Washington, DC***
XSolutions Corp (www.xlsolutions-corp.com) is proud to announce our Advanced R/Splus programming course taught by R Development Core Team Guru! www.xlsolutions-corp.com/Radv.htm *Washington, DC April 14th-15th, 2005 *Boston, MA-- TBD Early-bird discount ends March 15th! Ask for group discount and reserve your seat Now (payment due after the class) Email Sue Turner: [EMAIL PROTECTED] Phone: 206-686-1578 Course Outline: - Overview of R/S fundamentals: Syntax and Semantics - Class and Inheritance in R/S-Plus - Concepts, Construction and good use of language objects - Coercion and efficiency - Object-oriented programming in R and S-Plus - Advanced manipulation tools: Parse, Deparse, Substitute, etc. - How to fully take advantage of Vectorization - Generic and Method Functions; S4 (S-Plus 6) - Search path, databases and frames Visibility - Working with large objects - Handling Properly Recursion and iterative calculations - Managing loops; For (S-Plus) and for() loops - Consequences of Lazy Evaluation - Efficient Code practices for large computations - Memory management and Resource monitoring - Writing R/S-Plus functions to call compiled code - Writing and debugging compiled code for R/S-Plus system - Connecting R/S-Plus to External Data Sources - Understanding the structure of model fitting functions in R/S-Plus - Designing and Packaging efficiently a new model function It'll also deal with lots of S-Plus efficiency issues and any special topics from participants is welcome. Please let us know if you and your colleagues are interested in this class to take advantage of group discount. Over half of the seats in this class are currently reserved. Register now to secure your seat in this course! Cheers, Elvis Miller, PhD Manager Training. XLSolutions Corporation 206 686 1578 www.xlsolutions-corp.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] stopping a function
In Rgui, use the menu `Help' - `Console' and look at the third line from the bottom. Cheers, Andy From: Brian S Cade I've looked for this information in all the R help sources I could find and found nothing. Is it possible to use some function key to stop the execution of some R command without ending the R session (Windows, R 1.91)? I've several times started functions that for various reasons are not executing properly and it would be nice to stop them without killing the R session. I've been using taskmgr to end the R session but I then lose all objects created during the session. Brian Brian S. Cade U. S. Geological Survey Fort Collins Science Center 2150 Centre Ave., Bldg. C Fort Collins, CO 80526-8818 email: [EMAIL PROTECTED] tel: 970 226-9326 __ R-help@stat.math.ethz.ch 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 conctruct an inner grouping for nlme random statement?
JJ wrote: Hello. Im hoping someone can help with a grouping question related to the random= statement within the nlme function. How do you specify that some grouping levels are inner to others? I tried several things, given below. Lets say I have a data frame with five variables, resp, cov1, ran1, ran2, group1, and group 2. The formula is resp~cov1 + ran1 + ran2, where the ran are random variables. The data is of length 80, and there are 4 unique factors in group1 and 20 unique factors in group2. These are factors related to ran1 and ran2, respectively. The difficult part is that I want to estimate only 4 random variables for ran1|group1 and the full 20 for ran2|group2. I have tried many ways, and I cannot find a way to do this. Is there a way? Can someone suggest a code snippet? First I tried making the data frame a groupedData object, so that group2 is inner to group1, as it should be. Then I used the statement: random = as.formula(ran1+ran2~1). But this produced 20 estimates for both ran1 and ran2. I have also tried it without the data frame as a groupedData object, using the following: random = list(group1= c(ran1~1, group2=ran2~1)). But this gave only 4 estimates for ran2. I also tried: random = list(c(group1= ran1~1, group2=ran2~1)), but this just gave a parse error message. Any suggestions would be greatly appreciated. Is it even possible to do what I want to do? John I think we will need a bit more information before we are able to help. I'm not sure what you mean by ran1 and ran2 being random variables and group1 and group2 being factors related to ran1 and ran2. The lmer function in the lme4 package allows for specification of a mixed-effects model using grouping factors and model matrices. The model matrix determines the form of the random effects vector corresponding to each of the groups. The grouping factor determines the groups. In this specification I'm not sure what your ran1 and ran2 variables would be. It sounds as if you want a model that would be specified as fm1 - lmer(resp ~ cov1 + (1|group1) + (1|group2)) but I'm not sure. Feel free to correspond directly with me if you wish. __ R-help@stat.math.ethz.ch 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] BLAS or ATLAS?
I'm trying to tweak multi-threading on a multiprocessor box(Xeon). Reading appendix a.2.2 in the install and admin menu, I see that I can use the internal BLAS lib and/or a multi-thread version of ATLAS (Where to get?) Question: Intel claims their libraries are faster than ATLAS. Has anyone compiled R against the Intel BLAS? If so...how did you do it? http://developer.intel.com/software/products/mkl/features/lin_alg.htm And I'm in the market for a pre-built Atlas dual processor lib for Xeon. Know where I can find one? Thanks. -- Tom Colson Center for Earth Observation North Carolina State University Raleigh, NC 27695 (919) 515 3434 (919) 673 8023 [EMAIL PROTECTED] Online Calendar: http://www4.ncsu.edu/~tpcolson __ R-help@stat.math.ethz.ch 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] model.matrix for a factor effect with no intercept
Brian, many thanks for these helpful pointers: On 23 Feb 2005, at 15:45, Prof Brian Ripley wrote: MASS4 p.150 White Book p.38 Those are the only two reasonably comprehensive accounts that I am aware of (and they have only partial overlap). The underlying motivation is to span the _additional_ vector space covered by the term, the complement to what has gone before. Put another way, as each term is added, only enough columns are added to the model matrix to span the same space as if dummy coding had been used for that term and its predecessors. So think of this as a way to produce a parsimonious (usually full-rank) basis for the model space. Yes, indeed. My surprise was that *this* particular basis (dummy coding) was the one used here. I should have got a clue from what contrasts() does, and is documented to do, options(contrasts) $contrasts [1] contr.treatment contr.poly contrasts(a) .L .Q [1,] -7.071068e-01 0.4082483 [2,] -9.073800e-17 -0.8164966 [3,] 7.071068e-01 0.4082483 but when there's no intercept in the model the contrasts used appear to be contrasts(a, contrasts = FALSE) -1 0 1 -1 1 0 0 0 0 1 0 1 0 0 1 which are not the same as contr.poly(a, contrasts = FALSE) ^0^1 ^2 [1,] 1 -7.071068e-01 0.4082483 [2,] 1 -9.073800e-17 -0.8164966 [3,] 1 7.071068e-01 0.4082483 -- which is what I had naively expected to get in my model matrix. This is of course all a matter of convention. The present convention does seem a touch confusing though: the basis for the space spanned by a factor is determined by options(contrasts) or by the contrasts attribute of the factor or by the contrasts argument in the call, *except* when there's no intercept or other factor earlier in the model in which case all such settings are ignored (well, not *quite* ignored: they do get put in the contrasts attribute of the resultant model matrix). On the other hand, one good reason to use dummy coding for the first-encountered factor when there's no intercept is that the associated parameters are then often interpretable as group-specific intercepts. Would it be an improvement, though, if the contrasts attribute of the resultant model matrix contained contr.treatment in such cases instead of the name of a contrast function that was not actually used? Best wishes, David On Wed, 23 Feb 2005, David Firth wrote: I was surprised by this (in R 2.0.1): a - ordered(-1:1) a [1] -1 0 1 Levels: -1 0 1 model.matrix(~ a) (Intercept) a.La.Q 1 1 -7.071068e-01 0.4082483 2 1 -9.073800e-17 -0.8164966 3 1 7.071068e-01 0.4082483 attr(,assign) [1] 0 1 1 attr(,contrasts) attr(,contrasts)$a [1] contr.poly model.matrix(~ -1 + a) a-1 a0 a1 1 1 0 0 2 0 1 0 3 0 0 1 attr(,assign) [1] 1 1 1 attr(,contrasts) attr(,contrasts)$a [1] contr.poly Without the intercept, treatment contrasts seem to have been used (this despite the contr.poly in the contrasts attribute). It's not restricted to ordered factors. For example, if Helmert contrasts are used for nominal factors, the same sort of thing happens. I suppose it is a deliberate feature (perhaps to protect the user from accidentally fitting models that make no sense? or maybe some better reason?) -- is it explained somewhere? -- 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] BLAS or ATLAS?
What OS is this for? Several OSes (Windows, Linux, Solaris, FreeBSD ...) run on Xeons (and there are several sorts of Xeons, some even EM64T AMD-clones). On Wed, 23 Feb 2005, Tom Colson wrote: I'm trying to tweak multi-threading on a multiprocessor box(Xeon). Reading appendix a.2.2 in the install and admin menu, I see that I can use the internal BLAS lib and/or a multi-thread version of ATLAS (Where to get?) The last _is_ described in Appendix A.2.2, together with several other possibilities. Question: Intel claims their libraries are faster than ATLAS. Has anyone compiled R against the Intel BLAS? If so...how did you do it? http://developer.intel.com/software/products/mkl/features/lin_alg.htm That _is_ described in Appendix A.2.2. And I'm in the market for a pre-built Atlas dual processor lib for Xeon. Know where I can find one? -- 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] BLAS or ATLAS?
Fedora Core 3 I installed the Intel MLk, and tried : --with-blas=-lmkl -lguide -lpthread and got : External libraries:readline, BLAS(generic) thus I'm assuming I'm missing something when telling config where to look for the recently installed Intel BLAS? thanks for all replies thus far. On Wed, 2005-02-23 at 19:25 +, Prof Brian Ripley wrote: What OS is this for? Several OSes (Windows, Linux, Solaris, FreeBSD ...) run on Xeons (and there are several sorts of Xeons, some even EM64T AMD-clones). On Wed, 23 Feb 2005, Tom Colson wrote: I'm trying to tweak multi-threading on a multiprocessor box(Xeon). Reading appendix a.2.2 in the install and admin menu, I see that I can use the internal BLAS lib and/or a multi-thread version of ATLAS (Where to get?) The last _is_ described in Appendix A.2.2, together with several other possibilities. Question: Intel claims their libraries are faster than ATLAS. Has anyone compiled R against the Intel BLAS? If so...how did you do it? http://developer.intel.com/software/products/mkl/features/lin_alg.htm That _is_ described in Appendix A.2.2. And I'm in the market for a pre-built Atlas dual processor lib for Xeon. Know where I can find one? __ R-help@stat.math.ethz.ch 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] BLAS or ATLAS?
I believe that actually means configure was able to use the BLAS you specified (MKL in this case). You can scan through config.log to make sure. Andy From: Fedora Core 3 I installed the Intel MLk, and tried : --with-blas=-lmkl -lguide -lpthread and got : External libraries:readline, BLAS(generic) thus I'm assuming I'm missing something when telling config where to look for the recently installed Intel BLAS? thanks for all replies thus far. On Wed, 2005-02-23 at 19:25 +, Prof Brian Ripley wrote: What OS is this for? Several OSes (Windows, Linux, Solaris, FreeBSD ...) run on Xeons (and there are several sorts of Xeons, some even EM64T AMD-clones). On Wed, 23 Feb 2005, Tom Colson wrote: I'm trying to tweak multi-threading on a multiprocessor box(Xeon). Reading appendix a.2.2 in the install and admin menu, I see that I can use the internal BLAS lib and/or a multi-thread version of ATLAS (Where to get?) The last _is_ described in Appendix A.2.2, together with several other possibilities. Question: Intel claims their libraries are faster than ATLAS. Has anyone compiled R against the Intel BLAS? If so...how did you do it? http://developer.intel.com/software/products/m kl/features/lin_alg.htm That _is_ described in Appendix A.2.2. And I'm in the market for a pre-built Atlas dual processor lib for Xeon. Know where I can find one? __ R-help@stat.math.ethz.ch 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] BLAS or ATLAS?
In makeconf I get BLAS_LIBS = -lblas I believe that actually means configure was able to use the BLAS you specified (MKL in this case). You can scan through config.log to make sure. in config.logcan't seem to find any reference to the Intel BLAS. Thanks for all the replies. On Wed, 2005-02-23 at 13:45 -0600, Douglas Bates wrote: Tom Colson wrote: Fedora Core 3 I installed the Intel MLk, and tried : --with-blas=-lmkl -lguide -lpthread and got : External libraries:readline, BLAS(generic) thus I'm assuming I'm missing something when telling config where to look for the recently installed Intel BLAS? That's what you would expect to see in the summary section from configure. Look back in the configure output to see exactly what it reported for the blas. You may also check the Makeconf file that was generated to see what it defines for BLAS_LIBS __ R-help@stat.math.ethz.ch 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] bias of a boot statistic
The bias and standard error are calculated in the print.boot function and are not saved anywhere. It is very easy to calculate, however, as the mean of the replicate statistics minus the original statistic. apply(boot.out$t,2,mean)-boot.out$t0 HTH, Angelo On Wed, 23 Feb 2005, K. Steinmann wrote: Question: How can I get access to the bias value of a boot statistic? Details: Boot function: boot(data, statistic, R, sim=ordinary, stype=i, strata=rep(1,n), L=NULL, m=0, weights=NULL, ran.gen=function(d, p) d, mle=NULL, ...) When I create an object, containing the bootstrap statistic (object - boot ())I can call it and will get an output with t, bias and standarderror as follows: Bootstrap Statistics : original biasstd. error t1* 5.650.01 0.9134185 My question is now, where is the value of the bias stored? How can I get access to this value to do further caluculations? Thanks, K. St. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- -- | Angelo J. CantyEmail: [EMAIL PROTECTED] | | Mathematics and Statistics Phone: (905) 525-9140 x 27079 | | McMaster UniversityFax : (905) 522-0935 | | 1280 Main St. W. | | Hamilton ON L8S 4K1 | __ R-help@stat.math.ethz.ch 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] Solving systems of non-linear equations in R
Thank you, that was very helpful. My functions are in general monotonic, continous and differentiable(one exception sometimes encountered being y = min(a*x1,b*x2)) and do have a unique solution, if you specify the problem correctly. I have never worked with non-liner solving algoritms in my math courses, so maybe starting to setup algoritms in R may take the focus away from my thesis subject. On the other hand it could be fun :-D Maybe I should check out the optim()-function again. You can solve systems with more tha one equations with optim()? regards Spencer Graves wrote: A system of n equations in n unknowns has a unique solution if the n equations are linear and linearly independent. If the system is nonlinear, then one must characterize the nonlinearity before saying anything about whether a solution exists and if so how many solutions are there? Example 1: Solve sin(x)=0 for x. Answer: x = 2*n*pi, for n = any integer. Example 2: Solve sin(x) = 2 for x. Answer: If x must be a real number, then this equation has no solutions. Are your functions monotonic? Continuous? Differentiable? Without getting into pathologies like the Cantor function (e.g., http://www.cut-the-knot.org/do_you_know/cantor.shtml), my experience with a variety of practical problem like this suggests that it is best to recast the problem as one of minimizing, e.g., the sum of squared deviations from target. Moreover, I've had good luck transforming the parameter space to eliminate constraints -- or incorporating the constraints into the objective function and then solving the superficially unconstrained problem. If my functions have singularities where I might get 0/0 or Inf-Inf, for example, I use asymptotic expansions to approximate the function(s) near the singularities more accurately than can be achieved with any finite-precision arithmetic. With all of this, I'm confident that there are better algorithms than the different methods in optim, but I don't have not had the need to hunt for them. The methods in optim provide a reasonable range of options for the problems I've encountered. The R project has another advantage over a commercial software: You can see the source code. You can trace it step by step and find out where it does not work well for the specific problems you consider. If you're clever, you might be able to find a way to improve that algorithm and make it part of your thesis -- and get a publication on it in some statistical software journal. Where else can you so easily climb up and stand on the shoulders of giants? If you find a platform for innovation better than R, please let me know. hope this helps. spencer graves T Petersen wrote: No, this doesn't seem right. What I look for is something that could solve nonlinear systems with n unknowns and n equations. So there will be zero degrees of freedom, and statistical methods can't be the right way forward. Specifically I can see that the litterature mentions's Scarf's algoritm (Scarf 1967) and Merril's refinement of Scarf's algoritm in 1972, but there might be other algoritms too... Regards...TP yutaka hamaoka wrote: I believe library(systemfit) has nlsytemfit function. Yh T Petersen wrote: I'm about to write my thesis in economics and will need to setup and solve a system of non-linear equations. At our university we usually use GAMS for this, and though GAMS is a fine program, it bugs me a that I won't be able to run my code after I finish my thesis without buying a license for the program(about $3.500 :-(( ) So I've looked around for NL-stuff for R, but I can't find anything. The closest thing appears to be optim(), but it doesn't seem to allow constraints(as in fn = constant) or equations systems. So, anyone knows if there is a method in R that you can use for this purpose? regards __ R-help@stat.math.ethz.ch 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-help@stat.math.ethz.ch 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] BLAS or ATLAS?
Look for `sgemm' in config.log. I believe that's how configure test for usable BLAS. It looks like configure found another system installed BLAS. Andy From: Tom Colson In makeconf I get BLAS_LIBS = -lblas I believe that actually means configure was able to use the BLAS you specified (MKL in this case). You can scan through config.log to make sure. in config.logcan't seem to find any reference to the Intel BLAS. Thanks for all the replies. On Wed, 2005-02-23 at 13:45 -0600, Douglas Bates wrote: Tom Colson wrote: Fedora Core 3 I installed the Intel MLk, and tried : --with-blas=-lmkl -lguide -lpthread and got : External libraries:readline, BLAS(generic) thus I'm assuming I'm missing something when telling config where to look for the recently installed Intel BLAS? That's what you would expect to see in the summary section from configure. Look back in the configure output to see exactly what it reported for the blas. You may also check the Makeconf file that was generated to see what it defines for BLAS_LIBS __ R-help@stat.math.ethz.ch 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] how do I get means by factor?
I want the mean of a response variable for each level of a factor in a data.frame. How is this done? [[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 do I get means by factor?
tapply( x, f, mean ) maybe with tapply( x, f, mean, na.rm=T ) -- Bendix Carstensen Senior Statistician Steno Diabetes Center Niels Steensens Vej 2 DK-2820 Gentofte Denmark tel: +45 44 43 87 38 mob: +45 30 75 87 38 fax: +45 44 43 07 06 [EMAIL PROTECTED] www.biostat.ku.dk/~bxc -- -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Chris Sent: Wednesday, February 23, 2005 10:38 PM To: r-help@stat.math.ethz.ch Subject: [R] how do I get means by factor? I want the mean of a response variable for each level of a factor in a data.frame. How is this done? [[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] model.matrix for a factor effect with no intercept
David,
I am sorry, I did not read your question as `why dummy coding not poly
coding', but I still seemed to have answered it. BTW,
contr.poly(a, contrasts = FALSE)
should be an indicator matrix according to its help page, so I would not
have expected what you expected, and if it had been, you might not have
been surprised. [More on that below.]
I would say the algorithm was `code by the specified contrasts only if
dummy coding is not parsimonious'.
Your final point isn't the whole story. If I have 0+a+b+a:b, R uses the
contrasts for a in the interaction but not the main term. So the
attribute contrasts for 'a' represents
`the contrast used for those terms (possibly none) in which it needed
coding'.
Now, if there were none it might be better to say nothing, but that is not
actually calculated explicitly. I've just checked, and we do not document
on the help page what it means. I have just drafted
If there are any factors in terms in the model, there is an attribute
\code{contrasts}, a named list with an entry for each factor. This
specifies the contrasts that would be used in terms in which the
factor is coded by contrasts (in some terms dummy coding may be used),
either as a character vector naming a function or as a numeric matrix.
As far as I recall the contrasts attribute is returned to be tacked on a
fit and used later to reconstruct the model matrix, or perhaps a model
matrix for another formula and the same data. So I think it does need to
be the contrasts as specified at the time of the model.matrix call.
It is probably helpful to note that it is not specifically treatment
contrasts that are used but dummy coding: contr.sum and contr.helmert give
the same matrix.
It is a bug that contr.poly(x, contrasts=FALSE) does not behave as
documented (a bug shared with S, it seems). I do dimly recall this having
been raised in the past (probably by me), but cannot find any record as to
why the discrepancy remains. We should change the code or the
documentation, probably the latter.
Brian
On Wed, 23 Feb 2005, David Firth wrote:
Brian, many thanks for these helpful pointers:
On 23 Feb 2005, at 15:45, Prof Brian Ripley wrote:
MASS4 p.150
White Book p.38
Those are the only two reasonably comprehensive accounts that I am aware of
(and they have only partial overlap).
The underlying motivation is to span the _additional_ vector space covered
by the term, the complement to what has gone before. Put another way, as
each term is added, only enough columns are added to the model matrix to
span the same space as if dummy coding had been used for that term and its
predecessors. So think of this as a way to produce a parsimonious (usually
full-rank) basis for the model space.
Yes, indeed. My surprise was that *this* particular basis (dummy coding) was
the one used here. I should have got a clue from what contrasts() does, and
is documented to do,
options(contrasts)
$contrasts
[1] contr.treatment contr.poly
contrasts(a)
.L .Q
[1,] -7.071068e-01 0.4082483
[2,] -9.073800e-17 -0.8164966
[3,] 7.071068e-01 0.4082483
but when there's no intercept in the model the contrasts used appear to be
contrasts(a, contrasts = FALSE)
-1 0 1
-1 1 0 0
0 0 1 0
1 0 0 1
which are not the same as
contr.poly(a, contrasts = FALSE)
^0^1 ^2
[1,] 1 -7.071068e-01 0.4082483
[2,] 1 -9.073800e-17 -0.8164966
[3,] 1 7.071068e-01 0.4082483
-- which is what I had naively expected to get in my model matrix.
This is of course all a matter of convention. The present convention does
seem a touch confusing though: the basis for the space spanned by a factor is
determined by options(contrasts) or by the contrasts attribute of the
factor or by the contrasts argument in the call, *except* when there's no
intercept or other factor earlier in the model in which case all such
settings are ignored (well, not *quite* ignored: they do get put in the
contrasts attribute of the resultant model matrix).
On the other hand, one good reason to use dummy coding for the
first-encountered factor when there's no intercept is that the associated
parameters are then often interpretable as group-specific intercepts.
Would it be an improvement, though, if the contrasts attribute of the
resultant model matrix contained contr.treatment in such cases instead of
the name of a contrast function that was not actually used?
Best wishes,
David
On Wed, 23 Feb 2005, David Firth wrote:
I was surprised by this (in R 2.0.1):
a - ordered(-1:1)
a
[1] -1 0 1
Levels: -1 0 1
model.matrix(~ a)
(Intercept) a.La.Q
1 1 -7.071068e-01 0.4082483
2 1 -9.073800e-17 -0.8164966
3 1 7.071068e-01 0.4082483
attr(,assign)
[1] 0 1 1
attr(,contrasts)
attr(,contrasts)$a
[1] contr.poly
model.matrix(~ -1 + a)
a-1 a0 a1
1 1 0 0
2 0 1 0
3 0 0 1
attr(,assign)
[1] 1 1 1
attr(,contrasts)
attr(,contrasts)$a
[1]
[R] Pls unsubscribe me from this list.
Hi, I have difficulty unsubscribing from this list. I have already unsubscribed from the website mailing list. But the messages keep flooding my mailbox. Pls help! [[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] basic question about changing limits on generated plots
Is it possible to change the limits on plots that are already on the screen? In particular, is there any R equivalent to the sequence of matlab commands plot(1:10,1:10) hold on plot(2:12,5:15) I know I can use points and lines to add points and lines to plots, but the limits of the plot do not change when I do this. Looking at various examples, it seems that the answer is no, but I wanted to make sure. This seems to make exploratory data analysis somewhat more challenging. It also seems like it would be plausible to write a package on top of the standard graphics functions that keeps track of what you've done and automatically replots. Which makes me think, has someone already done this? Cheers, rif __ R-help@stat.math.ethz.ch 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] Pls unsubscribe me from this list.
Lim, Jasper H. N. [EMAIL PROTECTED] writes: Hi, I have difficulty unsubscribing from this list. I have already unsubscribed from the website mailing list. But the messages keep flooding my mailbox. Pls help! [[alternative HTML version deleted]] The mailing list manager is out of town for another couple of days. The other two thousand recipients of your message cannot help you. If you really cannot unsubscribe using the web interface, the list manager's address is found on the same page. Please notice that unsubscribing does not dequeue messaged that have already been sent to you so you may receive messages for some time afterwards, even though you are no longer on the list. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- O__ Peter Dalgaard Blegdamsvej 3 c/ /'_ --- Dept. of Biostatistics 2200 Cph. N (*) \(*) -- 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] Solving systems of non-linear equations in R
Write a driver function to compute the sum of squares of deviations from target. If the nonlinear equations are not pathological and a unique solution exists, optim will find it. If no solution exists, optim will find something close -- in terms of the sum of squared deviations. If the functions have discontinuities, multiple solutions, or possibly other pathologies, I'd try all 4 methods (after transforming parameters to eliminate constraints). You could also try to generate random starting points over a plausible region, evaluate your objective function at all the random starting points and feed the best one (or the best 10) to optim. How big is n? If it's bigger than, say, half a dozen, I suspect you have a variance component situation, in which case, you would be wise to use nlme; see Pinheiro and Bates (2000) Mixed-Effects Models in S and S-Plus (Springer). Or try Scilab (http://scilabsoft.inria.fr). hope this helps. spencer graves T Petersen wrote: Thank you, that was very helpful. My functions are in general monotonic, continous and differentiable(one exception sometimes encountered being y = min(a*x1,b*x2)) and do have a unique solution, if you specify the problem correctly. I have never worked with non-liner solving algoritms in my math courses, so maybe starting to setup algoritms in R may take the focus away from my thesis subject. On the other hand it could be fun :-D Maybe I should check out the optim()-function again. You can solve systems with more tha one equations with optim()? regards Spencer Graves wrote: A system of n equations in n unknowns has a unique solution if the n equations are linear and linearly independent. If the system is nonlinear, then one must characterize the nonlinearity before saying anything about whether a solution exists and if so how many solutions are there? Example 1: Solve sin(x)=0 for x. Answer: x = 2*n*pi, for n = any integer. Example 2: Solve sin(x) = 2 for x. Answer: If x must be a real number, then this equation has no solutions. Are your functions monotonic? Continuous? Differentiable? Without getting into pathologies like the Cantor function (e.g., http://www.cut-the-knot.org/do_you_know/cantor.shtml), my experience with a variety of practical problem like this suggests that it is best to recast the problem as one of minimizing, e.g., the sum of squared deviations from target. Moreover, I've had good luck transforming the parameter space to eliminate constraints -- or incorporating the constraints into the objective function and then solving the superficially unconstrained problem. If my functions have singularities where I might get 0/0 or Inf-Inf, for example, I use asymptotic expansions to approximate the function(s) near the singularities more accurately than can be achieved with any finite-precision arithmetic. With all of this, I'm confident that there are better algorithms than the different methods in optim, but I don't have not had the need to hunt for them. The methods in optim provide a reasonable range of options for the problems I've encountered. The R project has another advantage over a commercial software: You can see the source code. You can trace it step by step and find out where it does not work well for the specific problems you consider. If you're clever, you might be able to find a way to improve that algorithm and make it part of your thesis -- and get a publication on it in some statistical software journal. Where else can you so easily climb up and stand on the shoulders of giants? If you find a platform for innovation better than R, please let me know. hope this helps. spencer graves T Petersen wrote: No, this doesn't seem right. What I look for is something that could solve nonlinear systems with n unknowns and n equations. So there will be zero degrees of freedom, and statistical methods can't be the right way forward. Specifically I can see that the litterature mentions's Scarf's algoritm (Scarf 1967) and Merril's refinement of Scarf's algoritm in 1972, but there might be other algoritms too... Regards...TP yutaka hamaoka wrote: I believe library(systemfit) has nlsytemfit function. Yh T Petersen wrote: I'm about to write my thesis in economics and will need to setup and solve a system of non-linear equations. At our university we usually use GAMS for this, and though GAMS is a fine program, it bugs me a that I won't be able to run my code after I finish my thesis without buying a license for the program(about $3.500 :-(( ) So I've looked around for NL-stuff for R, but I can't find anything. The closest thing appears to be optim(), but it doesn't seem to allow constraints(as in fn = constant) or equations systems. So, anyone knows if there is a method in R that you can use for this purpose? regards __
[R] and [ESS] Starting ESS
Dear R People: I have finally seen the error of my ways and have decided to use ESS for R and S + stuff. However, I have a question right from the beginning. I'm somewhat confused by the installation instructions. Do I install XEMACS or ESS first, please? Windows XP R Version 2.0.1 (S + 6.2) Thanks so much! Sincerely, Laura Holt mailto: [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] Getting tick positions
While writing a function that includes placing grid lines at the same position as the axis ticks, I found that the axis* functions don't return anything. Thus I have had to copy the appropriate function, removing the call to axis() and adding a line to return the tick positions. Is there a more elegant way to determine the tick positions on an axis? Thanks. Jim (normally [EMAIL PROTECTED] - my modem is toast) __ R-help@stat.math.ethz.ch 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] Pls unsubscribe me from this list.
Two other possibilities come to my mind : 1) If memory serves, when you un-subscribe using the web interface, you will be sent an email asking for confirmation. Till you click on the link in that email, you will not be un-subscribed. Did you you receive any such mail ? If not, you might want to check your junk mail folder or spam settings. 2) You can disable the Mail delivery option under R-help subscription options in web interface. Regards, Adai On Thu, 2005-02-24 at 00:06 +0100, Peter Dalgaard wrote: Lim, Jasper H. N. [EMAIL PROTECTED] writes: Hi, I have difficulty unsubscribing from this list. I have already unsubscribed from the website mailing list. But the messages keep flooding my mailbox. Pls help! [[alternative HTML version deleted]] The mailing list manager is out of town for another couple of days. The other two thousand recipients of your message cannot help you. If you really cannot unsubscribe using the web interface, the list manager's address is found on the same page. Please notice that unsubscribing does not dequeue messaged that have already been sent to you so you may receive messages for some time afterwards, even though you are no longer on the list. __ R-help@stat.math.ethz.ch 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] and [ESS] Starting ESS
AFAIK, it does not matter as long as you got the correct paths set in your init.el file. Please see the excellent installation instructions from John Fox's http://socserv.mcmaster.ca/jfox/Books/Companion/ESS/ Regards, Adai On Wed, 2005-02-23 at 17:16 -0600, Laura Holt wrote: Dear R People: I have finally seen the error of my ways and have decided to use ESS for R and S + stuff. However, I have a question right from the beginning. I'm somewhat confused by the installation instructions. Do I install XEMACS or ESS first, please? Windows XP R Version 2.0.1 (S + 6.2) Thanks so much! Sincerely, Laura Holt mailto: [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] baseline shift / detrend
Hello, I am analyzing data of biological motion (mouth movements) over time (~ 10 mins). There is a considerable baseline shift which appears to be a linear trend. Is it possible to get rid of this baseline shift using R? Any help is appreciated! Peter Dr. Peter Soros Sunnybrook and Women's College Health Sciences Centre 2075 Bayview Avenue Room S646 Toronto, Ontario, M4N 3M5 Canada [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] colorbar for image
Klaus Ladner wrote: Dear colleges, Does anyone know, how to insert a color bar as used with filled.contour when using image? In the package plotrix, there is a function gradient.rect() that might do what you want. Jim This message was sent through MyMail http://www.mymail.com.au __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Getting tick positions
Which plotting function are you using ? I think most of plotting can accept xaxt=n which is the command to supress the x-axis. If this works, at least you do not have to redefine the function. Examples plot(1:10, xaxt=n) hist( rnorm(100), xaxt=n ) boxplot( rnorm(10), rnorm(10), rnorm(10), xaxt=n ) Reading help(par) will also shows you that xaxp might be useful but I have not managed to get this working. Maybe someone on the list can explain why the following does not work : plot(1:100, xaxp=c(x1=0,x2=100,n=20) ) par( xaxp=c(x1=0,x2=100,n=20) ); plot(1:100) Regards, Adai On Thu, 2005-02-24 at 10:08 +1100, fls wrote: While writing a function that includes placing grid lines at the same position as the axis ticks, I found that the axis* functions don't return anything. Thus I have had to copy the appropriate function, removing the call to axis() and adding a line to return the tick positions. Is there a more elegant way to determine the tick positions on an axis? Thanks. Jim (normally [EMAIL PROTECTED] - my modem is toast) __ R-help@stat.math.ethz.ch 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] Getting tick positions
?axis where you will find See Also: 'axTicks' returns the axis tick locations corresponding to 'at=NULL'; 'pretty' is more flexible for computing pretty tick coordinates and does _not_ depend on (nor adapt to) the coordinate system in use. Tom -Original Message- From: fls [mailto:[EMAIL PROTECTED] Sent: Thursday, 24 February 2005 7:08 AM To: r-help@stat.math.ethz.ch Subject: [R] Getting tick positions While writing a function that includes placing grid lines at the same position as the axis ticks, I found that the axis* functions don't return anything. Thus I have had to copy the appropriate function, removing the call to axis() and adding a line to return the tick positions. Is there a more elegant way to determine the tick positions on an axis? Thanks. Jim (normally [EMAIL PROTECTED] - my modem is toast) __ R-help@stat.math.ethz.ch 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
FYI from apple employee RE: [R] Memory error in Mac OS X Aqua GUI v1.01 with cluster
Please ignore this is you're not interested in what my friend who works for apple had to say regarding my issue which is explained in detail at the bottom of this exchange... me: The data from this simulation has to be processed at several stages of data processing so at some point I do have to have to make a dissimilarity matrix from this large matrix, as well as ones that I have filtered. appleguy: Your matrix is too big; adding more memory, or killing off other processes won't help. R is clearly not very smart when it comes to dealing with large chunks of data, and you're running into a fundamental upper limit on the amount of data that it can deal with. The only hope you have of dealing with this is to find a 64-bit system (many Unix boxes fit this description) and use the text or X-windows based R there. If you have a login on a Unix machine at work that has R installed, try your program there. In a few months, when Tiger is out and about, the G5 Macs will also be able to run 64-bit R. me: Yeah I was unsuccesful trying to download the unix version of r last night(which is also the windows version). Christopher said it needed the fortran compiler(?). According to takao the linux we use doesn't have any gui capability so I'm not sure I can make the dendograms I'm gonna need if I installed R on it. It Linux 2.6.7. and is some version of red hat(?). Oh well... appleguy: Being Linnex, it is probably running on a PC, so it won't be a 64-bit system anyway. What sort of timeframe do you have on this? I might be able to track down a 64-bit system, but it may take a little while. If time is tight, then you should ask around to see a) whether anyone has a 64-bit unix machine; i.e. a Sun, SGI or Alpha box, and b) whether they know how to generate the graphs without the GUI (I suspect the GUI is not required to generate them, only view them). Until then, you're going to need to reduce the size of your dataset I'm afraid. 8( OG message to rhelp: I'm sorry if the answer to my problem is buried in the archives. I have limited experience with R and I couldn't find a solution to my particular problem. I am running Mac OS X Aqua GUI v1.01 on a new G5 running os 10.3.8 with a 1.8Ghz processor and 1GB of sdram. I just downloaded bioconducter a week ago and I'm trying to cluster a matrix I created with a simulation with dimensions dim(nca35) [1] 1048112 with size object.size(nca352) [1] 1426204 I checked my ulimits variable on the unix terminal and it says it's unlimited as does mem.limits() nsize vsize NANA But I'm still getting errors like the following with funtions in the cluster package daisy(nca352, metric= euclidean, stand=FALSE)-dnca35 Error: cannot allocate vector of size 858213 Kb *** malloc: vm_allocate(size=878813184) failed (error code=3) *** malloc[599]: error: Can't allocate region if it helps i also checked gc() used (Mb) gc trigger (Mb) Ncells 448662 12.0 741108 19.8 Vcells 847630 6.5 135357901 1032.7 I tried the suggested unix command in the memory help doc but that doesn't work in the Aqua GUI. Can someone tell me how to change the Vcells? Although to the best of my understanding (which is limited) I shouldn't have to do that. Any suggestions would be greatly appreciated. thanks, betty -- Betty Gilbert [EMAIL PROTECTED] Taylor Lab Plant and Microbial Biology 321 Koshland Hall U.C. Berkeley Berkeley, Ca 94720 __ R-help@stat.math.ethz.ch 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] basic question about changing limits on generated plots
On Wed, 2005-02-23 at 17:42 -0500, rif wrote: Is it possible to change the limits on plots that are already on the screen? In particular, is there any R equivalent to the sequence of matlab commands plot(1:10,1:10) hold on plot(2:12,5:15) I know I can use points and lines to add points and lines to plots, but the limits of the plot do not change when I do this. Looking at various examples, it seems that the answer is no, but I wanted to make sure. This seems to make exploratory data analysis somewhat more challenging. It also seems like it would be plausible to write a package on top of the standard graphics functions that keeps track of what you've done and automatically replots. Which makes me think, has someone already done this? Cheers, rif I have not used Matlab, but I suspect that you want to use: par(new = TRUE) For example: plot(1:10,1:10) # Check plot region limits par(usr) [1] 0.64 10.36 0.64 10.36 # Set par(new) to not overwrite existing plot par(new = TRUE) plot(2:12,5:15) # Re-check plot region limits par(usr) [1] 1.6 12.4 4.6 15.4 See ?par for more information. Note also that par(usr) is not read-only, so you can explicitly change it when required: plot(1:10,1:10) par(usr) [1] 0.64 10.36 0.64 10.36 # Now reset the plot region limits par(usr = c(0, 20, 0, 20)) # Check it par(usr) [1] 0 20 0 20 You do not indicate what OS you are using, but Under Windows, there is an ability to record plots. See R Windows FAQ 4.2. Otherwise, save the R code that you use to generate the plot and either CP it or source() it or use ESS. HTH, Marc Schwartz __ R-help@stat.math.ethz.ch 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] survreg with gamma distribution: re-post
Dear r-help subscribers, A couple of weeks ago I sent the following message to the r-help mail list. It hasn't generated any response, and I could really use some help on this. Anyone able to help? Thanks again, Roger Dungan I am working on some survival analysis of some interval censored failure time data in R. I have done similar analysis before using PROC LIFEREG in SAS. In that instance, a gamma survival function was the optimum parametric model for describing the survival and hazard functions. I would like to be able to use a gamma function in R, but apparently the survival package does not support this distribution. I have been googling around for some help, and have found some threads to a similar question posted to the R-Help list in October last year. Because I am a bit of a survival analysis and R newbie, I didn't really understand the discussion thread. I've been working with a Weibull distribution, thus: leafsurv.weibull-survreg(Surv(minage, maxage, censorcode, type = interval)~1, dist = weib) And I guess I'd like to be able to do something that's the equivalent of leafsurv.gamma-survreg(Surv(minage, maxage, censorcode, type = interval)~1, dist = gamma) At least one of the R-help listserver comments mentioned using survreg.distributions to customise a gamma distribution, but I can't figure out how to make this work with the resources (intellectual and bibliographical!) that I have available. With thanks in advance for your help, Dr Roger Dungan School of Biological Sciences University of Cantebury Christchurch, New Zealand ph +64 3 366 7001 ext. 4848 fax +64 3 354 2590 __ R-help@stat.math.ethz.ch 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] basic question about changing limits on generated plots
On Wed, 2005-02-23 at 17:42 -0500, rif wrote: Is it possible to change the limits on plots that are already on the screen? In particular, is there any R equivalent to the sequence of matlab commands plot(1:10,1:10) hold on plot(2:12,5:15) rif I have not used Matlab, but I suspect that you want to use: par(new = TRUE) For example: plot(1:10,1:10) # Check plot region limits par(usr) [1] 0.64 10.36 0.64 10.36 # Set par(new) to not overwrite existing plot par(new = TRUE) plot(2:12,5:15) Marc Schwartz Marc, This does not do what the matlab code I posted does (the matlab code also works in the free program octave, if you want to try). The matlab code moves already plotted data within the window (replots it). When I first type plot(1:10,1:10), I see a graph with axis limits [1 10 1 10]. When I type hold on (to keep my original data), and execute plot(2:12,5:15), the plot I see is equivalent to the plot I'd have gotten if I'd originally specified axis limits [1 12 5 15]. By contrast, in the R code you sent, it's as if I'm superimposing two unrelated plots. Essentially, the underlying task is that I want to compare multiple functions, but I do not know good limits for the complete set of functions when I start. Being able to adjust the graph to show all the data I've plotted so far would be extremely useful for exploratory analysis. This is the mode I and colleagues generally use matlab and octave in. Does this question get asked all the time? It seems to be something that would come up a lot for people who switch from Matlab/Octave to R, but I searched the archives and didn't really see anything. Cheers, rif __ R-help@stat.math.ethz.ch 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 with html help system: help.start()
I am having problems installing R under Fedora Core 3. I installed R using YUM yum -install R. The basic installation appears to have worked, I can start and run R. When I try to start the HTML help system, I get the errors listed below. Any help solving my problem would be appreciated. Thanks John Sorkin R intalled using YUM under Red Hat Fedora Core 3. Computer: Dell Latitude C400. help.start()Making links in per-session dir .../usr/lib/R/share/sh/help-links.sh: line 9: /usr/lib/R/bin/mkinstalldirs: No such file or directory ln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/manual' to `/usr/lib/R/doc/manual': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/R.css': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/about.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/faq.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/index.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/left.jpg': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/logo.jpg': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/logosm.jpg': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/packages-head.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/packages.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/resources.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/right.jpg': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/template.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/thanks.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html' to `/usr/lib/R/doc/html/up.jpg': No such file or directorycp: cannot create regular file `/tmp/RtmpxI5110/.R/doc/html': No such file or directorycp: cannot create regular file `/tmp/RtmpxI5110/.R/doc/html/search': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html/search' to `/usr/lib/R/doc/html/search/SearchEngine-foot.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html/search' to `/usr/lib/R/doc/html/search/SearchEngine-head.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html/search' to `/usr/lib/R/doc/html/search/SearchEngine.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/doc/html/search' to `/usr/lib/R/doc/html/search/SearchObject.html': No such file or directoryln: creating symbolic link `/tmp/RtmpxI5110/.R/library' to `/usr/lib/R/doc/html/R.css': No such file or directoryIf /usr/bin/firefox is already running, it is *not* restarted, and you must switch to its window.Otherwise, be patient ...Warning message:cannot create HTML package index in: make.packages.html() Error: No running window foundJohn Sorkin M.D., Ph.D. Chief, Biostatistics and Informatics Baltimore VA Medical Center GRECC and University of Maryland School of Medicine Claude Pepper OAICUniversity of Maryland School of Medicine Division of Gerontology Baltimore VA Medical Center 10 North Greene Street GRECC (BT/18/GR) Baltimore, MD 21201-1524410-605-7119 - NOTE NEW EMAIL ADDRESS: [EMAIL PROTECTED] [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Adding up intervals by activity
Hello all, I have a data frame with the following columns: timestamp, activity, interval. Each row contains information about an activity that occured over a time interval. Timestamp is when the interval ended, Activity is a factor that represents the type of activity, and Interval is a double which represents an interval in time. The data frame is sorted by timestamp. For example: TimestampActivity Interval 2004-02-12 08:59:08 A 5 2004-02-12 09:23:03 A 7 2004-02-12 09:56:04 A 4 2004-02-12 10:39:00 B 5 2004-02-12 11:23:06 B 3 2004-02-12 12:56:01 A 4 2004-02-12 12:59:01 A 2 ... I would like to know how much time was spent on an activity before the activity changed (to compute time, I just care about the lengths of the intervals, not the values of the timestamps). For example, given the table above, I would like to get an output that looks like: Activity TotalTime A 16 B 8 A 6 Is this possible to do in R without resorting to a for loop or recursion (i.e. what's the R-ish way to do this?) I know I can use the diff function to identify which rows represent a change in activity, but is there anything I can do with that? Lorin -- Lorin Hochstein Graduate Research Assistant Experimental Software Engineering Group Computer Science Dept., University of Maryland email:[EMAIL PROTECTED] tel:301-405-2721 __ R-help@stat.math.ethz.ch 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] Adding up intervals by activity
Lorin
You could use rle(): Say your Activity and Interval data is in lorin, then:
tmp.rle - rle(as.vector(lorin[,1]))[[1]]
tapply(lorin[,2], rep(1:length(tmp.rle), tmp.rle), sum)
1 2 3
16 8 6
HTH
Peter Alspach
Lorin Hochstein [EMAIL PROTECTED] 24/02/05 16:18:43
Hello all,
I have a data frame with the following columns: timestamp, activity,
interval.
Each row contains information about an activity that occured over a time
interval. Timestamp is when the interval ended, Activity is a factor
that represents the type of activity, and Interval is a double which
represents an interval in time. The data frame is sorted by timestamp.
For example:
TimestampActivity Interval
2004-02-12 08:59:08 A 5
2004-02-12 09:23:03 A 7
2004-02-12 09:56:04 A 4
2004-02-12 10:39:00 B 5
2004-02-12 11:23:06 B 3
2004-02-12 12:56:01 A 4
2004-02-12 12:59:01 A 2
...
I would like to know how much time was spent on an activity before the
activity changed (to compute time, I just care about the lengths of the
intervals, not the values of the timestamps). For example, given the
table above, I would like to get an output that looks like:
Activity TotalTime
A 16
B 8
A 6
Is this possible to do in R without resorting to a for loop or recursion
(i.e. what's the R-ish way to do this?) I know I can use the diff
function to identify which rows represent a change in activity, but is
there anything I can do with that?
Lorin
--
Lorin Hochstein
Graduate Research Assistant
Experimental Software Engineering Group
Computer Science Dept., University of Maryland
email:[EMAIL PROTECTED] tel:301-405-2721
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
The contents of this e-mail are privileged and/or confidenti...{{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] a simple question
Dear All, I need to separate one column which is in a format like, can_region, into two columns like, can, and region. Do we have a function in R to do this? Thanks. Regards, Jin [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Row median of Date class variables in a data frame
I am trying to calculate the median of each row of a
data frame where the data frame consist of
columns of class Date.
Below are my test data and best attempt at using apply.
I didn't see a solution via Google or the Baron search
site.
I'd be grateful for any suggestions or solutions.
I'm using R 2.0.0 on Mac OS X.
Thank you,
Stephen Weigand
### Test data
date1 - c(1000, 2000, 3000,4000)
date2 - date1 + 100
date3 - date2 + 100
class(date1) - class(date2) - class(date3) - Date
test - data.frame(date1, date2, date3)
print(test)
### create a function for apply()
medDate - function(x){
obj - unclass(unlist(x))
med - median(obj, na.rm = TRUE)
med
class(med) - Date
med
}
medDate(test$date1) # works
medDate(test[1,]) # works
apply(test, 1, medDate) # gives error: 'need numeric data'
apply(test, 2, medDate) # gives error: 'need numeric data'
__
R-help@stat.math.ethz.ch 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] survreg with gamma distribution: re-post
Google for the gfcure package on statlib. I have not used it myself and am not entirely sure whether it works on R or just S. Roger Dungan roger.dungan at canterbury.ac.nz writes: : : Dear r-help subscribers, : : A couple of weeks ago I sent the following message to the r-help mail : list. It hasn't generated any response, and I could really use some help : on this. Anyone able to help? : : Thanks again, : : Roger Dungan : : : I am working on some survival analysis of some interval censored failure : time data in R. I have done similar analysis before using PROC LIFEREG : in SAS. In that instance, a gamma survival function was the optimum : parametric model for describing the survival and hazard functions. I : would like to be able to use a gamma function in R, but apparently the : survival package does not support this distribution. I have been : googling around for some help, and have found some threads to a similar : question posted to the R-Help list in October last year. Because I am a : bit of a survival analysis and R newbie, I didn't really understand the : discussion thread. : : I've been working with a Weibull distribution, thus: : : leafsurv.weibull-survreg(Surv(minage, maxage, censorcode, type = : interval)~1, dist = weib) : : And I guess I'd like to be able to do something that's the equivalent of : : leafsurv.gamma-survreg(Surv(minage, maxage, censorcode, type = : interval)~1, dist = gamma) : : At least one of the R-help listserver comments mentioned using : survreg.distributions to customise a gamma distribution, but I can't : figure out how to make this work with the resources (intellectual and : bibliographical!) that I have available. : : With thanks in advance for your help, : : Dr Roger Dungan : School of Biological Sciences : University of Cantebury : Christchurch, New Zealand : ph +64 3 366 7001 ext. 4848 : fax +64 3 354 2590 : : __ : R-help at stat.math.ethz.ch 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] a simple question
Does strsplit fit what you want? Tom -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Sent: Thursday, 24 February 2005 11:46 AM To: r-help@stat.math.ethz.ch Subject: [R] a simple question Dear All, I need to separate one column which is in a format like, can_region, into two columns like, can, and region. Do we have a function in R to do this? Thanks. Regards, Jin [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how do I get means by factor?
Chris grover2160 at metrocast.net writes: : : I want the mean of a response variable for each level of a factor in a : data.frame. How is this done? In addition to the tapply suggestion already given you can use summary in package Hmisc: summary(y ~ f) __ R-help@stat.math.ethz.ch 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] baseline shift / detrend
Peter Soros soros at mac.com writes: Hello, I am analyzing data of biological motion (mouth movements) over time (~ 10 mins). There is a considerable baseline shift which appears to be a linear trend. Is it possible to get rid of this baseline shift using R? Any help is appreciated! The following will detrend y: resid(lm(y ~ seq(y))) __ R-help@stat.math.ethz.ch 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] Row median of Date class variables in a data frame
Stephen D. Weigand weigand.stephen at charter.net writes:
:
: I am trying to calculate the median of each row of a
: data frame where the data frame consist of
: columns of class Date.
:
: Below are my test data and best attempt at using apply.
: I didn't see a solution via Google or the Baron search
: site.
:
: I'd be grateful for any suggestions or solutions.
: I'm using R 2.0.0 on Mac OS X.
:
: Thank you,
:
: Stephen Weigand
:
: ### Test data
:
: date1 - c(1000, 2000, 3000,4000)
: date2 - date1 + 100
: date3 - date2 + 100
:
: class(date1) - class(date2) - class(date3) - Date
:
: test - data.frame(date1, date2, date3)
:
: print(test)
:
: ### create a function for apply()
: medDate - function(x){
:obj - unclass(unlist(x))
:med - median(obj, na.rm = TRUE)
:med
:class(med) - Date
:med
: }
:
: medDate(test$date1) # works
: medDate(test[1,]) # works
:
: apply(test, 1, medDate) # gives error: 'need numeric data'
: apply(test, 2, medDate) # gives error: 'need numeric data'
:
Try this:
library(zoo)
as.Date(apply(data.matrix(test), 1, median)) # ignore warning
The only reason you need zoo is that it supplies as.Date methods
which convert numerics and integers to Date.
__
R-help@stat.math.ethz.ch 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] basic question about changing limits on generated plots
On Wed, 2005-02-23 at 21:14 -0500, rif wrote: snip Marc, This does not do what the matlab code I posted does (the matlab code also works in the free program octave, if you want to try). The matlab code moves already plotted data within the window (replots it). When I first type plot(1:10,1:10), I see a graph with axis limits [1 10 1 10]. When I type hold on (to keep my original data), and execute plot(2:12,5:15), the plot I see is equivalent to the plot I'd have gotten if I'd originally specified axis limits [1 12 5 15]. By contrast, in the R code you sent, it's as if I'm superimposing two unrelated plots. Essentially, the underlying task is that I want to compare multiple functions, but I do not know good limits for the complete set of functions when I start. Being able to adjust the graph to show all the data I've plotted so far would be extremely useful for exploratory analysis. This is the mode I and colleagues generally use matlab and octave in. Does this question get asked all the time? It seems to be something that would come up a lot for people who switch from Matlab/Octave to R, but I searched the archives and didn't really see anything. Cheers, rif A general statement: There are members of the R Community engaged in Octave, so there is some overlap, at least in terms of expertise with both tools. Perhaps they can offer some insight here. The good news is that I have Octave installed on my FC3 system, so I was able to get a feel for what you are referring to. A search of the R archives would suggest that there is not a direct parallel in terms of adding a new set of data to an existing plot, while having the entire coordinate system adjusted to the new data in a single step. There are references to what I suggested, the use of the 'add' argument in some plot functions and of course the use of points() and lines(). I would defer to others with more low-level knowledge of the standard R screen based plotting devices, but from my past review of code (both R and C) and reasonable knowledge, I am not sure that this can be done without some form of two step approach involving re-plotting the original data using an updated coordinate system based upon the new data and then overlaying the new data. Presumably at a low level, this is what Octave and Matlab do, since I noted that the plot device seems to be completely re-drawn upon the second plot call. Since Octave is open source, one can of course review the code to see what is truly happening. The key it seems would be to transparently save the original data as an object, re-plot it with the adjusted coordinate system and then add the new data. I would guess that with some thought, one could create some wrapper plotting functions or methods that would save the data object(s) in a plot environment so that each successive plot call layers the older data sets in turn and then add the newest data set. All this done in a coordinate system that is re-configured each time, based upon the maximum x and y axis ranges required for the multiple datasets that have been plotted to that point. Almost sounds like a plot stack, to use an assembly language metaphor. You might also want to look at the xgobi package on CRAN, which provides an interface to the XGobi system: http://www.research.att.com/areas/stat/xgobi/ or the ggobi system: http://www.ggobi.org/ which is accessible from R. These are both dynamic data visualization tools. HTH, Marc __ R-help@stat.math.ethz.ch 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] KalmanXXXX and deJong-Penzer statistic?
A question about: Kalman in R, time series and deJong-Penzer statistic - how to compute it using available artefacts of KalmanX? Background. in the paper http://www.lse.ac.uk/collections/statistics/documents/researchreport34.pdf 'Diagnosing Shocks in TIme Series', de Jong and Penzer construct a statistic (tau) which can be used to locate potential shocks. [p15, Theorem 6.1 and below]. They also state that all the components of that statistic (v_i, F_i, r_i, N_i) 'are computed with Kalman Filter Smoother applied to the null model'. Also, as I understand, that part has been implemented in one of the S packages , SsfPack, as the book on that states on p 531 'the standardized smoothed disturbances may be interpreted as t-statistics for impulse intervention variable in the transition and measurement equations.' and equations for the statistic are: eta_t / sqrt(var(eta_t), where eta_t and var have hats over them. The second equation is identical, with 'eta' replaced by 'epsilon'. On page 524 we also have: the smoothed disturbance estimates are the estimates of the measurement equation innovations epsilon and transition equation innovations eta based on all available information Y. ... the computation of hat(eta) and hat(epsilon) from the Kalman smoother algorithm is described in Durbin and Koopman chapter 7, 'Time series analysis by state space methods', OUP (2001) Local libraries do not have this book and it will take several weeks to get it. Assuming I will get the book: does the KalmanXXX set of functions produce all the necessary artefacts to compute this statistic either as per deJong-Penzer or as per SsfPack? Reading carefully through the manual I see that we have artefacts of states and normalized residuals (presumably of states - but how can I unscale them if I need them)? What about other stats? How to compute smoothed disturbance estimates? I am rather confused - that's my first approach to Kalman filter and state models. __ R-help@stat.math.ethz.ch 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] Row median of Date class variables in a data frame
Note that the argument to apply() is a matrix, not a data frame, so your
object has been coerced to a matrix. That is spelled out on the help page
for apply(), as it is frequent misuse.
You would be better off using a matrix of class Date in the first place if
you want row operations.
On Wed, 23 Feb 2005, Stephen D. Weigand wrote:
I am trying to calculate the median of each row of a
data frame where the data frame consist of
columns of class Date.
Below are my test data and best attempt at using apply.
I didn't see a solution via Google or the Baron search
site.
I'd be grateful for any suggestions or solutions.
I'm using R 2.0.0 on Mac OS X.
Thank you,
Stephen Weigand
### Test data
date1 - c(1000, 2000, 3000,4000)
date2 - date1 + 100
date3 - date2 + 100
class(date1) - class(date2) - class(date3) - Date
test - data.frame(date1, date2, date3)
print(test)
### create a function for apply()
medDate - function(x){
obj - unclass(unlist(x))
med - median(obj, na.rm = TRUE)
med
class(med) - Date
med
}
medDate(test$date1) # works
medDate(test[1,]) # works
apply(test, 1, medDate) # gives error: 'need numeric data'
apply(test, 2, medDate) # gives error: 'need numeric data'
--
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
[R] a question about function eval()
Hi,
I have a question about the usage of eval(). Wonder if any experienced user can
help me out of it.
I use eval() in the following function:
semireg.pwl -
function(coef.s=rnorm(1),coef.a=rnorm(1),knots.pos=knots.x,knots.ini.val=knots.val){
knotn - length(knots.pos)
def.par.env - sys.frame(1)
print(def.par.env)
print(environment(coef.s))
tg - eval( (parse(text=
paste(
function(coef.sex=coef.s,coef.age=coef.a,,
paste(knot.val.,1:knotn,=knots.ini.val[,1:knotn,],sep=,collapse=,),
),
{\n print(sys.frame());print(coef.sex)\n y - c(,
paste(knot.val.,1:knotn,sep=,collapse=,),
)\n,
-loglikelihood(coef.sex,coef.age,knots.pos,y)\n},
sep=,collapse=
)
)
)
,envir=def.par.env)
print(tg())
print(coef.s)
print(sys.frame(1))
mle(tg)
}
But when I ran semireg.pwl(), I got correct value for tg(), but an message
Error in eval(expr, envir, enclos) : Object coef.s not found
for mle(tg). I just don't know how to make the environment variable visible for
mle().
Thanks a lot!
Feng
[[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
