Re: [R] Underline in expression().
Sundar, Thanks. Unfortunately, I am looking for something that also works in the margins of the plot. John. Sundar Dorai-Raj wrote: John Janmaat wrote: Hello All, Is there an analogue to \underbar or the AMS math \underline in graphical math expressions? Thanks, John. Uwe Ligges posted a solution a couple of years ago. I don't know if there is anything built in yet. ?plotmath does not seem to say anything about underlining. http://finzi.psych.upenn.edu/R/Rhelp01/archive/7191.html plot(0:1, 0:1, type=n) underlined(0.5, 0.5, expression(widehat(x %*% y))) --sundar __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- = Dr. John Janmaat Department of Economics Acadia University Wolfville, Nova Scotia, Canada B4P 2R6 TEL: 902-585-1461 WWW: http://ace.acadiau.ca/~jjanmaat/ EMAIL: [EMAIL PROTECTED] June 10, 2004 to September 7, 2004 Dr. John Janmaat Environmental Econmics and Natural Resources Group Wageningen University. P.O. Box 8130 6700 EW, Wageningen The Netherlands. __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Best way to store negative indexes
StephaneDemurget [EMAIL PROTECTED] writes: Hi, I'm trying to figure out how to properly construct a graph of frequencies of a difference between 2 values, that is | i | - | j |. I'd like to know the best way to store the actual data because of course doing my_vector[i -j] will not work because the index could be negative. I know there's no hash table so what's the best solution, simplicity-wise ? Use a list of pair of values {index, i - j} ? Or can I somehow use negative index, perhaps using 0array ? Any help would be appreciated :) (i and j are integer vectors, right?) How about this: f - i - j f - factor(f,levels=min(f):max(f)) table(f) -- O__ Peter Dalgaard Blegdamsvej 3 c/ /'_ --- Dept. of Biostatistics 2200 Cph. N (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Another big data size problem
Federico Gherardini wrote: Hi all, I'm trying to read a 1220 * 2 table in R but I'm having lot of problems. Basically what it happens is that R.bin starts eating all my memory until it gets about 90%. At that point it locks itself in a uninterruptible sleep status (at least that's what top says) where it just sits there barely using the cpu at all but keeping its tons of memory. I've tried with read.table and scan but none of them did the trick. I've also tried some orrible hack like reading one line a time and gradually combining everything in a matrix using rbind... nope! It seems I can read up to 500 lines in a *decent* time but nothing more. The machine is a 3 GHz P4 with HT and 512 MB RAM running R-1.8.1. Will I have to write a little a C program myself to handle this thing or am I missing something? If your data is numeric, you will need roughly 1220 * 2 * 8 / 1024 / 1024 ~~ 200 MB just to store one copy in memory. If you need more than two copies, your machine with its 512MB will start to use swap space . Hence either use a machine with more memory, or don't use all the data at once in memory, e.g. by making use of a database. Uwe Ligges Thanks in advance for your help, fede __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Integration with adapt
Rodrigo == Rodrigo Drummond [EMAIL PROTECTED] on Tue, 27 Jul 2004 15:31:07 -0300 (BRT) writes: Rodrigo Hi all, I need to calculate a multidimensional Rodrigo integration on R. I am using the command Rodrigo adapt (from library adapt), although it's a package, not a library. Rodrigo sometimes I get the following error message: Rodrigo Ifail=2, lenwrk was too small. -- fix adapt() ! Rodrigo Check the returned relerr! in: adapt(3, linf, lsup, Rodrigo functn = Integrando1) If you could give as a *reproducible* example, we (the adapt authors) would have chance to do what the above message says, namely fix adapt() .. Rodrigo I guess it happens because the domain of Rodrigo integration is too small, maybe, maybe not. We need an example we can reproduce, see above. Rodrigo although I tried a change of variables to avoid Rodrigo this problem and it didnt help. The command adapt Rodrigo calls a fortran routine, but I dont know fortran Rodrigo enough to fix the problem. Martin Maechler, ETH Zurich __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] problem with the .Rhistory
Hello, The history function doesn't seems to work when I am working with R. I have a .Rhistory file in my working directory and I have tried savehistory(file = .Rhistory) or the command given in the history help : .Last - function() if(interactive()) try(savehistory(.Rhistory)) and I get the same error message Error in savehistory(file) : no history available to save Is there any other variable to set up? Nathalie Peyrard __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] problem with the .Rhistory
It depends on your platform. E.g. on Unix/Linux, ?savehistory says Details: This works under the 'readline' and GNOME interfaces, but not if 'readline' is not available (for example, in batch use). To be able to help you further, we would need quite precise details of your platform and which interface you are using. On Wed, 28 Jul 2004, Peyrard Nathalie wrote: The history function doesn't seems to work when I am working with R. I have a .Rhistory file in my working directory and I have tried savehistory(file = .Rhistory) or the command given in the history help : .Last - function() if(interactive()) try(savehistory(.Rhistory)) and I get the same error message Error in savehistory(file) : no history available to save -- 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 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Latex error about Schunk environment under Windows NT
Hi, I am running the Sweave function on a Rnw file. It produces a tex file with \begin{Schunk} \begin{Sinput} \end{Sinput} \end{Schunk} and when I try to compile the tex file, I get a Latex error : Environment Schunk undefined However, I wrote \usepackage{Sweave} at the beginning of my Rnw file And the file Sweave.sty does exist in my tex files. I am not experiencing the same trouble with my personal computer at home. However I did the same steps. Is it due to Windows NT (my home computer runs under XP) and shall I install another latex package? Anything would help. Thanks a lot in advance, Matthieu Cornec __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Another big data size problem
On Wed, 2004-07-28 at 03:10, Federico Gherardini wrote: Hi all, I'm trying to read a 1220 * 2 table in R but I'm having lot of problems. Basically what it happens is that R.bin starts eating all my memory until it gets about 90%. At that point it locks itself in a uninterruptible sleep status (at least that's what top says) where it just sits there barely using the cpu at all but keeping its tons of memory. I've tried with read.table and scan but none of them did the trick. I've also tried some orrible hack like reading one line a time and gradually combining everything in a matrix using rbind... nope! It seems I can read up to 500 lines in a *decent* time but nothing more. The machine is a 3 GHz P4 with HT and 512 MB RAM running R-1.8.1. Will I have to write a little a C program myself to handle this thing or am I missing something? Thanks in advance for your help, fede Hi, It looks like you're running linux !? if so it will be quite easy to create a table in MySQL, upload all the data into the database and access the data with RMySQL (it's _very_ fast). Probably there will be some operations that you can do on MySQL instead of eating memory in R. Regards EJ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Latex error about Schunk environment under Windows NT
Cornec Matthieu wrote: Hi, I am running the Sweave function on a Rnw file. It produces a tex file with \begin{Schunk} \begin{Sinput} \end{Sinput} \end{Schunk} and when I try to compile the tex file, I get a Latex error : Environment Schunk undefined However, I wrote \usepackage{Sweave} at the beginning of my Rnw file And the file Sweave.sty does exist in my tex files. I am not experiencing the same trouble with my personal computer at home. However I did the same steps. Is it due to Windows NT (my home computer runs under XP) and shall I install another latex package? Anything would help. Thanks a lot in advance, Matthieu Cornec __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html It is required to have the style file, of course. Calling Swaeve() from R should work (but there are some known MikTeX related problems) anyway. You can find the style file in folder .../share/texmf Uwe Ligges __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Fw: [R] Another big data size problem
On Wed, 28 Jul 2004 09:53:08 +0200 Uwe Ligges [EMAIL PROTECTED] wrote: If your data is numeric, you will need roughly 1220 * 2 * 8 / 1024 / 1024 ~~ 200 MB just to store one copy in memory. If you need more than two copies, your machine with its 512MB will start to use swap space . Hence either use a machine with more memory, or don't use all the data at once in memory, e.g. by making use of a database. Uwe Ligges Well I'd be happy if it used swap space instead of locking itself up! By the way I don't think that the problem is entirely related to memory consumption. I have written a little function that reads the data row by row and does a print each time, to monitor its functioning. Everything starts to crwal to an horrible slowness long before my memory is exhausted... i.e.: after about 100 lines. It seems like R has problems managing very large objects per se? By the way I'll try to upgrade to 1.9 and see what happens... Ernesto Jardim wrote: Hi, It looks like you're running linux !? if so it will be quite easy to create a table in MySQL, upload all the data into the database and access the data with RMySQL (it's _very_ fast). Probably there will be some operations that you can do on MySQL instead of eating memory in R. Regards EJ I'll give that a try. Thanks everybody for their time fede __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: Fw: [R] Another big data size problem
Federico Gherardini wrote: On Wed, 28 Jul 2004 09:53:08 +0200 Uwe Ligges [EMAIL PROTECTED] wrote: If your data is numeric, you will need roughly 1220 * 2 * 8 / 1024 / 1024 ~~ 200 MB just to store one copy in memory. If you need more than two copies, your machine with its 512MB will start to use swap space . Hence either use a machine with more memory, or don't use all the data at once in memory, e.g. by making use of a database. Uwe Ligges Well I'd be happy if it used swap space instead of locking itself up! By the way I don't think that the problem is entirely related to memory consumption. I have written a little function that reads the data row by row and does a print each time, to monitor its functioning. Everything starts to crwal to an horrible slowness long before my memory is exhausted... i.e.: after about 100 lines. It seems like R has problems managing very large objects per se? By the way I'll try to upgrade to 1.9 and see what happens... Well, using swap space takes much time. And what looks like hanging up is quite probably the use of swap space - you will see your hard disc LED flashing all the time! Are you sure that your memory was not exhausted? Note that it is better to initialize the object to full size before inserting -- rather than using rbind() and friends which is indeed slow since it need to re-allocate much memory for each step. Uwe Ernesto Jardim wrote: Hi, It looks like you're running linux !? if so it will be quite easy to create a table in MySQL, upload all the data into the database and access the data with RMySQL (it's _very_ fast). Probably there will be some operations that you can do on MySQL instead of eating memory in R. Regards EJ I'll give that a try. Thanks everybody for their time fede __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] problem with the .Rhistory
I recall having a similar problem a while back, and it turned out that it was due to command line editing not being available. This in turn was due to a problem with the readline package. I think we had a flakey version (obtained as a Sun package?) and when we got another version (from GNU?) and installed it --- and then rebuilt R --- command line editing (and savehistory()) worked fine. Sorry to be vague, but I don't have good records of what went on, and it was a couple of years ago. Check on command line editing first --- in R execute capabilities() and if ``cledit'' turns out to be FALSE, then readline is probably your problem. Peter Dalgaard may be able to give you better insight as to how to fix it. cheers, Rolf Turner [EMAIL PROTECTED] ===+===+===+===+===+===+===+===+===+===+===+===+===+===+===+===+===+===+=== Nathalie.Peyrard wrote: The history function doesn't seems to work when I am working with R. I have a .Rhistory file in my working directory and I have tried savehistory(file = .Rhistory) or the command given in the history help : .Last - function() if(interactive()) try(savehistory(.Rhistory)) and I get the same error message Error in savehistory(file) : no history available to save Is there any other variable to set up? __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Another big data size problem
Hi, i'm working with a ~ 250.000 * 150 data.frame and can share your problems - i've upgraded last weekend my notebook from 512MB - 1024MB, it's really better especially for load, write.table , mysqlReadTable, mysqlWriteTable, because machine begin caching if RAM is full. One example: With 512MB i get after some hours no success write a table to mysql. With 1024MB it does in some minutes. regards, christian Am Mittwoch, 28. Juli 2004 04:10 schrieb Federico Gherardini: Hi all, I'm trying to read a 1220 * 2 table in R but I'm having lot of problems. Basically what it happens is that R.bin starts eating all my memory until it gets about 90%. At that point it locks itself in a uninterruptible sleep status (at least that's what top says) where it just sits there barely using the cpu at all but keeping its tons of memory. I've tried with read.table and scan but none of them did the trick. I've also tried some orrible hack like reading one line a time and gradually combining everything in a matrix using rbind... nope! It seems I can read up to 500 lines in a *decent* time but nothing more. The machine is a 3 GHz P4 with HT and 512 MB RAM running R-1.8.1. Will I have to write a little a C program myself to handle this thing or am I missing something? Thanks in advance for your help, fede __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Another big data size problem
On Wed, 2004-07-28 at 12:40, Christian Schulz wrote: Hi, i'm working with a ~ 250.000 * 150 data.frame and can share your problems - i've upgraded last weekend my notebook from 512MB - 1024MB, it's really better especially for load, write.table , mysqlReadTable, mysqlWriteTable, because machine begin caching if RAM is full. One example: With 512MB i get after some hours no success write a table to mysql. With 1024MB it does in some minutes. Hi, When you're writing a table to MySQL you have to be carefull if the table is created by RMySQL. The fields definition may not be the most adequate and there will be no indexes in your table, which makes the queries _very_ slow. Regards EJ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Re: group definition for a bootstrap
Hi, There dose not seem to be any other replies so far, may be because the question is not quite clear. Do you want a separate estimation of the 'precision of the mean' for each age by quarter subsample? This reply is based on that interpretation. If not, you may need some sort of stratified bootstrap (see the article by Canty mentioned in the code below). Also, since you want to use the bootstrap for the mean, are the (age by quarter) subsample sizes small, or the subsamples heavily skewed? In any case it is better to use the package boot by Canty and Ripley than programming from scratch. A good introduction to this package is written by Angelo J. Canty in R Newsletter 2/3, December 2002. The following code first generates an aretificial (large) data frame dd with the same structure as your fram d. Then it performs bootstrap calculations in each age by quarter subsample. The output includes for each subsample bootstrap based bias and standard error for the mean, plots for checking discreteness problems of the bootstrap and normal and BCa confidence intervals. To use this code for your sample it should suffice to set dd - d and start at the line # prepares bootstrap procedure for the mean in the code. The code has been checked for the large artificial sample generated. It will probably break down if some subsamples are empty or nearly so. And there may be discreteness problems if some subsamples are small, although the bootstrap is fairly stable for means. Best, Tore Wentzel-Larsen # 'Group definition for a bootstrap' # generates artificial data frame dd to use (age and quarter as factors in case they are in your data): nn - floor(runif(1,1000,2000)) # random sample size dd - data.frame(factor(floor(runif(nn,0,6))), factor(floor(runif(nn,1,5))), rnorm(nn,15,3)) names(dd) - c('age','quarter','x') # prepares bootstrap procedure for the mean # (cf the article by Angelo J. Canty in R Newsletter 2/3, December 2002): library(boot) mean.w - function(xx,ww) sum(xx*ww) # means, with bootstrap estimates of bias and standard error, # for each age by quarter subsample: nboot - 1000 # number of bootstrap replications par(ask = TRUE) # prompts before next plot for (age in 0:5) { for (quarter in 1:4) { dd.sub - dd[dd$age==agedd$quarter==quarter,] # the subsample boot.subsample - boot(data=dd.sub$x, statistic=mean.w, R=nboot, stype='w') print(boot.subsample) print(c('age:', age, 'quarter:', quarter)) # for keeping track of the subsamples plot(boot.subsample,) # checks for discreteness problems # bootstrap confidence intervals (normal based and BCa): ci.subsample - boot.ci(boot.subsample, type=c('norm','bca'), conf=c(.95)) print(ci.subsample) } } par(ask = FALSE) # back to default setting for ask Hi, This is probably really simple, but I am clearly not R-minded, I have read the help files, and reread them, and I still can't work out what to do... I have a data frame (d) with 3 columns (age (0-5), quarter (1-4) and x). I want to estimate the precision of my mean x by age and quarter, so I want to carry out a bootstrap for each group. I am trying to do this within a loop, so I don't have to retype the whole thing out 20 times ... This is what I have done: n = length (d$x) N = 1000 stat = numeric(N) for (i in 1:N) { d$x2 = sample (d$x, n, replace=T) stat[i] = mean(d$x2) } I believe I should define the age and quarter groups in the line straight after for - using the split function, or should I use some variant of [d$age == 1 d$quarter ==1] in the sample definition? Please don't think I am looking for an easy answer - I have been puzzling for this for over a week already :( (I am using R1.9 under MS W2000) Thank you in advance Louize __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Simulation from a model fitted by survreg.
Dear list, I would like to simulate individual survival times from a model that has been fitted using the survreg procedure (library survival). Output shown below. My plan is to extract the shape and scale arguments for use with rweibull() since my error terms are assumed to be Weibull, but it does not make any sense. The mean survival time is easy to predict, but I would like to simulate individual survival times. I am probably missing something completely obvious. Any hints or advice are appreciated. Thanks Sixten summary(mod1) Call: survreg(formula = Surv(tid, study$first.event.death) ~ regim + age + stadium2, data = study, dist = weibull) Value Std. Error zp (Intercept) 11.6005 0.7539 15.387 2.01e-53 regimposto -0.1350 0.1558 -0.867 3.86e-01 age -0.0362 0.0102 -3.533 4.11e-04 stadium2ii -0.0526 0.2794 -0.188 8.51e-01 Log(scale) -0.5148 0.1116 -4.615 3.93e-06 Scale= 0.598 Weibull distribution Loglik(model)= -680.7 Loglik(intercept only)= -689.2 Chisq= 16.87 on 3 degrees of freedom, p= 0.00075 Number of Newton-Raphson Iterations: 8 n=1183 (4 observations deleted due to missing) version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major1 minor8.1 year 2003 month11 day 21 language R __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Writing polyline shapefiles
Stephane, dear R listers, I wonder if the idea of developping some functions in R writing shapefiles from a matrix of coordinates -eg poly2shape()- is still with us? I currently manage with a home-made function writing an ASCII grass file, then importing it into GRASS as vector file, then exporting the vector file as a shapefile (ouf!). The main issue to me is still to be capable to wite the *.shx and *.shp from coordinates (the *.dbf is easier to manage, even indirectly, and can for instance be completed simply via Excel). I tried to open *.shp, etc.. files to understand how they were structured and thus write some programme by myself, but unsuccessfully (the way those files are coded is not ascii) Maptools is fantastic for importation and data handling within R, but unfortunately cannot export to shapefiles. Any news? Patrick Giraudoux PS: ESRI has managed to distribute ArcTools (with ArcGIS) with no way to importe generate ASCII formats!!! One must get ArcINFO for that... No comment! __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Another big data size problem
On Wed, 28 Jul 2004 13:28:20 +0100 Ernesto Jardim [EMAIL PROTECTED] wrote: Hi, When you're writing a table to MySQL you have to be carefull if the table is created by RMySQL. The fields definition may not be the most adequate and there will be no indexes in your table, which makes the queries _very_ slow. So, if I understood correctly, if you want to use SQL you'll have to upload the table in SQL, directly from MySQL without using R at all, and then use RMySQL to read the elements in R? Uwe Ligges [EMAIL PROTECTED] wrote: Note that it is better to initialize the object to full size before inserting -- rather than using rbind() and friends which is indeed slow since it need to re-allocate much memory for each step. Do you mean something like this? tab - matrix(rep(0, 1227 * 2), 1227, 2, byrow = TRUE) for(i in 0:num.lines) tab[i + 1,] - scan(mytab, nlines = 1, what=PS, skip = i) The above doesn't get very far either... it seems that, once it has created the table, it becomes so slow that it's unusable. I'll have to try this with more RAM by the way. Cheers, fede __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] elegant matrix creation
Hello everybody. I am trying to reproduce a particular matrix in an elegant way. If I have jj1 - structure(c(1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2, 3,1,2,3,1,2,3,1,2,3,2,3,1,2,3,1,2,3,1,2,3, 1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,3,1,2,3,1, 2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1, 2),.Dim = as.integer(c(9,9))) [ image(jj1) is good here ] then I can get this with kronecker(matrix(1,3,1),kronecker(1+outer(0:2,0:2,+)%%3,matrix(1,1,3))) I want to reproduce the following matrices in an equivalent way: jj2 - matrix(c(1,2,3,1,2,3,1,2,3,2,3,1,2,3,1,2,3,1, 1,2,3,1,2,3,1,2,3,3,1,2,3,1,2,3,1,2,1,2,3,1,2, 3,1,2,3,3,1,2,3,1,2,3,1,2,2,3,1,2,3,1,2,3,1,3, 1,2,3,1,2,3,1,2,2,3,1,2,3,1,2,3,1),9,9) jj3 - structure(c(1,2,3,2,3,1,3,1,2,1,2,1,2,3,2,3,1, 3,1,3,1,2,1,2,3,2,3,2,3,1,3,1,2,1,2,3,2,3, 2,3,1,3,1,2,1,2,1,2,3,2,3,1,3,1,3,1,2,1,2, 3,2,3,1,3,1,3,1,2,1,2,3,2,3,2,3,1,3,1,2,1, 2),.Dim = as.integer(c(9,9))) [ note that jj1-jj3 each have precisely 3 occurrences of A, B, and C along each row, column and (broken) diagonal ]. Can anyone give me a nice elegant way of creating jj2 and jj3 please? -- Robin Hankin Uncertainty Analyst Southampton Oceanography Centre SO14 3ZH tel +44(0)23-8059-7743 [EMAIL PROTECTED] (edit in obvious way; spam precaution) [[alternative HTML version deleted]] __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Another big data size problem
On Wed, 2004-07-28 at 16:26, Federico Gherardini wrote: On Wed, 28 Jul 2004 13:28:20 +0100 Ernesto Jardim [EMAIL PROTECTED] wrote: Hi, When you're writing a table to MySQL you have to be carefull if the table is created by RMySQL. The fields definition may not be the most adequate and there will be no indexes in your table, which makes the queries _very_ slow. So, if I understood correctly, if you want to use SQL you'll have to upload the table in SQL, directly from MySQL without using R at all, and then use RMySQL to read the elements in R? If you can't get R to read the table then you can't use R/RMySQL ;) When creating tables in MySQL you may improve a lot the speed of scaning the table with the proper fields definition and indexing. See the MySQL manual, it is a very good reference and quite easy to read (at least the sections related with table design). The basis is to create indexes for the fields you're supose to restrict (fields to be used with the WHERE clause in SQL statements) considering the order of those fields. Regards EJ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Simulation from a model fitted by survreg.
Please check the documentation, e.g., Venables and Ripley (2002) Modern Applied Statistics with S, p. 360. The Weibull shape parameter is the reciprocal of the scale parameter 0.598 in your printout, so shape = 1/0.598 = 1.672; see also Meeker Escobar (1998) Statistical Methods for Reliability Data (Wiley). Does this answer the question? You can get coefficients with coef(mod1). Also, have you looked at attributes(summary(mod1))? If mod1 follows the old S3 standard, attributes may give you a list of names you can access via summary(mod1)$whateverthenameis (or via summary(mod1)[[whateverthenameis]]). If mod1 follows the new S4 standard, then getSlots(mod1) and getSlots(summary(mod1)) will give you the names of the slots and their classes, which can then be accessed via [EMAIL PROTECTED] Sorry, I don't have time to construct an example myself, but I've done this kind of thing many times. hope this helps. spencer graves Sixten Borg wrote: Dear list, I would like to simulate individual survival times from a model that has been fitted using the survreg procedure (library survival). Output shown below. My plan is to extract the shape and scale arguments for use with rweibull() since my error terms are assumed to be Weibull, but it does not make any sense. The mean survival time is easy to predict, but I would like to simulate individual survival times. I am probably missing something completely obvious. Any hints or advice are appreciated. Thanks Sixten summary(mod1) Call: survreg(formula = Surv(tid, study$first.event.death) ~ regim + age + stadium2, data = study, dist = weibull) Value Std. Error zp (Intercept) 11.6005 0.7539 15.387 2.01e-53 regimposto -0.1350 0.1558 -0.867 3.86e-01 age -0.0362 0.0102 -3.533 4.11e-04 stadium2ii -0.0526 0.2794 -0.188 8.51e-01 Log(scale) -0.5148 0.1116 -4.615 3.93e-06 Scale= 0.598 Weibull distribution Loglik(model)= -680.7 Loglik(intercept only)= -689.2 Chisq= 16.87 on 3 degrees of freedom, p= 0.00075 Number of Newton-Raphson Iterations: 8 n=1183 (4 observations deleted due to missing) version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major1 minor8.1 year 2003 month11 day 21 language R __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] using Rterm under cygwin, no possiblity to delete characters
Dear R-users, When I call Rterm from cygwin, I have no options but typing the exact syntax the first time. If I happen to hit the delete key (backspace), R dies when I press enter saying : Error: ... (error concerning the function on the last line of text) Execution halted Perhaps some of you have experienced this and found work arounds? One has to be pretty good to type without ever committing mistakes! By the way, I am running Cygwin (cygwin_NT-5.0 release 1.5.10(0.116/4/2), I reinstalled everything fresh last week) on Win 2k and R 1.9.1. Although some of you find cygwin inefficient in many ways, I don't have the option to migrate to Linux. Cheers, Thomas *** Le contenu de cet e-mail et de ses pièces jointes est destiné à l'usage exclusif du (des) destinataire(s) expressément désigné(s) comme tel(s). En cas de réception de cet e-mail par erreur, le signaler à son expéditeur et ne pas en divulguer le contenu. L'absence de virus a été vérifié à l'émission du message. Il convient néanmoins de vérifier l'absence de corruption à sa réception. The contents of this email and any attachments are confidential. They are intended for the named recipient(s) only. If you have received this email in error please notify the system manager or the sender immediately and do not disclose the contents to anyone or make copies. eSafe scanned this email for viruses, vandals and malicious content. *** __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] package.contents() deprecated, what replaces it?
Back on 18-June a colleague (Earl Glynn) asked: == I did a ?package.contents in version 1.9.0. The help file says this is deprecated, and says to see also Deprecated and Defunct, but what is the new function that replaces packge.contents? The function still seems to work, but it's just tagged as deprecated. efg == I have the same question, and also could not find the answer despite searching everywhere. I couldn't find a followup to the question. I've written an instruction manual introducing R to my institute's students, which includes the suggestion to use this command, e.g. For example, to see what's in the geostatistical package \texttt{gstat}: \begin{verbatim} package.contents(gstat) \end{verbatim} but now they see an off-putting deprecated message. Please advise. D G RossiterD G Rossiter Senior University Lecturer Department of Earth Systems Analysis (DESA) International Institute for Geo-Information Science and Earth Observation (ITC) Hengelosestraat 99 PO Box 6, 7500 AA Enschede, The Netherlands Phone: +31 (0)53 4874 499 Fax:+31 (0)53 4874 336 mailto:[EMAIL PROTECTED], Internet: http://www.itc.nl/personal/rossiter __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Re: Writing polyline shapefiles
On Wed, 28 Jul 2004, Patrick Giraudoux wrote: Stephane, dear R listers, I wonder if the idea of developping some functions in R writing shapefiles from a matrix of coordinates -eg poly2shape()- is still with us? I currently manage with a home-made function writing an ASCII grass file, then importing it into GRASS as vector file, then exporting the vector file as a shapefile (ouf!). The main issue to me is still to be capable to wite the *.shx and *.shp from coordinates (the *.dbf is easier to manage, even indirectly, and can for instance be completed simply via Excel). I tried to open *.shp, etc.. files to understand how they were structured and thus write some programme by myself, but unsuccessfully (the way those files are coded is not ascii) Maptools is fantastic for importation and data handling within R, but unfortunately cannot export to shapefiles. Any news? Nicholas Lewin-Koh wrote a simple sketch for type 1 (points) files (undocumented C function shpwrite). Lines are probably less interesting, so the polygon case remains. Is it important that the function pays close attention to the shapefile specification (ring direction, etc., lakes and islands) or can it assume that the user has everything under control? The ESRI specification is clear enough, and the shapelib code is OK. The longer term aim is to use the sourceforge sp package of spatial data classes in S4 form to read and write SpatialDataFramePolygons, but I guess doing a simpler function using S3 will get done sooner. Roger Patrick Giraudoux PS: ESRI has managed to distribute ArcTools (with ArcGIS) with no way to importe generate ASCII formats!!! One must get ArcINFO for that... No comment! -- 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] __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] elegant matrix creation
It's not obvious what pattern you want, but some variants of the following would work for jj1 and jj2: 1+outer(1:9, 1:9, +)%%3 In particular the following are equal to your jj1 and jj2: jj1. - 1+outer(0:8, rep(0:2, e=3), +)%%3 jj2. - 1+outer(0:8, c(1,2,1,3,1,3,2,3,2)-1, +)%%3 I couldn't figure out jj3, but this system may not work so easily for that. hope this helps. spencer graves Robin Hankin wrote: Hello everybody. I am trying to reproduce a particular matrix in an elegant way. If I have jj1 - structure(c(1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2, 3,1,2,3,1,2,3,1,2,3,2,3,1,2,3,1,2,3,1,2,3, 1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,3,1,2,3,1, 2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1, 2),.Dim = as.integer(c(9,9))) [ image(jj1) is good here ] then I can get this with kronecker(matrix(1,3,1),kronecker(1+outer(0:2,0:2,+)%%3,matrix(1,1,3))) I want to reproduce the following matrices in an equivalent way: jj2 - matrix(c(1,2,3,1,2,3,1,2,3,2,3,1,2,3,1,2,3,1, 1,2,3,1,2,3,1,2,3,3,1,2,3,1,2,3,1,2,1,2,3,1,2, 3,1,2,3,3,1,2,3,1,2,3,1,2,2,3,1,2,3,1,2,3,1,3, 1,2,3,1,2,3,1,2,2,3,1,2,3,1,2,3,1),9,9) jj3 - structure(c(1,2,3,2,3,1,3,1,2,1,2,1,2,3,2,3,1, 3,1,3,1,2,1,2,3,2,3,2,3,1,3,1,2,1,2,3,2,3, 2,3,1,3,1,2,1,2,1,2,3,2,3,1,3,1,3,1,2,1,2, 3,2,3,1,3,1,3,1,2,1,2,3,2,3,2,3,1,3,1,2,1, 2),.Dim = as.integer(c(9,9))) [ note that jj1-jj3 each have precisely 3 occurrences of A, B, and C along each row, column and (broken) diagonal ]. Can anyone give me a nice elegant way of creating jj2 and jj3 please? __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] using Rterm under cygwin, no possiblity to delete characters
You have to call Rterm.exe *from a command shell*. `cygwin' is not a command shell, but a collection of tools that provides several such. So we can only guess at what you are using. Rterms runs satisfactorily in many shells (I use tcsh, others use Cygwin bash ...). This is all described in the README and rw-FAQ files, for example. On Wed, 28 Jul 2004, Dewez Thomas wrote: Dear R-users, When I call Rterm from cygwin, I have no options but typing the exact syntax the first time. If I happen to hit the delete key (backspace), R dies when I press enter saying : backspace and delete are separate keys, so which did you mean? Error: ... (error concerning the function on the last line of text) Execution halted Perhaps some of you have experienced this and found work arounds? One has to be pretty good to type without ever committing mistakes! By the way, I am running Cygwin (cygwin_NT-5.0 release 1.5.10(0.116/4/2), I reinstalled everything fresh last week) on Win 2k and R 1.9.1. Although some of you find cygwin inefficient in many ways, I don't have the option to migrate to Linux. The alternative is to use a Windows-native shell, such as tcsh.exe. -- 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 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] package.contents() deprecated, what replaces it?
package.contents(gstat) is really ugly, and not something we would ever recommend for end users. E.g. library(help=gstat) has always done a better job. The information package.contents regurgitates is itself deprecated, to be replaced by metadata (in package subdir Meta) that other tools, e.g. library(help=), can format nicely. Had namespaces be available back then, package.contents would never have been user-visible. On Wed, 28 Jul 2004, David Rossiter wrote: Back on 18-June a colleague (Earl Glynn) asked: == I did a ?package.contents in version 1.9.0. The help file says this is deprecated, and says to see also Deprecated and Defunct, but what is the new function that replaces packge.contents? The function still seems to work, but it's just tagged as deprecated. efg == I have the same question, and also could not find the answer despite searching everywhere. I couldn't find a followup to the question. I've written an instruction manual introducing R to my institute's students, which includes the suggestion to use this command, e.g. For example, to see what's in the geostatistical package \texttt{gstat}: \begin{verbatim} package.contents(gstat) \end{verbatim} but now they see an off-putting deprecated message. Please advise. D G RossiterD G Rossiter Senior University Lecturer Department of Earth Systems Analysis (DESA) International Institute for Geo-Information Science and Earth Observation (ITC) Hengelosestraat 99 PO Box 6, 7500 AA Enschede, The Netherlands Phone:+31 (0)53 4874 499 Fax: +31 (0)53 4874 336 mailto:[EMAIL PROTECTED], Internet: http://www.itc.nl/personal/rossiter __ [EMAIL PROTECTED] mailing list https://www.stat.math.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 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] elegant matrix creation
Spencer thank you! just what I wanted. I have a feeling that jj3 is qualitatively different from the other two. Nevertheless, I'm sure it'll crack sooner or later! best wishes rksh At 07:34 am -0700 28/07/04, Spencer Graves wrote: It's not obvious what pattern you want, but some variants of the following would work for jj1 and jj2: 1+outer(1:9, 1:9, +)%%3 In particular the following are equal to your jj1 and jj2: jj1. - 1+outer(0:8, rep(0:2, e=3), +)%%3 jj2. - 1+outer(0:8, c(1,2,1,3,1,3,2,3,2)-1, +)%%3 I couldn't figure out jj3, but this system may not work so easily for that. hope this helps. spencer graves Robin Hankin wrote: Hello everybody. I am trying to reproduce a particular matrix in an elegant way. If I have jj1 - structure(c(1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2, 3,1,2,3,1,2,3,1,2,3,2,3,1,2,3,1,2,3,1,2,3, 1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,3,1,2,3,1, 2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1, 2),.Dim = as.integer(c(9,9))) [ image(jj1) is good here ] then I can get this with kronecker(matrix(1,3,1),kronecker(1+outer(0:2,0:2,+)%%3,matrix(1,1,3))) I want to reproduce the following matrices in an equivalent way: jj2 - matrix(c(1,2,3,1,2,3,1,2,3,2,3,1,2,3,1,2,3,1, 1,2,3,1,2,3,1,2,3,3,1,2,3,1,2,3,1,2,1,2,3,1,2, 3,1,2,3,3,1,2,3,1,2,3,1,2,2,3,1,2,3,1,2,3,1,3, 1,2,3,1,2,3,1,2,2,3,1,2,3,1,2,3,1),9,9) jj3 - structure(c(1,2,3,2,3,1,3,1,2,1,2,1,2,3,2,3,1, 3,1,3,1,2,1,2,3,2,3,2,3,1,3,1,2,1,2,3,2,3, 2,3,1,3,1,2,1,2,1,2,3,2,3,1,3,1,3,1,2,1,2, 3,2,3,1,3,1,3,1,2,1,2,3,2,3,2,3,1,3,1,2,1, 2),.Dim = as.integer(c(9,9))) [ note that jj1-jj3 each have precisely 3 occurrences of A, B, and C along each row, column and (broken) diagonal ]. Can anyone give me a nice elegant way of creating jj2 and jj3 please? -- Robin Hankin Uncertainty Analyst Southampton Oceanography Centre SO14 3ZH tel +44(0)23-8059-7743 [EMAIL PROTECTED] (edit in obvious way; spam precaution) __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] as(obj,matrix)
Hi! Here a simple example. setClass(myclass ,representation(info=character) ,contains=matrix ) rownames(dd)-c(a,b) tt-new(myclass,dd) #the source of pain. as(tt,matrix)-matrix(1,3,3) Error: length of dimnames [1] not equal to array extent Is there a different way to do what I would like to do (I would like to change the @.Data and all its attributes in the object)? Eryk. __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] a bug with LAPACK ? non orthogonal vectors obtained with eigen of a symmetric matrix
Hello, I have send send this message one week ago but I have receive no answer. Perhaps, some of RListers were in holidays and do not read my message. I try again.. My problem is that I obtained non orthonormal eigenvectors with some matrices with LAPACK while EISPACK seems to provide good results. Is there some restrictions with the use of LAPACK ? Is it a bug ? I did not find the answer. Here is my experiment: I have obtained strange results using eigen on a symmetric matrix: # this function perform a double centering of a matrix (xij-rowmean(i)-colmean(j)+meantot) dbcenter=function(mat){ rmean=apply(mat,1,mean) cmean=apply(mat,2,mean) newmat=sweep(mat,1,rmean,-) newmat=sweep(newmat,2,cmean,-) newmat=newmat+mean(mat) newmat} # i use spdep package to create a spatial contiguity matrix library(spdep) x=dbcenter(nb2mat(cell2nb(3,3),style=B)) #compute eigenvalues of a 9 by 9 matrix res=eigen(x) # some eigenvalues are equal to 0 eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) # I remove the corresponding eigenvectors res0=res$vec[,-which(eq0)] # then I compute the Froebenius norm with the identity matrix sum((diag(1,ncol(res0))-crossprod(res0))^2) [1] 1.515139e-30 # The results are correct, # then I do it again with a biggest matrix(100 by 100) x=dbcenter(nb2mat(cell2nb(10,10),style=B)) res=eigen(x) eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) res0=res$vec[,-which(eq0)] sum((diag(1,ncol(res0))-crossprod(res0))^2) [1] 3.986387 I have try the same with res=eigen(x,EISPACK=T): x=dbcenter(nb2mat(cell2nb(10,10),style=B)) res=eigen(x,EISPACK=T) eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) res0=res$vec[,-which(eq0)] sum((diag(1,ncol(res0))-crossprod(res0))^2) [1] 1.315542e-27 So I wonder I there is a bug in the LAPACK algorithm or if there are some things that I have not understood. Note that my matrix has some pairs of equal eigenvalues. Thanks in advance. I have continue my experiments in changing the size of my matrix : (3^2 by 3^2, 4^2 by 4^2... 20^2 by 20^2) EISPACK is always correct but LINPACK provide very strange results: for(i in 3:20){ + x=dbcenter(nb2mat(cell2nb(i,i),style=B)) + res=eigen(x,EIS=T) + eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) + res0=res$vec[,-which(eq0)] + print(sum((diag(1,ncol(res0))-crossprod(res0))^2)) + } [1] 7.939371e-30 [1] 2.268788e-29 [1] 9.237286e-29 [1] 1.806393e-28 [1] 3.24619e-28 [1] 5.239195e-28 [1] 9.78079e-28 [1] 1.315542e-27 [1] 1.838600e-27 [1] 3.114150e-27 [1] 5.499297e-27 [1] 5.471782e-27 [1] 1.075098e-26 [1] 1.534822e-26 [1] 1.771326e-26 [1] 2.342404e-26 [1] 3.462522e-26 [1] 4.310143e-26 for(i in 3:20){ + x=dbcenter(nb2mat(cell2nb(i,i),style=B)) + res=eigen(x) + eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) + res0=res$vec[,-which(eq0)] + print(sum((diag(1,ncol(res0))-crossprod(res0))^2)) + } [1] 1.515139e-30 [1] 1.054286e-27 [1] 9.553017e-29 [1] 2.263455e-28 [1] 5.641993e-27 [1] 4.442088e-26 [1] 3.996714 [1] 3.986387 [1] 3.996545 [1] 7.396718 [1] NaN [1] 7.980621 [1] 7.996769 [1] 3.984399 [1] NaN [1] NaN [1] NaN [1] NaN Note that I have do the same with random number and never find this kind of problems R.Version() $platform [1] i386-pc-mingw32 $arch [1] i386 $os [1] mingw32 $system [1] i386, mingw32 $status [1] $major [1] 1 $minor [1] 9.1 $year [1] 2004 $month [1] 06 $day [1] 21 $language [1] R Stéphane DRAY -- Département des Sciences Biologiques Université de Montréal, C.P. 6128, succursale centre-ville Montréal, Québec H3C 3J7, Canada Tel : (514) 343-6111 poste 1233 Fax : (514) 343-2293 E-mail : [EMAIL PROTECTED] -- Web http://www.steph280.freesurf.fr/ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Parallel Functions on AIX
Dear R Development Team, Does the latest version R-1.9.1 provide parallel functions if R is running on Multiple CPUs Unix platform (IBM AIX e-server)? Kexiao [[alternative HTML version deleted]] __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] covariate selection in cox model (counting process)
On Wed, 28 Jul 2004, Mayeul KAUFFMANN wrote: No, I mean recurrent events. With counting process notation but no recurrent revents the partial likelihood is still valid, and the approach of treating it as a real likelihood for AIC (and presumably BIC) makes sense. Roughly speaking, you can't tell there is dependence until you see multiple events. Thanks a lot, I got it (well, I hope so)! I've read in several places that events in the Andersen-Gill model must be conditionnaly independent, which is sometimes more precisely written as conditionnaly independent given the covariates or even more precisely: the Andersen-Gill (AG) model assumes that each [individual] has a multi-event counting process with independent increments. The observed increments must be conditionally independent given the history of all observable information up to the event times. (http://www.stat.umu.se/egna/danardono/licdd.pdf) More precisely still, for the criterion function in coxph() to be a partial likelihood the estimating function must be a martingale. This is actually a slightly weaker assumption than independent increments. The proportional rates model doesn't require this assumption, and is also sometimes called the Andersen-Gill model. The criterion function isn't a likelihood but it still gives valid estimators. Then, there is still another option. In fact, I already modelled explicitely the influence of past events with a proximity of last event covariate, assuming the dependence on the last event decreases at a constant rate (for instance, the proximity covariate varies from 1 to 0.5 in the first 10 years after an event, then from 0.5 to 0.25 in the next ten years, etc). With a well chosen modelisation of the dependence effect, the events become conditionnaly independent, I do not need a +cluster(id) term, and I can use fit$loglik to make a covariate selection based on BIC, right? If you can get the conditional independence (martingaleness) then, yes, BIC is fine. One way to check might be to see how similar the standard errors are with and without the cluster(id) term. -thomas Thanks a lot again for your time. Mayeul KAUFFMANN Univ. Pierre Mendes France Grenoble - France PS: I wrongly concluded from the R statement (Note: the likelihood ratio and score tests assume independence of observations within a cluster, the Wald and robust score tests do not). that it meant independence between two consecutive observations (without any event). It made sense to me because when only one covariate changes for a given individual, and with a small change, there is a new observation, with a risk very simlar to the risk for the previous observation. But there is still independence with respect to the question of recurrent event. Maybe the warning should be rewritten saying assume *conditionnal* independence of *events* (given the covariates) __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Thomas Lumley Assoc. Professor, Biostatistics [EMAIL PROTECTED] University of Washington, Seattle __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] Another big data size problem
-Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Federico Gherardini Sent: Wednesday, July 28, 2004 5:26 PM To: [EMAIL PROTECTED] Subject: Re: [R] Another big data size problem On Wed, 28 Jul 2004 13:28:20 +0100 Ernesto Jardim [EMAIL PROTECTED] wrote: Hi, When you're writing a table to MySQL you have to be carefull if the table is created by RMySQL. The fields definition may not be the most adequate and there will be no indexes in your table, which makes the queries _very_ slow. So, if I understood correctly, if you want to use SQL you'll have to upload the table in SQL, directly from MySQL without using R at all, and then use RMySQL to read the elements in R? Uwe Ligges [EMAIL PROTECTED] wrote: Note that it is better to initialize the object to full size before inserting -- rather than using rbind() and friends which is indeed slow since it need to re-allocate much memory for each step. Do you mean something like this? tab - matrix(rep(0, 1227 * 2), 1227, 2, byrow = TRUE) for(i in 0:num.lines) tab[i + 1,] - scan(file=fh, nlines=1, what=PS, skip = i) It is better to open a file connection, keep it open during the loop and the close it afterwards. Something like tab - matrix(rep(0, 1227 * 2), 1227, 2, byrow = TRUE) fh - file(filename, open=r); for(i in 0:num.lines) tab[i + 1,] - scan(file=fh, nlines=1); close(fh); As you have done it, the file is opened once in each iteration of the loop, scan() starts reading from the beginning, parse all lines to skip 'i' lines, and the reads one line. This is done num.lines+1 times! Anyway, I think you also should read the help for scan(). What do you want with argument 'what=PS'? PS is not a valid data type; 'what' does not specify a name of field/column to be read. The above doesn't get very far either... it seems that, once it has created the table, it becomes so slow that it's unusable. I'll have to try this with more RAM by the way. My suggestions to you are that try read.table() with specified data type for the columns using vector argument 'colClasses'. This way you can help R by specify that, say, column 3 is an integer (have the memory of a double), and that column 6-10 are doubles. Unfortunately you can tell read.table() to skip some of the columns that you are not interested in, which in your case to help you out a lot. To do this, you have to use scan(), which read.table() uses internally. In scan() 'what' works similar to 'colClasses' *and* if you specify 'what' as a 'list' you can tell scan() to skip some columns by setting its 'what' value to NULL, e.g. what=list(integer, integer, NULL, double, character). I think you can get pretty far doing this! Cheers, fede Good luck! Henrik Bengtsson __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] using Rterm under cygwin, no possiblity to delete characters
Hi. I used to run R under Cygwin for quite some time, but after, I think, R v1.8.1 or so, the problems started to show up. There has been some discussion about (recent) problems running R under Cygwin. Once problem was for instance that hitting Ctrl-C in R running Cygwin would bring you into a concurrent forked R and shell environment. Search the r-help archive for that discussion. Anyway, the final take home message, which was underlined by B. Ripley, was that R does *not* support Cygwin. For command line R, I run R in a specially setup Command prompt. It works fine indeed; Below is the bat-file that I use to initiate a new Command prompt to run Rterm in. If it does not run out of the box for you, it should not be hard to modify. Cheers Henrik Bengtsson @echo off rem ## rem # Usage: RCMDprompt.bat [path] rem # rem # This script opens a MS-DOS prompt with a enviroment variables rem # set such that R can be ran and packages can be build. rem # If 'path' is given, the working directory will be set accordingly. rem # rem # NOTE: This scripts works even if Cygwin is installed. HOWEVER, you rem # can not have any Cygwin applications running (not even a shell or rem # XEmacs for Cygwin) at the same time you try to run RCMD. rem # rem # Requires: rem # To build and install packages two things must be installed, i.e. rem # exists in the PATH. First, the Rtools compilation [1,2] must exists. rem # The path (R_TOOLS) to it is set below. Second, Perl (must not be rem # Cygwin/Perl) must also exists. The path to it is set below. rem # rem # Reference: rem # [1] http://www.stats.ox.ac.uk/pub/Rtools/ rem # [2] http://www.murdoch-sutherland.com/Rtools/ rem # [3] Duncan Murdoch, Using MiKTeX with R for Windows, 2004. rem # http://www.murdoch-sutherland.com/Rtools/miktex.html rem # rem # Henrik Bengtsson, [EMAIL PROTECTED], March-June 2004. rem ## rem # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - rem # 1. Global environment variables rem # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - rem # Short version of PROGRAMFILES, e.g. 'C:\Progra~1' instead of rem # 'C:\Program Files\', which contains spaces that are BAD for R Co. rem # rem # On Windows XP and NT, variable substitution using 'for' can be used rem # to do this automatically from %ProgramFiles%. For details see rem # http://www.microsoft.com/windowsxp/home/using/productdoc/en/for.asp for %%v in (%ProgramFiles%) do set PROGRAMFILES_SHORT=%%~sv rem # Set the main R directory set R_ROOT=%PROGRAMFILES_SHORT%\R rem # Set the R_HOME directory set R_HOME=%R_ROOT%\rw1090 rem # Set the HOME directory. This is the directory where R looks rem # for the .Rprofile and .Renviron files. See ?Startup. set HOME=%UserProfile% rem # Set TMPDIR to a temporary directory set TMPDIR=%TEMP% rem # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - rem # 2. Setup the PATH rem # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - rem # Clear the PATH (making the main Cygwin installation invisible) path ; rem # Set the LaTeX directory rem # http://www.miktex.org/ path %SystemDrive%\texmf\miktex\bin;%PATH% rem # There are recent issues with MikTeX v2.4 and R. The problem occurs rem # because e-TeX is used instead of TeX. See [3] for details. rem # Try latex -version. If you see MikTeX-e-TeX then you should read rem # the instructions at [3]. rem # Set the Microsoft HTML Help Compiler directory rem # http://msdn.microsoft.com/library/tools/htmlhelp/chm/HH1Start.htm path %ProgramFiles%\HTML Help Workshop;%PATH% rem # Set the Perl directory rem # http://www.activestate.com/Products/ActivePerl/Download.html path %SystemDrive%\Perl\bin;%PATH% rem # Set the Rtools [1] directory rem # http://www.stats.ox.ac.uk/pub/Rtools/ path %R_ROOT%\Rtools;%PATH% rem # Set the R_HOME directory path %R_HOME%\bin;%PATH% rem # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - rem # 3. Start the MSDOS prompt in the given directory rem # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - rem # Change directory according to argument 1 cd /D %1 rem # Start the MSDOS commando prompt %SystemRoot%\system32\cmd.exe -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Dewez Thomas Sent: Wednesday, July 28, 2004 4:24 PM To: '[EMAIL PROTECTED]' Subject: [R] using Rterm under cygwin, no possiblity to delete characters Dear R-users, When I call Rterm from cygwin, I have no options but typing the exact syntax the first time. If I happen to hit the delete key (backspace), R dies when I press enter saying : Error: ... (error concerning the function on the last line of text) Execution halted Perhaps some of you have experienced this and found work arounds? One has to be
Re: [R] Simulation from a model fitted by survreg.
On Wed, 28 Jul 2004, Sixten Borg wrote: Dear list, I would like to simulate individual survival times from a model that has been fitted using the survreg procedure (library survival). Output shown below. My plan is to extract the shape and scale arguments for use with rweibull() since my error terms are assumed to be Weibull, but it does not make any sense. The mean survival time is easy to predict, but I would like to simulate individual survival times. It's parametrized differently from rweibull. The linear predictor for the model is the logarithm of the scale, and the scale parameter from the model is the reciprocal of the shape. As with the gamma family you always need to check the parametrisation for weibull distributions. -thomas I am probably missing something completely obvious. Any hints or advice are appreciated. Thanks Sixten summary(mod1) Call: survreg(formula = Surv(tid, study$first.event.death) ~ regim + age + stadium2, data = study, dist = weibull) Value Std. Error zp (Intercept) 11.6005 0.7539 15.387 2.01e-53 regimposto -0.1350 0.1558 -0.867 3.86e-01 age -0.0362 0.0102 -3.533 4.11e-04 stadium2ii -0.0526 0.2794 -0.188 8.51e-01 Log(scale) -0.5148 0.1116 -4.615 3.93e-06 Scale= 0.598 Weibull distribution Loglik(model)= -680.7 Loglik(intercept only)= -689.2 Chisq= 16.87 on 3 degrees of freedom, p= 0.00075 Number of Newton-Raphson Iterations: 8 n=1183 (4 observations deleted due to missing) version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major1 minor8.1 year 2003 month11 day 21 language R __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Thomas Lumley Assoc. Professor, Biostatistics [EMAIL PROTECTED] University of Washington, Seattle __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Gaussian frailty leads to segmentation fault
Dear R gurus, for a simulation concerning study effects and historical controls in survival analysis, i would like to experiment with a gaussian frailty model. The simulated scenario consists of a randomized trial (treatment and placebo) and historical controls (only placebo). So the simulated data frames consist of four columns $time, $cens, $study, $treat. $time, $cens are the usual survival data. For the binary thretment indicator we have $treat == 0 or 1, if $study == 1, $treat == 1 if $study 1 Typical parameters for my simulations are: sample sizes (per arm): between 100 and 200 number of historical studies: between 7 and 15 hazard ratio treatment/placebo: between 0.7 and 1 variance of the study effekt: between 0 and 0.3 Depending on the sample sizes, the following call sometimes leads to a segmentation fault: coxph(Surv(time,cens) ~ as.factor(treatment) + frailty(study, distribution=gaussian), data=data) I noticed, that this segmentation fault occures most frequently, if the number of randomized treatment patients is higher than the number of randomized placebo patients, and the number of historical studies is large. There seems to be no problem, if there are at least as many randomized placebo patients as treated patients. Unfortunately, this is not the situation i want to investigate (historical controls should be used to decrease the number of treated patients). Is there a way to circumwent this problem? Christian P.S. Is it allowed, to attach gzipped sample data sets in this mailing list? __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] How to add an object to an RData file ?
Hi, I've saved an RData file with save and now I want to add a new object to this file. At the moment I do: attach(file.RData) save(list=c(new,obj, ls(pos=2)), file=file.RData, compress=T) detach() Is there a quicker method that just add the object to the file ? Thanks EJ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] as(obj,matrix)
The definition of dd is. dd-matrix(0,2,2) *** REPLY SEPARATOR *** On 7/28/2004 at 4:12 PM Adaikalavan Ramasamy wrote: On Wed, 2004-07-28 at 15:55, Wolski wrote: Hi! Here a simple example. setClass(myclass ,representation(info=character) ,contains=matrix ) rownames(dd)-c(a,b) Er, I don't think you have defined 'dd' in the example. tt-new(myclass,dd) #the source of pain. as(tt,matrix)-matrix(1,3,3) Error: length of dimnames [1] not equal to array extent Is there a different way to do what I would like to do (I would like to change the @.Data and all its attributes in the object)? Eryk. __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Dipl. bio-chem. Eryk Witold Wolski@MPI-Moleculare Genetic Ihnestrasse 63-73 14195 Berlin 'v' tel: 0049-30-83875219 / \ mail: [EMAIL PROTECTED]---W-Whttp://www.molgen.mpg.de/~wolski __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Modelling compound logistic growth curves
Motivated by the discovery of 'loglet analysis' (http://phe.rockefeller.edu/LogletLab/) that allows one to decompose growth curves into a series of logistic equations, I attempted to do the same thing in R. SIMULATED DATA Time - 1:200 pop.size - SSlogis(Time,10,20,5) + SSlogis(Time,20,100,20) + rnorm(length(Time)) MY ANALYSIS results - nls(size ~ SSlogis(Time, Asym1, xmid1, scal1) + SSlogis(Time, Asym2, xmid2, scal2), start = list(Asym1=5, xmid1=15, scal1=30, Asym2=25, xmid2=67, scal2=25)) THE RESULT I get the error message: Error in nls(size ~ SSlogis(Time, Asym1, xmid1, scal1) + SSlogis(Time, : step factor 0.000488281 reduced below `minFactor' of 0.000976563 Assistance in doing this analysis would be much appreciated. Dan Bebber Dr. Daniel P. Bebber Department of Plant Sciences University of Oxford South Parks Road Oxford OX1 3RB Tel. 01865 275060 Web. http://www.forestecology.co.uk/ Data, data, data! he cried impatiently. I can't make bricks without clay - Sherlock Holmes, The Adventure of the Copper Beeches, 1892 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Parallel Functions on AIX
The first thing to ascertain is that R will actually build on such a machine -- there have been a lot of reports of failure on AIX, and no recent reports of success. If it does, R itself is single-threaded (but can make use of multi-threaded BLAS) but packages such as Rmpi, snow, RScaLAPACK provide parallel facilities (and are in the FAQ). On Wed, 28 Jul 2004, Liao, Kexiao wrote: Dear R Development Team, You actually wrote to R-help! Does the latest version R-1.9.1 provide parallel functions if R is running on Multiple CPUs Unix platform (IBM AIX e-server)? -- 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 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] a bug with LAPACK ? non orthogonal vectors obtained with eigen of a symmetric matrix
This is interesting. I can reproduce your results but cannot come up with an explanation. However, using svd(LINPACK = FALSE) seems to work every time. Might you consider trying that instead? -roger Stephane DRAY wrote: Hello, I have send send this message one week ago but I have receive no answer. Perhaps, some of RListers were in holidays and do not read my message. I try again.. My problem is that I obtained non orthonormal eigenvectors with some matrices with LAPACK while EISPACK seems to provide good results. Is there some restrictions with the use of LAPACK ? Is it a bug ? I did not find the answer. Here is my experiment: I have obtained strange results using eigen on a symmetric matrix: # this function perform a double centering of a matrix (xij-rowmean(i)-colmean(j)+meantot) dbcenter=function(mat){ rmean=apply(mat,1,mean) cmean=apply(mat,2,mean) newmat=sweep(mat,1,rmean,-) newmat=sweep(newmat,2,cmean,-) newmat=newmat+mean(mat) newmat} # i use spdep package to create a spatial contiguity matrix library(spdep) x=dbcenter(nb2mat(cell2nb(3,3),style=B)) #compute eigenvalues of a 9 by 9 matrix res=eigen(x) # some eigenvalues are equal to 0 eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) # I remove the corresponding eigenvectors res0=res$vec[,-which(eq0)] # then I compute the Froebenius norm with the identity matrix sum((diag(1,ncol(res0))-crossprod(res0))^2) [1] 1.515139e-30 # The results are correct, # then I do it again with a biggest matrix(100 by 100) x=dbcenter(nb2mat(cell2nb(10,10),style=B)) res=eigen(x) eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) res0=res$vec[,-which(eq0)] sum((diag(1,ncol(res0))-crossprod(res0))^2) [1] 3.986387 I have try the same with res=eigen(x,EISPACK=T): x=dbcenter(nb2mat(cell2nb(10,10),style=B)) res=eigen(x,EISPACK=T) eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) res0=res$vec[,-which(eq0)] sum((diag(1,ncol(res0))-crossprod(res0))^2) [1] 1.315542e-27 So I wonder I there is a bug in the LAPACK algorithm or if there are some things that I have not understood. Note that my matrix has some pairs of equal eigenvalues. Thanks in advance. I have continue my experiments in changing the size of my matrix : (3^2 by 3^2, 4^2 by 4^2... 20^2 by 20^2) EISPACK is always correct but LINPACK provide very strange results: for(i in 3:20){ + x=dbcenter(nb2mat(cell2nb(i,i),style=B)) + res=eigen(x,EIS=T) + eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) + res0=res$vec[,-which(eq0)] + print(sum((diag(1,ncol(res0))-crossprod(res0))^2)) + } [1] 7.939371e-30 [1] 2.268788e-29 [1] 9.237286e-29 [1] 1.806393e-28 [1] 3.24619e-28 [1] 5.239195e-28 [1] 9.78079e-28 [1] 1.315542e-27 [1] 1.838600e-27 [1] 3.114150e-27 [1] 5.499297e-27 [1] 5.471782e-27 [1] 1.075098e-26 [1] 1.534822e-26 [1] 1.771326e-26 [1] 2.342404e-26 [1] 3.462522e-26 [1] 4.310143e-26 for(i in 3:20){ + x=dbcenter(nb2mat(cell2nb(i,i),style=B)) + res=eigen(x) + eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) + res0=res$vec[,-which(eq0)] + print(sum((diag(1,ncol(res0))-crossprod(res0))^2)) + } [1] 1.515139e-30 [1] 1.054286e-27 [1] 9.553017e-29 [1] 2.263455e-28 [1] 5.641993e-27 [1] 4.442088e-26 [1] 3.996714 [1] 3.986387 [1] 3.996545 [1] 7.396718 [1] NaN [1] 7.980621 [1] 7.996769 [1] 3.984399 [1] NaN [1] NaN [1] NaN [1] NaN Note that I have do the same with random number and never find this kind of problems R.Version() $platform [1] i386-pc-mingw32 $arch [1] i386 $os [1] mingw32 $system [1] i386, mingw32 $status [1] $major [1] 1 $minor [1] 9.1 $year [1] 2004 $month [1] 06 $day [1] 21 $language [1] R Stéphane DRAY -- Département des Sciences Biologiques Université de Montréal, C.P. 6128, succursale centre-ville Montréal, Québec H3C 3J7, Canada Tel : (514) 343-6111 poste 1233 Fax : (514) 343-2293 E-mail : [EMAIL PROTECTED] -- Web http://www.steph280.freesurf.fr/ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Underline in expression().
John Janmaat wrote: Sundar, Thanks. Unfortunately, I am looking for something that also works in the margins of the plot. As a workaround, what about frac() (well, vertical justification is not that perfect in this case)? plot(1:10, xlab = ) title(xlab = expression(frac(y == X * beta + e)), line = 4.5) I'll probably take a look how to add that feature during August. Uwe Ligges John. Sundar Dorai-Raj wrote: John Janmaat wrote: Hello All, Is there an analogue to \underbar or the AMS math \underline in graphical math expressions? Thanks, John. Uwe Ligges posted a solution a couple of years ago. I don't know if there is anything built in yet. ?plotmath does not seem to say anything about underlining. http://finzi.psych.upenn.edu/R/Rhelp01/archive/7191.html plot(0:1, 0:1, type=n) underlined(0.5, 0.5, expression(widehat(x %*% y))) --sundar __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] as(obj,matrix)
Wolski wrote: Hi! Here a simple example. setClass(myclass ,representation(info=character) ,contains=matrix ) rownames(dd)-c(a,b) tt-new(myclass,dd) #the source of pain. as(tt,matrix)-matrix(1,3,3) Error: length of dimnames [1] not equal to array extent Is there a different way to do what I would like to do (I would like to change the @.Data and all its attributes in the object)? It's not particularly a problem with as(). Your class just doesn't behave as you expect. matrix is not a formal class with slots. (It isn't even an S3 class in R; attr(x,class) is NULL.) So you cannot expect classes extending matrix to know what a matrix is supposed to be. Part of the problem is that matrix objects sometimes have dimnames and sometimes don't. And there is basic code in R that assumes or applies constraints on the dim or dimnames. In S-Plus, matrix is a formal class, always having a slot for dimnames. R has not gone that route, at least not yet. It may be possible to define a new class, Matrix, say, that looks like a matrix to old-style code but has a formal definition. But the details are likely to be tricky, and it's definitely a topic for r-devel, not r-help. John Chambers Eryk. __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- John M. Chambers [EMAIL PROTECTED] Bell Labs, Lucent Technologiesoffice: (908)582-2681 700 Mountain Avenue, Room 2C-282 fax:(908)582-3340 Murray Hill, NJ 07974web: http://www.cs.bell-labs.com/~jmc __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Question ?
HI, I would like use header file in the R-HOME/scr/include. What should I do please ? JL __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] elegant matrix creation
Not sure if this qualifies as elegant or not but it (1) does allow one to generate all three matrices using the same scheme, (2) is simple requiring only a single one line function, (3) reduces the number of numbers you must specify from 81 to 18 per matrix and (4) gives some kminimal insight into the patterns. The key observation is that each row of each matrix is a cyclically shifted version of the first row of that same matrix. Thus given the first row and a vector of shifts we can reconstruct the remaining rows. Define a function which shifts its vector argument v by shift positions to the left producing a one row matrix. If shift is a vector the each row of the result corresponds to one shift in vector shift. shiftL - function(v, shift) outer(shift,seq(along=v)-1, function(i,j)v[(i+j)%%length(v)+1]) # now run shiftL using the first row of each matrix and the shift vector # for each matrix jj1 - shiftL(c(1,1,1,2,2,2,3,3,3),c(0,3,6,0,3,6,0,3,6)) jj2 - shiftL(c(1,2,1,3,1,3,2,3,2),c(0,6,3,0,6,3,0,6,3)) jj3 - shiftL(c(1,1,1,2,2,2,3,3,3),c(0,4,8,3,7,2,6,1,5)) or turning the input vectors into expressions themselves: jj1 - shiftL( rep(1:3,c(3,3,3)), rep(c(0,3,6),3) ) jj2 - shiftL( c(shiftL(c(1,3,2),c(0,2,0))), rep(c(0,6,3),3) ) jj3 - shiftL( rep(1:3,c(3,3,3)), seq(0,32,4) %% 9 ) Robin Hankin rksh at soc.soton.ac.uk writes: : : Hello everybody. : : I am trying to reproduce a particular matrix in an elegant way. If I : have : : jj1 - : structure(c(1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2, : 3,1,2,3,1,2,3,1,2,3,2,3,1,2,3,1,2,3,1,2,3, : 1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,3,1,2,3,1, : 2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1, : 2),.Dim = as.integer(c(9,9))) : : [ image(jj1) is good here ] then I can get this with : : kronecker(matrix(1,3,1),kronecker(1+outer(0:2,0:2,+)%%3,matrix(1,1,3))) : : I want to reproduce the following matrices in an equivalent way: : : jj2 - matrix(c(1,2,3,1,2,3,1,2,3,2,3,1,2,3,1,2,3,1, : 1,2,3,1,2,3,1,2,3,3,1,2,3,1,2,3,1,2,1,2,3,1,2, : 3,1,2,3,3,1,2,3,1,2,3,1,2,2,3,1,2,3,1,2,3,1,3, : 1,2,3,1,2,3,1,2,2,3,1,2,3,1,2,3,1),9,9) : : jj3 - structure(c(1,2,3,2,3,1,3,1,2,1,2,1,2,3,2,3,1, : 3,1,3,1,2,1,2,3,2,3,2,3,1,3,1,2,1,2,3,2,3, : 2,3,1,3,1,2,1,2,1,2,3,2,3,1,3,1,3,1,2,1,2, : 3,2,3,1,3,1,3,1,2,1,2,3,2,3,2,3,1,3,1,2,1, 2),.Dim = : as.integer(c(9,9))) : : [ note that jj1-jj3 each have precisely 3 occurrences of A, B, and C : along each row, column and (broken) diagonal ]. : : Can anyone give me a nice elegant way of creating jj2 and jj3 please? : __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] How to add an object to an RData file ?
I can't think of a faster way. -roger Ernesto Jardim wrote: Hi, I've saved an RData file with save and now I want to add a new object to this file. At the moment I do: attach(file.RData) save(list=c(new,obj, ls(pos=2)), file=file.RData, compress=T) detach() Is there a quicker method that just add the object to the file ? Thanks EJ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Another big data size problem
Thanks for your suggestions, On Wed, 28 Jul 2004 17:04:57 +0200 Henrik Bengtsson [EMAIL PROTECTED] wrote: Anyway, I think you also should read the help for scan(). What do you want with argument 'what=PS'? PS is not a valid data type; 'what' does not specify a name of field/column to be read. From the scan help page... what: the type of 'what' gives the type of data to be read. It seems to me that I had read somewhere (maybe on the mailing list archives), that 'what' was supposed to be a sort of example of the kind of data you had to read... so I put a character string because I wanted the data to be read as character because I had a column of factors. I know this is not a great way to do it (better have a matrix made up of nubers only, instead of having to subsequently convert columns of character strings to number) but I wanted to do a quick test without having to rearrange my file. Cheers, fede __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] elegant matrix creation
[Sorry if this gets posted twice but I am having more gmane posting problems.] Not sure if this qualifies as elegant or not but it does (1) bring all three matrices under a single scheme, (2) reduce the number of numbers from 81 to 18 per matrix, (3) requires only a single one line utility function, (4) is simple and (5) gives some minimal insight into the patterns. The key thing to note is that each row of each matrix is a cyclic shift of the first row of that matrix. Define a shift function which shifts its vector argument v by shift positions to the left creating a one row matrix. If shift is a vector it creates a matrix with one row per shift. shiftL - function(v, shift) outer(shift,seq(along=v)-1, function(i,j)v[(i+j)%%length(v)+1]) jj1a - shiftL(c(1,1,1,2,2,2,3,3,3),c(0,3,6,0,3,6,0,3,6)) jj2a - shiftL(c(1,2,1,3,1,3,2,3,2),c(0,6,3,0,6,3,0,6,3)) jj3a - shiftL(c(1,1,1,2,2,2,3,3,3),c(0,4,8,3,7,2,6,1,5)) # or finding expressions for the two args in each case: jj1b - shiftL( rep(1:3,c(3,3,3)), rep(c(0,3,6),3) ) jj2b - shiftL( c(shiftL(c(1,3,2),c(0,2,0))), rep(c(0,6,3),3) ) jj3b - shiftL( rep(1:3,c(3,3,3)), seq(0,32,4) %% 9 ) __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] DESCRIPTION.in
Kevin To work around this problem I generate the DESCRIPTION and DESCRIPTION.in files from my Makefile. (BTW, using make to control the building process is pretty useful.) If I am building a single package I generate DESCRIPTION and if I am generating a bundle I generate DESCRIPTION.in. If you want more details on my whole Makefile system then let me know. Paul Gilbert Kevin Bartz wrote: Hello R world! I'm building a bundle of four packages, but I don't always want to build the whole bundle. Usually I just want to tweak one function in one of the packages and rebuild just that package. As such, I have DESCRIPTION and DESCRIPTION.in files sitting in all the package folders. Unfortunately, there's a pesky line in R CMD build that hacks away DESCRIPTION whenever it sees DESCRIPTION.in. This is okay for building the whole bundle, but it makes it a major pain to build any of the packages individually. I am root, so for now I have commented the offending line in /usr/local/lib/R/bin/build. Is this the proper solution or am I overlooking something? Thanks for any help you can provide, Kevin __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] elegant matrix creation
Robin Hankin wrote: Hello everybody. I am trying to reproduce a particular matrix in an elegant way. If I have jj1 - structure(c(1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2, 3,1,2,3,1,2,3,1,2,3,2,3,1,2,3,1,2,3,1,2,3, 1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,3,1,2,3,1, 2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1, 2),.Dim = as.integer(c(9,9))) [ image(jj1) is good here ] then I can get this with kronecker(matrix(1,3,1),kronecker(1+outer(0:2,0:2,+)%%3,matrix(1,1,3))) I want to reproduce the following matrices in an equivalent way: jj2 - matrix(c(1,2,3,1,2,3,1,2,3,2,3,1,2,3,1,2,3,1, 1,2,3,1,2,3,1,2,3,3,1,2,3,1,2,3,1,2,1,2,3,1,2, 3,1,2,3,3,1,2,3,1,2,3,1,2,2,3,1,2,3,1,2,3,1,3, 1,2,3,1,2,3,1,2,2,3,1,2,3,1,2,3,1),9,9) jj3 - structure(c(1,2,3,2,3,1,3,1,2,1,2,1,2,3,2,3,1, 3,1,3,1,2,1,2,3,2,3,2,3,1,3,1,2,1,2,3,2,3, 2,3,1,3,1,2,1,2,1,2,3,2,3,1,3,1,3,1,2,1,2, 3,2,3,1,3,1,3,1,2,1,2,3,2,3,2,3,1,3,1,2,1, 2),.Dim = as.integer(c(9,9))) some musings. You have not told us hove jj3 arises naturally, but its structure seems to have something to do with magic squares. So library(magic) # on CRAN is.magic(jj3) # TRUE This is not in accordance with my concept of magic squares, since an n*n magic square should have all the numbers from 1 to n^2 exactly once. But all rowSums and colSums are equal. Note that if we view jj3 as a 3*3 block matrix, with each block a 3*3 matrix, then all the blockx are generated rowwise the following way: Let x in 1:3 be the generator, and + be sum modulo 3, but we take the rep 3 and not 0. Then we have x xx x+1 x+1 x+2 x+2 xx and the 3*3 matrix of the generators can be generated in the following way: xx - ((magic(3) %% 3)+1)[,3:1] xx [,1] [,2] [,3] [1,]123 [2,]231 [3,]312 ++ - function(x,a) { if ( (x+a)%%3 == 0 ) 3 else (x+a) %% 3 } makeBlock - function(x) { +matrix( c( rep(x,3), rep(++(x,1), 2), rep(++(x,2), 2), + rep(x,2) ) , 3,3, byrow=TRUE) } makeBlock(1) [,1] [,2] [,3] [1,]111 [2,]223 [3,]311 ans - matrix (NA,9,9) for (i in 1:3) for (j in 1:3) { + ans[3*(i-1)+1:3, 3*(j-1)+1:3] - makeBlock(xx[i,j]) } ans [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [1,]111222333 [2,]223331112 [3,]311122233 [4,]222333111 [5,]331112223 [6,]122233311 [7,]333111222 [8,]112223331 [9,]233311122 ans ==jj3 [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [1,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [2,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [3,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [4,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [5,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [6,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [7,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [8,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [9,] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE Kjetil Halvorsen [ note that jj1-jj3 each have precisely 3 occurrences of A, B, and C along each row, column and (broken) diagonal ]. Can anyone give me a nice elegant way of creating jj2 and jj3 please? __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Question ?
Jean-Luc Allard wrote: HI, I would like use header file in the R-HOME/scr/include. What should I do please ? JL Just by using include and the filename. If you use R CMD SHLIB to create the library or the library is compiled during R CMD INSTALL, the path is included in the search path for header files. Uwe Ligges __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] covariate selection in cox model (counting process)
If you can get the conditional independence (martingaleness) then, yes, BIC is fine. One way to check might be to see how similar the standard errors are with and without the cluster(id) term. (Thank you again !, Thomas.) At first look, the values seemed very similar (see below, case 2). However, to check this without being too subjective, and without a specific test, I needed other values to assess the size of the differences: what is similar, what is not? == = CASE 1 I first estimated the model without modeling dependence: Call: coxph(formula = Surv(start, stop, status) ~ cluster(ccode) + pop + pib + pib2 + crois + instab.x1 + instab.autres, data = xstep) coef exp(coef) se(coef) robust se z p pop0.3606 1.434 0.09780.1182 3.05 2.3e-03 pib -0.5947 0.552 0.19520.1828 -3.25 1.1e-03 pib2 -0.4104 0.663 0.14520.1270 -3.23 1.2e-03 crois -0.0592 0.943 0.02450.0240 -2.46 1.4e-02 instab.x1 2.2059 9.079 0.46920.4097 5.38 7.3e-08 instab.autres 0.9550 2.599 0.47000.4936 1.93 5.3e-02 Likelihood ratio test=74 on 6 df, p=6.2e-14 n= 7286 There seems to be a strong linear relationship between standard errors (se, or naive se) and robust se. summary(lm(sqrt(diag(cox1$var))~ sqrt(diag(cox1$naive.var)) -1)) Coefficients: Estimate Std. Error t value Pr(|t|) sqrt(diag(cox1$naive.var)) 0.961030.04064 23.65 2.52e-06 *** Multiple R-Squared: 0.9911, Adjusted R-squared: 0.9894 == = CASE 2 Then I added a variable (pxcw) measuring the proximity of the previous event (1pxcw0) n= 7286 coef exp(coef) se(coef) robust se z p pxcw 0.9063 2.475 0.42670.4349 2.08 3.7e-02 pop0.3001 1.350 0.10410.1295 2.32 2.0e-02 pib -0.5485 0.578 0.20140.1799 -3.05 2.3e-03 pib2 -0.4033 0.668 0.14500.1152 -3.50 4.6e-04 crois -0.0541 0.947 0.02360.0227 -2.38 1.7e-02 instab.x1 1.9649 7.134 0.48390.4753 4.13 3.6e-05 instab.autres 0.8498 2.339 0.46930.4594 1.85 6.4e-02 Likelihood ratio test=78.3 on 7 df, p=3.04e-14 n= 7286 Estimate Std. Error t value Pr(|t|) sqrt(diag(cox1$naive.var)) 0.983970.02199 44.74 8.35e-09 *** Multiple R-Squared: 0.997, Adjusted R-squared: 0.9965 The naive standard errors (se) seem closer to the robust se than they were when not modeling for dependence. 0.98397 is very close to one, R^2 grew, etc. The dependence is high (risk is multiplied by 2.475 the day after an event) but conditional independence (given covariates) seems hard to reject. == = CASE 3 Finally, I compared these results with those without repeated events (which gives a smaller dataset). A country is removed as soon as we observe its first event. (robust se is still computed, even if naive se should in fact be used here to compute the pvalue) coxph(formula = Surv(start, stop, status) ~ cluster(ccode) + pop + pib + pib2 + crois + instab.x1 + instab.autres, data = xstep[no.previous.event, ]) coef exp(coef) se(coef) robust se z p pop0.4236 1.528 0.10300.1157 3.66 2.5e-04 pib -0.7821 0.457 0.20720.1931 -4.05 5.1e-05 pib2 -0.3069 0.736 0.14770.1254 -2.45 1.4e-02 crois -0.0432 0.958 0.02810.0258 -1.67 9.5e-02 instab.x1 1.9925 7.334 0.53210.3578 5.57 2.6e-08 instab.autres 1.3571 3.885 0.54280.5623 2.41 1.6e-02 Likelihood ratio test=66.7 on 6 df, p=1.99e-12 n=5971 (2466 observations deleted due to missing) summary(lm(sqrt(diag(cox1$var))~ sqrt(diag(cox1$naive.var)) -1)) Estimate Std. Error t value Pr(|t|) sqrt(diag(cox1$naive.var)) 0.866820.07826 11.08 0.000104 *** Residual standard error: 0.06328 on 5 degrees of freedom Multiple R-Squared: 0.9608, Adjusted R-squared: 0.953 There seems to be no evidence that robust se is more different from se in case 2 than in case 3 (and case 1). It even seems closer. I conclude that conditional independence (martingaleness) cannot be rejected in CASE 2, when modeling the dependence between events with a covariate. Mayeul KAUFFMANN Univ. Pierre Mendes France Grenoble - France Then, there is still another option. In fact, I already modelled explicitely the influence of past events with a proximity of last event covariate, assuming the dependence on the last event decreases at a constant rate (for instance, the proximity covariate varies from 1 to 0.5 in the first 10 years after an event, then from 0.5 to 0.25 in the next ten
RE: [R] elegant matrix creation
In reading this over I noticed that the two input vectors to shiftL could be put in a parameterized form that reduces the input to two integers: shiftL - function(v, shift) outer(shift,seq(along=v)-1, function(i,j)v[(i+j)%%length(v)+1]) jj - function(i,j) shiftL( seq(0,by=i,len=9) %% 9 %/% 3 + 1, seq(0, len=9, by=j) %% 9) jj1 - jj(1,3) jj2 - jj(5,6) jj3 - jj(1,4) From: Gabor Grothendieck [EMAIL PROTECTED] [Sorry if this gets posted twice but I am having more gmane posting problems.] Not sure if this qualifies as elegant or not but it does (1) bring all three matrices under a single scheme, (2) reduce the number of numbers from 81 to 18 per matrix, (3) requires only a single one line utility function, (4) is simple and (5) gives some minimal insight into the patterns. The key thing to note is that each row of each matrix is a cyclic shift of the first row of that matrix. Define a shift function which shifts its vector argument v by shift positions to the left creating a one row matrix. If shift is a vector it creates a matrix with one row per shift. shiftL - function(v, shift) outer(shift,seq(along=v)-1, function(i,j)v[(i+j)%%length(v)+1]) jj1a - shiftL(c(1,1,1,2,2,2,3,3,3),c(0,3,6,0,3,6,0,3,6)) jj2a - shiftL(c(1,2,1,3,1,3,2,3,2),c(0,6,3,0,6,3,0,6,3)) jj3a - shiftL(c(1,1,1,2,2,2,3,3,3),c(0,4,8,3,7,2,6,1,5)) # or finding expressions for the two args in each case: jj1b - shiftL( rep(1:3,c(3,3,3)), rep(c(0,3,6),3) ) jj2b - shiftL( c(shiftL(c(1,3,2),c(0,2,0))), rep(c(0,6,3),3) ) jj3b - shiftL( rep(1:3,c(3,3,3)), seq(0,32,4) %% 9 ) __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] covariate selection in cox model (counting process)
On Wed, 28 Jul 2004, Mayeul KAUFFMANN wrote: If you can get the conditional independence (martingaleness) then, yes, BIC is fine. One way to check might be to see how similar the standard errors are with and without the cluster(id) term. (Thank you again !, Thomas.) At first look, the values seemed very similar (see below, case 2). However, to check this without being too subjective, and without a specific test, I needed other values to assess the size of the differences: what is similar, what is not? I think the econometricians have theory for this (comparing the whole covariance matrices). -thomas == = CASE 1 I first estimated the model without modeling dependence: Call: coxph(formula = Surv(start, stop, status) ~ cluster(ccode) + pop + pib + pib2 + crois + instab.x1 + instab.autres, data = xstep) coef exp(coef) se(coef) robust se z p pop0.3606 1.434 0.09780.1182 3.05 2.3e-03 pib -0.5947 0.552 0.19520.1828 -3.25 1.1e-03 pib2 -0.4104 0.663 0.14520.1270 -3.23 1.2e-03 crois -0.0592 0.943 0.02450.0240 -2.46 1.4e-02 instab.x1 2.2059 9.079 0.46920.4097 5.38 7.3e-08 instab.autres 0.9550 2.599 0.47000.4936 1.93 5.3e-02 Likelihood ratio test=74 on 6 df, p=6.2e-14 n= 7286 There seems to be a strong linear relationship between standard errors (se, or naive se) and robust se. summary(lm(sqrt(diag(cox1$var))~ sqrt(diag(cox1$naive.var)) -1)) Coefficients: Estimate Std. Error t value Pr(|t|) sqrt(diag(cox1$naive.var)) 0.961030.04064 23.65 2.52e-06 *** Multiple R-Squared: 0.9911, Adjusted R-squared: 0.9894 == = CASE 2 Then I added a variable (pxcw) measuring the proximity of the previous event (1pxcw0) n= 7286 coef exp(coef) se(coef) robust se z p pxcw 0.9063 2.475 0.42670.4349 2.08 3.7e-02 pop0.3001 1.350 0.10410.1295 2.32 2.0e-02 pib -0.5485 0.578 0.20140.1799 -3.05 2.3e-03 pib2 -0.4033 0.668 0.14500.1152 -3.50 4.6e-04 crois -0.0541 0.947 0.02360.0227 -2.38 1.7e-02 instab.x1 1.9649 7.134 0.48390.4753 4.13 3.6e-05 instab.autres 0.8498 2.339 0.46930.4594 1.85 6.4e-02 Likelihood ratio test=78.3 on 7 df, p=3.04e-14 n= 7286 Estimate Std. Error t value Pr(|t|) sqrt(diag(cox1$naive.var)) 0.983970.02199 44.74 8.35e-09 *** Multiple R-Squared: 0.997, Adjusted R-squared: 0.9965 The naive standard errors (se) seem closer to the robust se than they were when not modeling for dependence. 0.98397 is very close to one, R^2 grew, etc. The dependence is high (risk is multiplied by 2.475 the day after an event) but conditional independence (given covariates) seems hard to reject. == = CASE 3 Finally, I compared these results with those without repeated events (which gives a smaller dataset). A country is removed as soon as we observe its first event. (robust se is still computed, even if naive se should in fact be used here to compute the pvalue) coxph(formula = Surv(start, stop, status) ~ cluster(ccode) + pop + pib + pib2 + crois + instab.x1 + instab.autres, data = xstep[no.previous.event, ]) coef exp(coef) se(coef) robust se z p pop0.4236 1.528 0.10300.1157 3.66 2.5e-04 pib -0.7821 0.457 0.20720.1931 -4.05 5.1e-05 pib2 -0.3069 0.736 0.14770.1254 -2.45 1.4e-02 crois -0.0432 0.958 0.02810.0258 -1.67 9.5e-02 instab.x1 1.9925 7.334 0.53210.3578 5.57 2.6e-08 instab.autres 1.3571 3.885 0.54280.5623 2.41 1.6e-02 Likelihood ratio test=66.7 on 6 df, p=1.99e-12 n=5971 (2466 observations deleted due to missing) summary(lm(sqrt(diag(cox1$var))~ sqrt(diag(cox1$naive.var)) -1)) Estimate Std. Error t value Pr(|t|) sqrt(diag(cox1$naive.var)) 0.866820.07826 11.08 0.000104 *** Residual standard error: 0.06328 on 5 degrees of freedom Multiple R-Squared: 0.9608, Adjusted R-squared: 0.953 There seems to be no evidence that robust se is more different from se in case 2 than in case 3 (and case 1). It even seems closer. I conclude that conditional independence (martingaleness) cannot be rejected in CASE 2, when modeling the dependence between events with a covariate. Mayeul KAUFFMANN Univ. Pierre Mendes France Grenoble - France Then, there is still another option. In fact, I already modelled explicitely the influence of past events with a
Re: [R] Integration with 'adapt'
Thank you for your answer to my question. Here is a reproducible example of the problem: A-2 B--1 C-1 linf-c(-10,-1,0.0003) lsup-c(10,1,0.0004) Integrand1-function(v) {exp(-1*(A*v[1]^2-B*v[1]+C+(((v[1]-v[2])^2)/(2*v[3]^2/(K*v[3])} Const-adapt(3,linf,lsup,functn=Integrand1)$value Warning message: Ifail=2, lenwrk was too small. -- fix adapt() ! Check the returned relerr! in: adapt(3, linf, lsup, functn = Integrand1) The problem is related to the small range of the third variable, (0.0003-0.0004), so I tried the change of variables z2-z/(max(z)-min(z)): rang=0.0004-0.0003 linf[3]-linf[3]/rang lsup[3]-lsup[3]/rang Integrand1-function(v) {exp(-1*(A*v[1]^2-B*v[1]+C+(((v[1]-v[2])^2)/(2*(rang^2)*(v[3]^2)/(K*v[3])} #With the change of variables, the constant rang appears in the function. Const-adapt(3,linf,lsup,functn=Integrand1)$value Warning message: Ifail=2, lenwrk was too small. -- fix adapt() ! Check the returned relerr! in: adapt(3, linf, lsup, functn = Integrand1) And the problem persists, probably because of the small value of rang, what can be checked below: rang-.01 Integrand1-function(v) {exp(-1*(A*v[1]^2-B*v[1]+C+(((v[1]-v[2])^2)/(2*(rang^2)*(v[3]^2)/(K*v[3])} #The same function as before, with a bigger value for rang Const-adapt(3,linf,lsup,functn=Integrand1)$value Const [1] 8.766637e-05 How can I evaluate the integral, even with small values of rang? Thanks a lot for your help. Best regards. Rodrigo Drummond Rodrigo == Rodrigo Drummond [EMAIL PROTECTED] on Tue, 27 Jul 2004 15:31:07 -0300 (BRT) writes: Rodrigo Hi all, I need to calculate a multidimensional Rodrigo integration on R. I am using the command Rodrigo adapt (from library adapt), although it's a package, not a library. Rodrigo sometimes I get the following error message: Rodrigo Ifail=2, lenwrk was too small. -- fix adapt() ! Rodrigo Check the returned relerr! in: adapt(3, linf, lsup, Rodrigo functn = Integrando1) If you could give as a *reproducible* example, we (the adapt authors) would have chance to do what the above message says, namely fix adapt() .. Rodrigo I guess it happens because the domain of Rodrigo integration is too small, maybe, maybe not. We need an example we can reproduce, see above. Rodrigo although I tried a change of variables to avoid Rodrigo this problem and it didnt help. The command adapt Rodrigo calls a fortran routine, but I dont know fortran Rodrigo enough to fix the problem. Martin Maechler, ETH Zurich Rodrigo D. Drummond Laboratorio Genoma Funcional Centro de Biologia Molecular e Eng. Genetica Universidade Estadual de Campinas Caixa Postal 6010 13083-875 - Campinas - SP - Brasil Tel: xx-19-3788-1119 Fax: xx-19-3788-1089 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] elegant matrix creation
Here is yet one more simplification. Note that seq9 in this one is not the same as seq9 in the last one. This one looks much more symmetric in the arguments and is shorter: seq9 - function(i) seq(0, len = 9, by = i) jj - function(i,j) outer( seq9(i), seq9(j), + ) %% 9 %/% 3 + 1 jj1. - jj(3,1) jj2. - jj(3,5) jj3. - jj(4,1) __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Gaussian frailty leads to segmentation fault
We really need a reproducible example to find segmentation faults. Can you make one? -thomas On Wed, 28 Jul 2004, Christian Lederer wrote: Dear R gurus, for a simulation concerning study effects and historical controls in survival analysis, i would like to experiment with a gaussian frailty model. The simulated scenario consists of a randomized trial (treatment and placebo) and historical controls (only placebo). So the simulated data frames consist of four columns $time, $cens, $study, $treat. $time, $cens are the usual survival data. For the binary thretment indicator we have $treat == 0 or 1, if $study == 1, $treat == 1 if $study 1 Typical parameters for my simulations are: sample sizes (per arm): between 100 and 200 number of historical studies: between 7 and 15 hazard ratio treatment/placebo: between 0.7 and 1 variance of the study effekt: between 0 and 0.3 Depending on the sample sizes, the following call sometimes leads to a segmentation fault: coxph(Surv(time,cens) ~ as.factor(treatment) + frailty(study, distribution=gaussian), data=data) I noticed, that this segmentation fault occures most frequently, if the number of randomized treatment patients is higher than the number of randomized placebo patients, and the number of historical studies is large. There seems to be no problem, if there are at least as many randomized placebo patients as treated patients. Unfortunately, this is not the situation i want to investigate (historical controls should be used to decrease the number of treated patients). Is there a way to circumwent this problem? Christian P.S. Is it allowed, to attach gzipped sample data sets in this mailing list? __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Thomas Lumley Assoc. Professor, Biostatistics [EMAIL PROTECTED] University of Washington, Seattle __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] How to add an object to an RData file ?
I guess it would be nice if save() has an `append' option... Andy From: Roger D. Peng I can't think of a faster way. -roger Ernesto Jardim wrote: Hi, I've saved an RData file with save and now I want to add a new object to this file. At the moment I do: attach(file.RData) save(list=c(new,obj, ls(pos=2)), file=file.RData, compress=T) detach() Is there a quicker method that just add the object to the file ? Thanks EJ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] a bug with LAPACK ? non orthogonal vectors obtained with eigen of a symmetric matrix
At 11:25 28/07/2004, Roger D. Peng wrote: This is interesting. I can reproduce your results but cannot come up with an explanation. However, using svd(LINPACK = FALSE) seems to work every time. Might you consider trying that instead? -roger I have try to do it with svd but here the problem is that the matrices under study have some equal singular values and in this case, I do not have u=v for a symmetric matrix (I supose that a linear combination of u corresponding to equal singular values allow to obtain the v for the same singular values). i.e. if d2=d3=d4 u2=av2+bv3+cv4 ... and not necessarily u2=v2 But here, I found another problem: #This is for svd with LINPACK. It works, u vectors are orthogonal (same for v) or(i in 3:20){ x=dbcenter(nb2mat(cell2nb(i,i),style=B)) res=svd(x,LIN=T) eq0 - apply(as.matrix(res$d),1,function(x) identical(all.equal(x, 0), TRUE)) resu=res$u[,-which(eq0)] resv=res$v[,-which(eq0)] print(paste(Grid size,i)) print(paste(u ortho ?,all.equal(diag(1,ncol(resu)),t(resu)%*%resu))) print(paste(v ortho ?,all.equal(diag(1,ncol(resv)),t(resv)%*%resv))) print(paste(u == v ?,all.equal(abs(resv),abs(resu } #But without LINPACK: for(i in 3:20){ x=dbcenter(nb2mat(cell2nb(i,i),style=B)) res=svd(x,LIN=F) eq0 - apply(as.matrix(res$d),1,function(x) identical(all.equal(x, 0), TRUE)) resu=res$u[,-which(eq0)] resv=res$v[,-which(eq0)] print(paste(Grid size,i)) print(paste(u ortho ?,all.equal(diag(1,ncol(resu)),t(resu)%*%resu))) print(paste(v ortho ?,all.equal(diag(1,ncol(resv)),t(resv)%*%resv))) print(paste(u == v ?,all.equal(abs(resv),abs(resu } All is fine, except for i=14 where vectors u (and v) are not orthogonal [1] Grid size 14 [1] u ortho ? Mean relative difference: 269.2564 [1] v ortho ? Mean relative difference: 2097.224 [1] u == v ? Mean relative difference: 0.7612647 I wonder if it is a general problem, a problem due to the structure of my matrix (I don't think so) , a windoze problem or a bug in my head !!! Sincerely. Stephane DRAY wrote: Hello, I have send send this message one week ago but I have receive no answer. Perhaps, some of RListers were in holidays and do not read my message. I try again.. My problem is that I obtained non orthonormal eigenvectors with some matrices with LAPACK while EISPACK seems to provide good results. Is there some restrictions with the use of LAPACK ? Is it a bug ? I did not find the answer. Here is my experiment: I have obtained strange results using eigen on a symmetric matrix: # this function perform a double centering of a matrix (xij-rowmean(i)-colmean(j)+meantot) dbcenter=function(mat){ rmean=apply(mat,1,mean) cmean=apply(mat,2,mean) newmat=sweep(mat,1,rmean,-) newmat=sweep(newmat,2,cmean,-) newmat=newmat+mean(mat) newmat} # i use spdep package to create a spatial contiguity matrix library(spdep) x=dbcenter(nb2mat(cell2nb(3,3),style=B)) #compute eigenvalues of a 9 by 9 matrix res=eigen(x) # some eigenvalues are equal to 0 eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) # I remove the corresponding eigenvectors res0=res$vec[,-which(eq0)] # then I compute the Froebenius norm with the identity matrix sum((diag(1,ncol(res0))-crossprod(res0))^2) [1] 1.515139e-30 # The results are correct, # then I do it again with a biggest matrix(100 by 100) x=dbcenter(nb2mat(cell2nb(10,10),style=B)) res=eigen(x) eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) res0=res$vec[,-which(eq0)] sum((diag(1,ncol(res0))-crossprod(res0))^2) [1] 3.986387 I have try the same with res=eigen(x,EISPACK=T): x=dbcenter(nb2mat(cell2nb(10,10),style=B)) res=eigen(x,EISPACK=T) eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) res0=res$vec[,-which(eq0)] sum((diag(1,ncol(res0))-crossprod(res0))^2) [1] 1.315542e-27 So I wonder I there is a bug in the LAPACK algorithm or if there are some things that I have not understood. Note that my matrix has some pairs of equal eigenvalues. Thanks in advance. I have continue my experiments in changing the size of my matrix : (3^2 by 3^2, 4^2 by 4^2... 20^2 by 20^2) EISPACK is always correct but LINPACK provide very strange results: for(i in 3:20){ + x=dbcenter(nb2mat(cell2nb(i,i),style=B)) + res=eigen(x,EIS=T) + eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) + res0=res$vec[,-which(eq0)] + print(sum((diag(1,ncol(res0))-crossprod(res0))^2)) + } [1] 7.939371e-30 [1] 2.268788e-29 [1] 9.237286e-29 [1] 1.806393e-28 [1] 3.24619e-28 [1] 5.239195e-28 [1] 9.78079e-28 [1] 1.315542e-27 [1] 1.838600e-27 [1] 3.114150e-27 [1] 5.499297e-27 [1] 5.471782e-27 [1] 1.075098e-26 [1] 1.534822e-26 [1] 1.771326e-26 [1] 2.342404e-26 [1] 3.462522e-26 [1] 4.310143e-26 for(i in 3:20){ + x=dbcenter(nb2mat(cell2nb(i,i),style=B)) + res=eigen(x) + eq0 - apply(as.matrix(res$values),1,function(x) identical(all.equal(x, 0), TRUE)) + res0=res$vec[,-which(eq0)] +
RE: [R] How to add an object to an RData file ?
On Wed, 28 Jul 2004, Liaw, Andy wrote: I guess it would be nice if save() has an `append' option... It would still have to read it in and out again, and then we would have all the issues about multiple objects with the same name. -thomas Andy From: Roger D. Peng I can't think of a faster way. -roger Ernesto Jardim wrote: Hi, I've saved an RData file with save and now I want to add a new object to this file. At the moment I do: attach(file.RData) save(list=c(new,obj, ls(pos=2)), file=file.RData, compress=T) detach() Is there a quicker method that just add the object to the file ? Thanks EJ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Thomas Lumley Assoc. Professor, Biostatistics [EMAIL PROTECTED] University of Washington, Seattle __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] R/S-plus Course***In Atlanta, GA***R/Splus Fundamentals and Programming Techniques, August 19-20, 2004
XLSolutions Corporation (www.xlsolutions-corp.com) is proud to announce June 2004 2-day R/S-plus Fundamentals and Programming Techniques. Atlanta, GA -- August, 19-20 Interested in our R/Splus Advanced Programming course? Please email us! Reserve your seat now at the early bird rates! Payment due AFTER the class. Course Description: This two-day beginner to intermediate R/S-plus course focuses on a broad spectrum of topics, from reading raw data to a comparison of R and S. We will learn the essentials of data manipulation, graphical visualization and R/S-plus programming. We will explore statistical data analysis tools,including graphics with data sets. How to enhance your plots. We will perform basic statistics and fit linear regression models. Participants are encouraged to bring data for interactive sessions With the following outline: - An Overview of R and S - Data Manipulation and Graphics - Using Lattice Graphics - A Comparison of R and S-Plus - How can R Complement SAS? - Writing Functions - Avoiding Loops - Vectorization - Statistical Modeling - Project Management - Techniques for Effective use of R and S - Enhancing Plots - Using High-level Plotting Functions - Building and Distributing Packages (libraries) Email us for group discounts. Email Sue Turner: [EMAIL PROTECTED] Phone: 206-686-1578 Visit us: www.xlsolutions-corp.com/training.htm Please let us know if you and your colleagues are interested in this classto take advantage of group discount. Register now to secure your seat! Interested in R/Splus Advanced course? email us. Cheers, Elvis Miller, PhD Manager Training. XLSolutions Corporation 206 686 1578 www.xlsolutions-corp.com [EMAIL PROTECTED] __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] How to add an object to an RData file ?
I don't think it is that simple. There are several save formats, some historical that no one wants to touch. So suppose we confined it to the current format. Then one would need a pass over the saved file to see what (e.g. namespaces) should be shared and to find the existing references, so I don't think it would actually save anything much. (If you want to work out if shared can be skipped, the code is in serialize.c ) On Wed, 28 Jul 2004, Liaw, Andy wrote: I guess it would be nice if save() has an `append' option... Andy From: Roger D. Peng I can't think of a faster way. -roger Ernesto Jardim wrote: Hi, I've saved an RData file with save and now I want to add a new object to this file. At the moment I do: attach(file.RData) save(list=c(new,obj, ls(pos=2)), file=file.RData, compress=T) detach() Is there a quicker method that just add the object to the file ? Thanks EJ __ [EMAIL PROTECTED] mailing list https://www.stat.math.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 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] lattice variable by page
Hi, Using lattice's xyplot, is it possible to specify a variable to group plots by page. For example, if I have xyplot(y~x|A*B*C) could I get a page created for each unique value in variable C ? I am hoping to avoid having pages with the same strip above each plot on a page. Thanks, M Matt Pocernich NCAR - Research Applications Program 303-497-8312 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] where is average shifted histogram?
Hello! In the book Modern Applied Statistics with S (4th ed), section 5.6 the concept of the average shifted histogram or ASH is mentionend. Also it is mentioned in the same section The code used is in the scripts for this chapter (from figure caption 5.8, analysis of the geyser duration data). *However*, I have trouble finding the code for that function! Admittedly, I am a newby to R. Where exactly do I find the scripts to that chapter? Commands like library(MASS), or help.search(average shifted histogram) did not help. Secondly: Is the fact that I have trouble finding a routine for the average shifted histogram method a hint, that the method is outdated, and instead method X or Y should be used? The problem I am trying to solve is to smoothe a two-dimensional data set, the y-values change very rapidly between one and zero. The function loess in R seems to do a very reasonable job for that problem, but I would like to use a second and different method. Any help is highly appreciated! Thanks, Eckart Bindewald __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] lattice variable by page
Matt Pocernich wrote: Hi, Using lattice's xyplot, is it possible to specify a variable to group plots by page. For example, if I have xyplot(y~x|A*B*C) could I get a page created for each unique value in variable C ? I am hoping to avoid having pages with the same strip above each plot on a page. Thanks, M Others may have a better solution, but when I've done this in the past I used a for loop with subset: # assuming C is a factor to use nlevels for(ci in nlevels(C)) { pl - xyplot(y ~ x | A * B, subset = C == ci, main = paste(C ==, ci)) print(pl) } HTH, --sundar __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] lattice variable by page
On Wednesday 28 July 2004 16:36, Matt Pocernich wrote: Hi, Using lattice's xyplot, is it possible to specify a variable to group plots by page. For example, if I have xyplot(y~x|A*B*C) could I get a page created for each unique value in variable C ? I am hoping to avoid having pages with the same strip above each plot on a page. If I understand you correctly (which I'm not sure I do), then no, not directly. But you can do the following (assuming C is a factor): for (l in levels(C)) { print(xyplot(y ~ x | A * B, subset = C == l)) } Is that what you are looking for? Deepayan __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Transparent backgrounds in png files
The result I'm aiming to achieve is a bitmap that can be imported into a PowerPoint file that shows what's behind the lines of the plot. There's a way in PowerPoint that almost works. By choosing a colour to set as transparent, what is behind the graphic is indeed visible, but it's at the expense of losing line and text definition. I notice there have been discussions about transparent backgrounds mostly with lattice plots. The problem most people seem to have was getting a blue background when they wanted a white one. Mine is the reverse (and I'm using standard graphics, not Lattice). I'm trying to get a transparent background but it always comes out white. Setting bg = transparent, I've tried using a bitmap device to create a png file. I've also tried creating a postscript file and converting it to a PNG file using the Gimp. I've always used a resolution of 300 dpi in bitmaps since the default is far too low. Other ideas appreciated. version _ platform i686-pc-linux-gnu arch i686 os linux-gnu system i686, linux-gnu status major1 minor9.1 year 2004 month06 day 21 language R Thanks. -- Patrick Connolly HortResearch Mt Albert Auckland New Zealand Ph: +64-9 815 4200 x 7188 ~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~ I have the world`s largest collection of seashells. I keep it on all the beaches of the world ... Perhaps you`ve seen it. ---Steven Wright ~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] lattice variable by page
Sundar Dorai-Raj wrote: Matt Pocernich wrote: Hi, Using lattice's xyplot, is it possible to specify a variable to group plots by page. For example, if I have xyplot(y~x|A*B*C) could I get a page created for each unique value in variable C ? I am hoping to avoid having pages with the same strip above each plot on a page. Thanks, M Others may have a better solution, but when I've done this in the past I used a for loop with subset: # assuming C is a factor to use nlevels for(ci in nlevels(C)) { pl - xyplot(y ~ x | A * B, subset = C == ci, main = paste(C ==, ci)) print(pl) } HTH, --sundar Sorry made a typo: nlevels should be levels. --sundar __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] R with perl
Just a quick, slightly off-topic, question: Could people using RSperl and/or Statistics::R comment on the relative merits of each? I would like to work mainly from perl, using R as a computation engine. I have gotten Statistics::R to work, but the interface is entirely string-based, it seems, which is a drawback and will probably require temporary files for communication. I haven't gotten RSperl to work from perl (on MacOS 10.3.4), but I have to admit I haven't gotten too deep into things. Any insight is appreciated. Thanks, Sean [[alternative HTML version deleted]] __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Gaussian frailty leads to segmentation fault
Dear Thomas, attached you find a data frame which produces the error. I am using survival 2.11-5 under R 1.9.1-1 and 1.9.0-1. By the way, if i randomly omit 50% of the data, i usually get no crash, but a warning message like this: Inner loop failed to coverge for iterations 1 2 3 in: coxpenal.fit(X, Y, strats, offset, init = init, control, weights = weights, Maybe, the model is not appropriate for this kind of data. But on the other hand, as soon the treatment group (study == 1, treatment == 1) is smaller than the randomized placebo group (study == 1, treatment == 0), the warnings disappear. and the model gives reasonable results in my first simulations with normally distributed study effects. Christian Thomas Lumley wrote: We really need a reproducible example to find segmentation faults. Can you make one? -thomas On Wed, 28 Jul 2004, Christian Lederer wrote: Dear R gurus, for a simulation concerning study effects and historical controls in survival analysis, i would like to experiment with a gaussian frailty model. The simulated scenario consists of a randomized trial (treatment and placebo) and historical controls (only placebo). So the simulated data frames consist of four columns $time, $cens, $study, $treat. $time, $cens are the usual survival data. For the binary thretment indicator we have $treat == 0 or 1, if $study == 1, $treat == 1 if $study 1 Typical parameters for my simulations are: sample sizes (per arm): between 100 and 200 number of historical studies: between 7 and 15 hazard ratio treatment/placebo: between 0.7 and 1 variance of the study effekt: between 0 and 0.3 Depending on the sample sizes, the following call sometimes leads to a segmentation fault: coxph(Surv(time,cens) ~ as.factor(treatment) + frailty(study, distribution=gaussian), data=data) I noticed, that this segmentation fault occures most frequently, if the number of randomized treatment patients is higher than the number of randomized placebo patients, and the number of historical studies is large. There seems to be no problem, if there are at least as many randomized placebo patients as treated patients. Unfortunately, this is not the situation i want to investigate (historical controls should be used to decrease the number of treated patients). Is there a way to circumwent this problem? Christian P.S. Is it allowed, to attach gzipped sample data sets in this mailing list? __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Thomas Lumley Assoc. Professor, Biostatistics [EMAIL PROTECTED] University of Washington, Seattle __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Re: [BioC] normalisation for universal reference in 2 channel arrays
Thanks for the suggestions, one option i was thinking was to center the universal channel and sample channels by normalising the medians of each column of the matrix on each channel. But I don't know if it is appropriate to do a loess after this? Does the quantile option actually do this? Peter At 06:59 PM 7/28/2004, you wrote: At 06:02 AM 29/07/2004, Peter Wilkinson wrote: I would like to know some alternative to normalization for 2 channel experiments against universal, where samples may have been hybridized in seperate batches where the universal RNA has changed lot (should not have happened but it did). What is the same between the batches is that what looks like up-regulated compared to the universal iis upr-egulated, and what is down-regulated looks down-regulated. The difference is that the down-regulated (and not in the up), so so much more down-regulated in one of the batches. It looks like to be that the universal has more mRNA abundance in one batch over the other. Gquantile won't help because it doesn't change the M-values. (It is intended for use with single channel analyses.) You could try 'quantile' or 'scale' normalization (not both) but there are no magic bullets in a situation like this. If you use 'scale normalization, you should always do within-array normalization first. If you use 'quantile' normalization, the within-array normalization step is optional. Gordon so ... I have samples that can be divided into 2 classes: 0,1, and within each class I have samples that have been run at different times. I would like to treat my universal channel uniformly across all samples (assuming that my universal changed lot), and then adjust the Sample (Red) channel to that. is the normalizeBetweenArrays with the method=Gquantile option the right option for this? What is the complete work-flow for this case? And after I have normalized within the Arrays, can I go on to scale option for normalizing between arrays. Peter __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] Parallel Functions on AIX
From: Brian Ripley The first thing to ascertain is that R will actually build on such a machine -- there have been a lot of reports of failure on AIX, and no recent reports of success. According to our AIX sysadmin, R-1.9.1 actually compiles a lot more cleanly than R-1.8.1 and prior. I take it the recent reports of difficulties has to do with (possibly outdated) GCC. We have no problem with xlc/xlf. The trouble we're having is getting R to link against readline in 64-bit. Haven't figure out how to do that yet. Life w/o readline is rather trying... Andy If it does, R itself is single-threaded (but can make use of multi-threaded BLAS) but packages such as Rmpi, snow, RScaLAPACK provide parallel facilities (and are in the FAQ). On Wed, 28 Jul 2004, Liao, Kexiao wrote: Dear R Development Team, You actually wrote to R-help! Does the latest version R-1.9.1 provide parallel functions if R is running on Multiple CPUs Unix platform (IBM AIX e-server)? -- 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 __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] RE: [S] tree function in R language
1. Could it be that your computer is behind a firewall? If so, try reading the R for Windows FAQ. 2. Please ask R-related question on R-help instead of S-news. Andy From: cheng wu Hi, Andy Thank you for your answer. Why I can't load CRAN packages? the error message is: {a - CRAN.packages() + install.packages(select.list(a[,1],,TRUE), .libPaths()[1], available=a)} trying URL `http://cran.r-project.org/bin/windows/contrib/PACKAGES' unable to connect to 'cran.r-project.org'. Error in download.file(url = paste(contriburl, PACKAGES, sep = /), : cannot open URL `http://cran.r-project.org/bin/windows/contrib/PACKAGES' From: Chushu Gu [EMAIL PROTECTED] To: cheng wu [EMAIL PROTECTED] Subject: Fw: [S] tree function in R language Date: Wed, 28 Jul 2004 09:14:48 -0400 - Original Message - From: Liaw, Andy [EMAIL PROTECTED] To: 'chushu Gu' [EMAIL PROTECTED]; [EMAIL PROTECTED] Sent: Tuesday, July 27, 2004 11:22 PM Subject: RE: [S] tree function in R language Have you read the (latest edition of the) book for which the package you are using supports? There are differences in S-PLUS and R (and the 4th edition of MASS supports both, thus ought to tell you this particular difference between the two). tree() in S-PLUS is written originally by Clark and Pregibon. If you want that functionality in R, you need to load the `tree' package (available on CRAN), which is an independent implementation by one of the co-authors of MASS. Another hint: Look in the `scripts' subdirectory of where the `MASS' package is installed. Andy From: chushu Gu Hi all, I am using R 1.8.1, I have the following code, library(MASS) data(iris) ir.tr - tree(Species ~., iris) ir.tr summary(ir.tr) I got the following message: Error: couldn't find function tree I don't the reason, as I already load the library MASS. Could anyone tell me the possible reasons? Thanks, Chushu Gu _ Share a single photo or an entire photo slide show right inside MSN(r) Messenger. http://join.msn.com/?pgmarket=en-capage=byoa/premxAPID=1994; DI=1034SU=http://hotmail.com/encaHL=Market_MSNIS_Taglines Start enjoying all the benefits of MSN(r) Premium right now and get the first two months FREE*. This message was distributed by [EMAIL PROTECTED] To unsubscribe send e-mail to [EMAIL PROTECTED] with -- Notice: This e-mail message, together with any attachments, contains information of Merck Co., Inc. (One Merck Drive, Whitehouse Station, New Jersey, USA 08889), and/or its affiliates (which may be known outside the United States as Merck Frosst, Merck Sharp Dohme or MSD and in Japan, as Banyu) that may be confidential, proprietary copyrighted and/or legally privileged. It is intended solely for the use of the individual or entity named on this message. If you are not the intended recipient, and have received this message in error, please notify us immediately by reply e-mail and then delete it from your system. -- __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] RE: [S] tree function in R language
Shouldn't the URL be (for R 1.8.1 on Windows): http://cran.r-project.org/bin/windows/contrib/1.8/PACKAGES There is no URL as listed below, which is presumably why the error message. Was options()$CRAN changed improperly or is there some other Windows specific issue that is escaping me at the moment? BTW, you should upgrade to R 1.9.1, as you are two versions behind at this point. HTH, Marc Schwartz On Wed, 2004-07-28 at 23:08, Liaw, Andy wrote: 1. Could it be that your computer is behind a firewall? If so, try reading the R for Windows FAQ. 2. Please ask R-related question on R-help instead of S-news. Andy From: cheng wu Hi, Andy Thank you for your answer. Why I can't load CRAN packages? the error message is: {a - CRAN.packages() + install.packages(select.list(a[,1],,TRUE), .libPaths()[1], available=a)} trying URL `http://cran.r-project.org/bin/windows/contrib/PACKAGES' unable to connect to 'cran.r-project.org'. Error in download.file(url = paste(contriburl, PACKAGES, sep = /), : cannot open URL `http://cran.r-project.org/bin/windows/contrib/PACKAGES' From: Chushu Gu [EMAIL PROTECTED] To: cheng wu [EMAIL PROTECTED] Subject: Fw: [S] tree function in R language Date: Wed, 28 Jul 2004 09:14:48 -0400 - Original Message - From: Liaw, Andy [EMAIL PROTECTED] To: 'chushu Gu' [EMAIL PROTECTED]; [EMAIL PROTECTED] Sent: Tuesday, July 27, 2004 11:22 PM Subject: RE: [S] tree function in R language Have you read the (latest edition of the) book for which the package you are using supports? There are differences in S-PLUS and R (and the 4th edition of MASS supports both, thus ought to tell you this particular difference between the two). tree() in S-PLUS is written originally by Clark and Pregibon. If you want that functionality in R, you need to load the `tree' package (available on CRAN), which is an independent implementation by one of the co-authors of MASS. Another hint: Look in the `scripts' subdirectory of where the `MASS' package is installed. Andy From: chushu Gu Hi all, I am using R 1.8.1, I have the following code, library(MASS) data(iris) ir.tr - tree(Species ~., iris) ir.tr summary(ir.tr) I got the following message: Error: couldn't find function tree I don't the reason, as I already load the library MASS. Could anyone tell me the possible reasons? Thanks, Chushu Gu __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] attach files on start up
Hi, I was wondering if anyone could please tell me if there's a way to attach the same files each time R starts up? Thanks! Kathryn Jones __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html