Re: [R] For each element in vector do...
Hi Or mabe easier in this case b - a+1-(a==0) HTH Petr On 17 Jan 2006 at 5:31, tom wright wrote: From: tom wright [EMAIL PROTECTED] To: Andrej Kastrin [EMAIL PROTECTED] Date sent: Tue, 17 Jan 2006 05:31:53 -0500 Copies to: r-help r-help@stat.math.ethz.ch Subject:Re: [R] For each element in vector do... a-c(0,1,2,3,0,4,5) b-vector(length=length(a)) b[a0]-a[a0]+1 b[a=0]-a[a=0] b [1] 0 2 3 4 0 5 6 On Tue, 2006-17-01 at 15:30 +0100, Andrej Kastrin wrote: Dear R useRs, I have a vector with positive and negative numbers: A=c(0,1,2,3,0,4,5) Now if i-th element in vector A is 0, then i-th element in vector B is a+1 else i-th element in vector b=a (or 0) vector A: 0 1 2 3 0 4 5 vector B: 0 2 3 4 0 5 6 What's the right way to do this. I still have some problems with for and if statements... Cheers, Andrej __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Petr Pikal [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Data frame index?
try: DF2 - as.data.frame(matrix(vec, nr=nrow(DF),nc=ncol(DF))== matrix(1:ncol(DF),nr=nrow(DF),nc=ncol(DF),byrow=T)) DF3 - data.frame(mapply(function(z,x,y) { x[y] - 0 ; x }, names(DF), DF, DF2, SIMPLIFY=F)) but there must be an easier way... Kenneth Cabrera a écrit : Hi, R users: I have a data.frame (not a matrix), I got a vector with the same length as the number of records (rows) of the data frame, and each element of that vector is the column number (in a specific range of columns) of the corresponding record that I must set to zero. How can I do this without a for loop? Thank you for your help. Kenneth __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Data frame index?
Hi eg. your data frame has 35 rows and 6 columns a-sample(1:6, 35, replace=T) b-1:35 vec-rep(0,35*6) vec[a+6*(b-1)]-1 This shall do the replacement your.d.f[matrix(vec,35,6, byrow=T)==1] - 0 But I am not sure if it is quicker than a loop. HTH Petr On 18 Jan 2006 at 2:35, Kenneth Cabrera wrote: Date sent: Wed, 18 Jan 2006 02:35:35 -0500 From: Kenneth Cabrera [EMAIL PROTECTED] To: r-help@stat.math.ethz.ch Subject:[R] Data frame index? Hi, R users: I have a data.frame (not a matrix), I got a vector with the same length as the number of records (rows) of the data frame, and each element of that vector is the column number (in a specific range of columns) of the corresponding record that I must set to zero. How can I do this without a for loop? Thank you for your help. Kenneth Petr Pikal [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] compiling R on powerpc-ibm-aix5.1.0.0
r-help, I am trying to compile R on an powerpc-IBM-AIX5.1.0.0 machine, Is R suitable to be used in this system? The R Installation and Administration document mentioned rs6000-ibm-aix not powerpc-IBM-AIX5.1.0.0 . When I tried to compile R.2.2.0 in powerpc-IBM-AIX5.1.0.0 using the following steps: (1) ./configure There is an error message :configure: error: --with-readline=yes (default) and headers/libs are not available So I used ./configure --with-readline=no --prefix=/export/home/hypan/R/bin instead, and configure is through. The following is information when finished configure R is now configured for powerpc-ibm-aix5.1.0.0 Source directory: . Installation directory:/export/home/hypan/R/bin C compiler:gcc -mno-fp-in-toc -g -O2 C++ compiler: g++ -g -O2 Fortran compiler: g77 -g -O2 Interfaces supported: X11, tcltk External libraries: Additional capabilities: PNG, JPEG, MBCS, NLS Options enabled: R profiling Recommended packages: yes configure: WARNING: you cannot build info or html versions of the R manuals configure: WARNING: I could not determine a browser configure: WARNING: I could not determine a PDF viewer *** Warning: the GNU linker, at least up to release 2.9.1, is reported *** to be unable to reliably create shared libraries on AIX. *** Therefore, libtool is disabling shared libraries support. If you *** really care for shared libraries, you may want to modify your PATH *** so that a non-GNU linker is found, and then restart. (2) make make is failed, the error message is: gcc -Wl,-bdynamic -Wl,-bE:../../etc/R.exp -Wl,-bM:SRE -L/usr/local/lib -o R.bin Rmain.o CConverters.o CommandLineArgs.o Rdynload.o Renviron.o RNG.o apply.o arithmetic.o apse.o array.o attrib.o base.o bind.o builtin.o character.o coerce.o colors.o complex.o connections.o context.o cov.o cum.o dcf.o datetime.o debug.o deparse.o deriv.o dotcode.o dounzip.o dstruct.o duplicate.o engine.o envir.o errors.o eval.o format.o fourier.o gevents.o gram.o gram-ex.o graphics.o identical.o internet.o iosupport.o lapack.o list.o logic.o main.o mapply.o match.o memory.o model.o names.o objects.o optim.o optimize.o options.o par.o paste.o pcre.o platform.o plot.o plot3d.o plotmath.o print.o printarray.o printvector.o printutils.o qsort.o random.o regex.o registration.o relop.o saveload.o scan.o seq.o serialize.o size.o sort.o source.o split.o sprintf.o startup.o subassign.o subscript.o subset.o summary.o sysutils.o unique.o util.o version.o vfonts.o xxxpr.o ../unix/libunix.a ../appl/libap! pl.a ../nmath/libnmath.a -lg2c -lm -lgcc_s /usr/local/lib/gcc-lib/powerpc-ibm-aix5.1.0.0/3.3/libgcc.a -lg /lib/crt0.o ../extra/zlib/libz.a ../extra/bzip2/libbz2.a ../extra/pcre/libpcre.a ../extra/intl/libintl.a -ldl -lm -lc -liconv /usr/local/lib/gcc-lib/powerpc-ibm-aix5.1.0.0/3.3/../../../../powerpc-ibm-aix5.1.0.0/bin/ld: target dynamic not found collect2: ld returned 1 exit status make: 1254-004 The error code from the last command is 1. Stop. make: 1254-004 The error code from the last command is 2. Stop. make: 1254-004 The error code from the last command is 1. Stop. make: 1254-004 The error code from the last command is 1. Stop. your help will be greatly appreciated. Thanks! Haiyan = = = = = = = = = = = = = = = = = = = = Haiyan Pan [EMAIL PROTECTED] Tel: 021-64363311-123 Shanghai Center for Bioinformatics Technology Floor 12th,100# QinZhou Road Shanghai,China,200235 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Data frame index?
you could try something like the following: dat - data.frame(matrix(rnorm(200), 20, 10)) index - sample(10, 20, TRUE) ### mat.ind - matrix(FALSE, nrow(dat), length(dat)) mat.ind[cbind(seq(along = index), index)] - TRUE dat[mat.ind] - 0 index dat I hope it helps. Best, Dimitris Dimitris Rizopoulos Ph.D. Student Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven, Belgium Tel: +32/(0)16/336899 Fax: +32/(0)16/337015 Web: http://www.med.kuleuven.be/biostat/ http://www.student.kuleuven.be/~m0390867/dimitris.htm - Original Message - From: Kenneth Cabrera [EMAIL PROTECTED] To: r-help@stat.math.ethz.ch Sent: Wednesday, January 18, 2006 8:35 AM Subject: [R] Data frame index? Hi, R users: I have a data.frame (not a matrix), I got a vector with the same length as the number of records (rows) of the data frame, and each element of that vector is the column number (in a specific range of columns) of the corresponding record that I must set to zero. How can I do this without a for loop? Thank you for your help. Kenneth __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Own Color Palette
Own Color Palette Dear all,?xml:namespace prefix = o ns = urn:schemas-microsoft-com:office:office /o:p/o:p I would like to generate a contour-plot according to a master plot. The problem is that the rainbow-palette included in R does not answer this purpose. I need a darker blue, no turquoise, relatively less green, more yellow and more red. Haw can I adjust the rainbow? Alternatively: How can I generate my own palette with at least 100 colors with smooth transitions?o:p/o:p Can anybody help me?o:p/o:p Thanks a loto:p/o:p Roberto:p/o:p o:p/o:p o:p/o:p [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Own Color Palette
On Wed, 18 Jan 2006, Robert Michael Rausch wrote: ?colorRampPalette Own Color Palette Dear all,?xml:namespace prefix = o ns = urn:schemas-microsoft-com:office:office /o:p/o:p I would like to generate a contour-plot according to a master plot. The problem is that the rainbow-palette included in R does not answer this purpose. I need a darker blue, no turquoise, relatively less green, more yellow and more red. Haw can I adjust the rainbow? Alternatively: How can I generate my own palette with at least 100 colors with smooth transitions?o:p/o:p Can anybody help me?o:p/o:p Thanks a loto:p/o:p Roberto:p/o:p o:p/o:p o:p/o:p [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Roger Bivand Economic Geography Section, Department of Economics, Norwegian School of Economics and Business Administration, Helleveien 30, N-5045 Bergen, Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43 e-mail: [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Help with mixed effects models
Dear R-users I have problems using lme The model i want to fit can be viewed as a two-level bivariate model Two-level bivariate: bivariate (S coded as -1,T coded as 1) endpoint within trial OR It can equivalently be considered as a three-level model.Three-level: endpoint within patient, patient within trial. My code tries to model the levels through a RANDOM statement and a within group correlation structure. -- Now then, i used the following R code: bm - lme (outcome~ -1 + as.factor(endpoint)+ as.factor(endpoint):trt, data=datt, random=~-1 + as.factor(endpoint) + as.factor(endpoint):trt |as.factor(Trial), corr = corSymm(form~-1+as.factor(endpoint)|trial/subject)) I beleive the fixed effects part of the code is okay. My intention for the random effects part is to estimate an intercept and treatment effect for each endpoint at the trial level. The correlation structure should produce a within correlation matrix for the enpoints at the subject level. Thus the random effects matrix is 4 by 4 and the within correlation matrix is 2 by 2 When i run the code in R, i get the following error message Error in Initialize.corSymm(X[[2]], ...) : Initial value for corSymm parameters of wrong dimension I hope someone will correct my codes . Kind regards Pryseley - Ring in the New Year with Photo Calendars. Add photos, events, holidays, whatever. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Current state of support for BUGS access for Linux users?
Hi all. Sorry for any confusion. At one point, close to a year ago I think, I put in the effort to get the bugs() function in R to work with OpneBUGS. But at the time, I couldn't get OpenBUGS to do much. After some effort I got it to work on some simple examples but I couldn't get it to work on some of my other examples (e.g., a mulitlevel logistic regression). In the meantime, Uwe and Sibylle set up R2WinBUGS. I still think that BRugs would be improved by an R front-end (basically, the bugs() function adapted to BRugs). But I haven't tried to get OpenBUGS working for awhile so I don't know if it currently fits the models I'd like to fit. It certainly wouldn't be difficult to extend the bugs() fuction to do the appropriate BRugs calls. In a different direction, Jouni Kerman wrote a function to take output from chains of iterative simulation--not necessarily from Bugs, they could come from mcmcsamp() or just your own customized Gibbs or Metropolis program--and convert them into bugs objects, which is nice because then you can automatically monitor convergence and look at the inferences using the built-in plot() and print() options in R2WinBUGS. We'll try to put this in R2WinBUGS soon. A key part of a bugs object is the way a vector or matrix of parameters with the same name is structured as a vector or matrix, not simply as part of a long vector of parameters (this can be clearly seen in the right half of the plot() output for a hiearchical model fit in Bugs, and also can be seen if you fit a bugs model and then do attach.bugs() and look at simulations of parameter vectors or arrays. Gregor Gorjanc wrote: Paul Johnson wrote: Do you mean to say that you have actually made OpenBUGS run with R2WinBUGS in Linux? No, I did not say this. Gelman's page seems to state that OpenBUGS support is brought in from BRugs, which is still Windows-only. Well, Gelman changed his site a bit. Few days (weeks?) ago there was a red text at the top of [1] stating that bugs.R was accomodated to work with OpenBUGS. Since bugs.R was base for R2WinBUGS, I conclude that BRugs was not involved here. Since then Gelman was active with his pages and things changed/improved. I do not know what he has done to handle also OpenBUGS. Perhaps Gelman and R2WinBUGS maintainers could tell us and I hope that upstream (by Gelman) changes will find way into R2WinBUGS package. Re R packages: - R2WinBUGS is compatible with WinBUGS-1.4.x only, its newest version can speak with WinBUGS under wine thanks to user contributions. But it still depends on WinBUGS-1.4.x, hence Windows only (considering wine as Windows). However Andrew Gelman, has added also support[1] for OpenBUGS so this might be also a good news for Linux if wine is used. Changelog can be found at [2]. [1]http://www.stat.columbia.edu/~gelman/bugsR/ [2]http://www.stat.columbia.edu/~gelman/bugsR/bugs.R -- Andrew Gelman Professor, Department of Statistics Professor, Department of Political Science [EMAIL PROTECTED] www.stat.columbia.edu/~gelman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] compiling R on powerpc-ibm-aix5.1.0.0
Hi, As the message pointed out, the current linker cannot produce share library. However gcc used option -bdynamic, is that for share library? I guess. Need to read its gcc manual. It might be caused by very slight difference of the gcc compiler in your AIX. Please seek for Zhao Wei's help. He has very good experience in operating UNIX-like system. If the gcc in your AIX is with low version, not suitable, how to do? hmmm ... anyway, let zhaowei have a look firstly. Regards, Lu -Original Message- From: Haiyan Pan [mailto:[EMAIL PROTECTED] Sent: 18 January 2006 08:51 To: r-help Cc: Dan Yang; Lu Qiang Subject: compiling R on powerpc-ibm-aix5.1.0.0 r-help, I am trying to compile R on an powerpc-IBM-AIX5.1.0.0 machine, Is R suitable to be used in this system? The R Installation and Administration document mentioned rs6000-ibm-aix not powerpc-IBM-AIX5.1.0.0 . When I tried to compile R.2.2.0 in powerpc-IBM-AIX5.1.0.0 using the following steps: (1) ./configure There is an error message :configure: error: --with-readline=yes (default) and headers/libs are not available So I used ./configure --with-readline=no --prefix=/export/home/hypan/R/bin instead, and configure is through. The following is information when finished configure R is now configured for powerpc-ibm-aix5.1.0.0 Source directory: . Installation directory:/export/home/hypan/R/bin C compiler:gcc -mno-fp-in-toc -g -O2 C++ compiler: g++ -g -O2 Fortran compiler: g77 -g -O2 Interfaces supported: X11, tcltk External libraries: Additional capabilities: PNG, JPEG, MBCS, NLS Options enabled: R profiling Recommended packages: yes configure: WARNING: you cannot build info or html versions of the R manuals configure: WARNING: I could not determine a browser configure: WARNING: I could not determine a PDF viewer *** Warning: the GNU linker, at least up to release 2.9.1, is reported *** to be unable to reliably create shared libraries on AIX. *** Therefore, libtool is disabling shared libraries support. If you *** really care for shared libraries, you may want to modify your PATH *** so that a non-GNU linker is found, and then restart. (2) make make is failed, the error message is: gcc -Wl,-bdynamic -Wl,-bE:../../etc/R.exp -Wl,-bM:SRE -L/usr/local/lib -o R.bin Rmain.o CConverters.o CommandLineArgs.o Rdynload.o Renviron.o RNG.o apply.o arithmetic.o apse.o array.o attrib.o base.o bind.o builtin.o character.o coerce.o colors.o complex.o connections.o context.o cov.o cum.o dcf.o datetime.o debug.o deparse.o deriv.o dotcode.o dounzip.o dstruct.o duplicate.o engine.o envir.o errors.o eval.o format.o fourier.o gevents.o gram.o gram-ex.o graphics.o identical.o internet.o iosupport.o lapack.o list.o logic.o main.o mapply.o match.o memory.o model.o names.o objects.o optim.o optimize.o options.o par.o paste.o pcre.o platform.o plot.o plot3d.o plotmath.o print.o printarray.o printvector.o printutils.o qsort.o random.o regex.o registration.o relop.o saveload.o scan.o seq.o serialize.o size.o sort.o source.o split.o sprintf.o startup.o subassign.o subscript.o subset.o summary.o sysutils.o unique.o util.o version.o vfonts.o xxxpr.o ../unix/libunix.a ../appl/libap! pl.a ../nmath/libnmath.a -lg2c -lm -lgcc_s /usr/local/lib/gcc-lib/powerpc-ibm-aix5.1.0.0/3.3/libgcc.a -lg /lib/crt0.o ../extra/zlib/libz.a ../extra/bzip2/libbz2.a ../extra/pcre/libpcre.a ../extra/intl/libintl.a -ldl -lm -lc -liconv /usr/local/lib/gcc-lib/powerpc-ibm-aix5.1.0.0/3.3/../../../../powerpc-ibm-aix5.1.0.0/bin/ld: target dynamic not found collect2: ld returned 1 exit status make: 1254-004 The error code from the last command is 1. Stop. make: 1254-004 The error code from the last command is 2. Stop. make: 1254-004 The error code from the last command is 1. Stop. make: 1254-004 The error code from the last command is 1. Stop. your help will be greatly appreciated. Thanks! Haiyan = = = = = = = = = = = = = = = = = = = = Haiyan Pan [EMAIL PROTECTED] Tel: 021-64363311-123 Shanghai Center for Bioinformatics Technology Floor 12th,100# QinZhou Road Shanghai,China,200235 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] y axis text truncated
I have been trying to find which par settings can help me avoid truncated text at the bottom of the y axis in a mosaic plot (created when I plot a result of a 2d xtabs) without much success. Using las=1 has helped but the text (the 500+ level) is still cropped. I get the same result on XP/2.2.0 and FC4/2.2.1. Any tips would be appreciated. # dput(foo.df) foo.df = structure(list(vol1 = structure(c(1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6 ), .Label = c(100, 101-250, 251-500, 501-750, 751-1000, 1000+), class = factor), vol2 = structure(c(1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5), .Label = c(20, 20-50, 50-100, 100-500, 500+), class = factor), Freq = c(4, 3, 0, 0, 2, 0, 4, 3, 6, 4, 1, 2, 1, 3, 3, 4, 5, 2, 3, 1, 3, 2, 2, 12, 0, 0, 1, 0, 2, 4)), .Names = c(vol1, vol2, Freq), row.names = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30), class = data.frame) xtabs(Freq ~ vol1 + vol2, data=foo.df) vol2 vol1 20 20-50 50-100 100-500 500+ 100 4 4 1 30 101-2503 3 3 10 251-5000 6 3 31 501-7500 4 4 20 751-1000 2 1 5 22 1000+ 0 2 2 124 plot(xtabs(Freq ~ vol1 + vol2, data=foo.df)) plot(xtabs(Freq ~ vol1 + vol2, data=foo.df), las=1) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Current state of support for BUGS access for Linux users?
Andrew Gelman wrote: Hi all. Sorry for any confusion. At one point, close to a year ago I think, I put in the effort to get the bugs() function in R to work with OpneBUGS. But at the time, I couldn't get OpenBUGS to do much. After some effort I got it to work on some simple examples but I couldn't get it to work on some of my other examples (e.g., a mulitlevel logistic regression). In the meantime, Uwe and Sibylle set up R2WinBUGS. I still think that BRugs would be improved by an R front-end (basically, the bugs() function adapted to BRugs). But I haven't tried to get Andrew, let me add a few details: The BRugs equivalent to your bugs() function is called BRugsFit() - it runs the whole stuff fo you with one function call. Some plots can be produced afterwards with a second call. For Linux, what we really want is to get the OpenBUGS library (i.e. BRugs.so) run in Linux natively. All other stuff are ugly workarounds. Since it looks not that promising to get a quick solution these days, it might make sense to communicate to some WinBUGS/OpenBUGS running under wine. Then, R2WinBUGS could be used for WinBUGS 1.4.x under wine and for OpenBUGS something similar to - communication structure of R2WinBUGS *and* - interface to the OpenBUGS command language is desirable. In principle, one could look at BRugs and take the relevant information originally passed through .C() calls. Instead of calling .C(), print the commands into text files and start OpenBUGS so that it processed the generated command files. Unfortunately, one cannot run interactive simulations that way and is limited to the functionality R2WinBUGS already provides. OpenBUGS working for awhile so I don't know if it currently fits the models I'd like to fit. It certainly wouldn't be difficult to extend the bugs() fuction to do the appropriate BRugs calls. In a different direction, Jouni Kerman wrote a function to take output from chains of iterative simulation--not necessarily from Bugs, they could come from mcmcsamp() or just your own customized Gibbs or Metropolis program--and convert them into bugs objects, which is nice because then you can automatically monitor convergence and look at the inferences using the built-in plot() and print() options in R2WinBUGS. We'll try to put this in R2WinBUGS soon. In which case one could use the coda package. A key part of a bugs object is the way a vector or matrix of parameters with the same name is structured as a vector or matrix, not simply as part of a long vector of parameters (this can be clearly seen in the right half of the plot() output for a hiearchical model fit in Bugs, and also can be seen if you fit a bugs model and then do attach.bugs() and look at simulations of parameter vectors or arrays. Hope this is the case in all representations of similar objects. Anyway, among build-in plot facilities, BRugs supports conversion to coda's mcmc objects. Uwe Gregor Gorjanc wrote: Paul Johnson wrote: Do you mean to say that you have actually made OpenBUGS run with R2WinBUGS in Linux? No, I did not say this. Gelman's page seems to state that OpenBUGS support is brought in from BRugs, which is still Windows-only. Well, Gelman changed his site a bit. Few days (weeks?) ago there was a red text at the top of [1] stating that bugs.R was accomodated to work with OpenBUGS. Since bugs.R was base for R2WinBUGS, I conclude that BRugs was not involved here. Since then Gelman was active with his pages and things changed/improved. I do not know what he has done to handle also OpenBUGS. Perhaps Gelman and R2WinBUGS maintainers could tell us and I hope that upstream (by Gelman) changes will find way into R2WinBUGS package. Re R packages: - R2WinBUGS is compatible with WinBUGS-1.4.x only, its newest version can speak with WinBUGS under wine thanks to user contributions. But it still depends on WinBUGS-1.4.x, hence Windows only (considering wine as Windows). However Andrew Gelman, has added also support[1] for OpenBUGS so this might be also a good news for Linux if wine is used. Changelog can be found at [2]. [1]http://www.stat.columbia.edu/~gelman/bugsR/ [2]http://www.stat.columbia.edu/~gelman/bugsR/bugs.R __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Current state of support for BUGS access for Linux users?
Uwe, Does coda put parameters into arrays based on their names? I had the impression that coda just considered parameters as one very long vector (as in the summary output inside the Bugs window), hence it would not make those nice graphs on the right side of the output from plot() applied to a Bugs object. But if coda is the standard, then I suppose what's needed is to adapt Jouni's function to also convert coda objects into bugs objects? Then any iterative simulation object could be displayed and accessed by parameters, not simply scalar components of parameters. (Or maybe coda can do more than I realize, I don't know.) Andrew -- Uwe Ligges wrote: In which case one could use the coda package. A key part of a bugs object is the way a vector or matrix of parameters with the same name is structured as a vector or matrix, not simply as part of a long vector of parameters (this can be clearly seen in the right half of the plot() output for a hiearchical model fit in Bugs, and also can be seen if you fit a bugs model and then do attach.bugs() and look at simulations of parameter vectors or arrays. Hope this is the case in all representations of similar objects. Anyway, among build-in plot facilities, BRugs supports conversion to coda's mcmc objects. Uwe -- Andrew Gelman Professor, Department of Statistics Professor, Department of Political Science [EMAIL PROTECTED] www.stat.columbia.edu/~gelman __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] y axis text truncated
Use par(xpd=TRUE). On Wed, 18 Jan 2006, paul sorenson wrote: I have been trying to find which par settings can help me avoid truncated text at the bottom of the y axis in a mosaic plot (created when I plot a result of a 2d xtabs) without much success. Using las=1 has helped but the text (the 500+ level) is still cropped. I get the same result on XP/2.2.0 and FC4/2.2.1. Any tips would be appreciated. # dput(foo.df) foo.df = structure(list(vol1 = structure(c(1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6 ), .Label = c(100, 101-250, 251-500, 501-750, 751-1000, 1000+), class = factor), vol2 = structure(c(1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5), .Label = c(20, 20-50, 50-100, 100-500, 500+), class = factor), Freq = c(4, 3, 0, 0, 2, 0, 4, 3, 6, 4, 1, 2, 1, 3, 3, 4, 5, 2, 3, 1, 3, 2, 2, 12, 0, 0, 1, 0, 2, 4)), .Names = c(vol1, vol2, Freq), row.names = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30), class = data.frame) xtabs(Freq ~ vol1 + vol2, data=foo.df) vol2 vol1 20 20-50 50-100 100-500 500+ 100 4 4 1 30 101-2503 3 3 10 251-5000 6 3 31 501-7500 4 4 20 751-1000 2 1 5 22 1000+ 0 2 2 124 plot(xtabs(Freq ~ vol1 + vol2, data=foo.df)) plot(xtabs(Freq ~ vol1 + vol2, data=foo.df), las=1) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] some EPS rotated in journal preview
I am trying to send a manuscript to a journal. One of the figures build by R is in the right orientation and 4 are rotated clockwise 90 deg in the preview. I used the right click save to PS option and I used the command line postscript(c:/temp/fig04.eps,bg=transparent,onefile = TRUE ,pointsize=20,paper = letter,height=8,width=8,horizontal=FALSE,family = Helvetica, font = Helvetica) I treed Horizontal=TRUE Ghostsript show the rotated image but not the preview from the journal. :-( Is there anything to change that the - unknown system - of the journal will be forced to display the image in the right direction? Regards Knut __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Possible improvement in lm
Folks, I do a series of regressions (one for each quarter in the dataset) and then go and extract the residuals from each stored lm object that is returned as follows: vResiduals - as.vector(unlist(resid(lQuarterlyRegressions[[i]]))); Here lQuarterlyRegressions is a vector of objects returned by lm(). Next, I may go find outliers using identify() on a plot or do some other analysis which tells me which row of the quarterly data I need to take a closer look at. However, if I try to match some point in one of the quarters that I have with its residual, then I have to drop the points from my current Data which have NA's for either the explanatory variables or the explained, so that the vector or residuals and the data have the same indexes. This lead to some serious confusion/bugs for me, and I am wondering if it might not be better for lm to put an NA into those rows where the point was dropped because of NA's in the explanatory or explained variables (currently it just returns nothing at that index). Ofcourse, there might be some arguments against this idea, and I would be interested to hear them. Thank you for your time and attention, -- Vivek Satsangi Student, Rochester, NY USA __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Data frame index?
On 1/18/2006 2:35 AM, Kenneth Cabrera wrote: Hi, R users: I have a data.frame (not a matrix), I got a vector with the same length as the number of records (rows) of the data frame, and each element of that vector is the column number (in a specific range of columns) of the corresponding record that I must set to zero. How can I do this without a for loop? It sounds as though you've found that you can use two-column matrix indexing on a data frame for reading but not assigning. You create a matrix where the first column is the row number, and the second column is the column number. Then indexing by that selects those particular elements in order. For instance, if you have named your vector of columns cols, you'd do my.data.frame[ cbind(1:rows, cols) ] - 0 Here's an example: df x y 1 1 a 2 1 a 3 1 a 4 1 a 5 1 a 6 1 a 7 1 a 8 1 a 9 1 a 10 1 a df[cbind(1:4,c(1,2,1,2))] [1] 1 a 1 a But df[cbind(1:4,c(1,2,1,2))] - 0 Error in [-.data.frame(`*tmp*`, cbind(1:4, c(1, 2, 1, 2)), value = 0) : only logical matrix subscripts are allowed in replacement To get around this, construct the logical matrix using this method, then use it as an index: mat - matrix(FALSE, 10, 2) mat[cbind(1:4,c(1,2,1,2))] - TRUE df[mat] - 0 Warning message: invalid factor level, NAs generated in: [-.factor(`*tmp*`, thisvar, value = 0) df xy 1 0a 2 1 NA 3 0a 4 1 NA 5 1a 6 1a 7 1a 8 1a 9 1a 10 1a If your columns are all numeric, you won't get the warning I got. Duncan Murdoch __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] R Wiki and R-sig-wikii
Hello all, This is to announce the creation of R-sig-wiki, a new R SIG (Special Interest Group) mailing list dedicated to the elaboration and maintenance of a R Wiki. You can subscribe at: https://stat.ethz.ch/mailman/listinfo/r-sig-wiki. There is currently a prototype for a new R Wiki at http://www.sciviews.org/_rgui/wiki (temporary address). The main idea is to offer a site where users could collaborate in writting various kind of documentation for R. It mainly targets R beginners, but is left open for more advanced sections. Any R user interested in the setting up of this Wiki is warmly welcome to participate and to subscribe to R-sig-wiki. For the others, we will send an announcement for the final R Wiki on this list when it will be ready. Hereunder is a (rather long) summary of the discussion we had so far on this topic. Best, Philippe Grosjean == There is a proposition to create a R Wiki, stimulated essentially by two facts: 1) The traffic in R-Help is very high, with many trivial questions asked repeatedly. Obviously, searching in R-Help archives is not so obvious for some R users. Perhaps another presentation, like plain HTML pages would be fine. Since the building of these HTML pages should be a collaborative work, a Wiki seems to be a possible solution (recall that a Wiki is essentially a simple way to collaborate on writting Web pages; see Wikipedia definition at http://en.wikipedia.org/wiki/Wiki). 2) Some threads in R-Help are very long and difficult to follow in the form of a succession of emails. Again, a more structured presentation, allowed by HTML / Wiki pages is suggested as a possible solution. There is already one attempt to build a Wiki by Detlef Steuer at: http://fawn.unibw-hamburg.de/cgi-bin/Rwiki.pl. Despite not much effort was put in this wiki, several people collaborated to it with time, but it appears that it was below the minimum required to make it fly as it should. Being unaware of that Wiki, I proposed recently (beginning of January 2006) a prototype of another R Wiki, just to explore if and how it could answer those two problems on R-Help. The prototype is at http://www.sciviews.org/_rgui/wiki (temporary address). There is a couple of other Wikis dedicated to R floating around, like: http://www.okada.jp.org/RWiki/ (in Japanese) and http://learnserver.csd.univie.ac.at/rcomwiki dedicated to R(D)COM and RExcel essentially. After a discussion, the R-Core Team decided to give support to one or several intiatives to make a R Wiki, with a big concern about the quality of information in the Wiki and how to keep it in phase with the rapid development of R. Here is the mail send by Martin Maechler: Martin Maechler wrote: We've had a small review time within R-core on this topic, amd would like to state the following: The R-core team welcomes proposals to develop an R-wiki. - We would consider linking a very small number of Wikis (ideally one) from www.r-project.org and offering an address in the r-project.org domain (such as 'wiki.r-project.org'). - The core team has no support time to offer, and would be looking for a medium-term commitment from a maintainer team for the Wiki(s). - Suggestions for the R documentation would best be filtered through the Wiki maintainers, who could e.g. supply suggested patches during the alpha phase of an R release. -- Our main concerns have been about ensuring the quality of such extra documentation projects, hence the 2nd point above. Several of our more general, not mainly R, experiences have been of outdated web pages which are continued to be used as reference when their advice has long been superseded. I think it's very important to try ensuring that this won't happen with an R Wiki. === After that announcement, several people started to discuss the structure and content of the Wiki (Tony Plate and Ben Bolker took the initiative to propose one structure on the experimental R Wiki at: http://www.sciviews.org/_rgui/wiki/doku.php?id=varia:organization_discussion, and I propose to leave it open until the end of February. At that time, we will implement the proposed structure as a proof-of-concept. Another topic is about how to allow discussion on the R documentation in the Wiki. This topic was fed by Frank Harrell Jr, Gabor Grothendieck and others and it leads to a trial on: http://www.sciviews.org/_rgui/wiki/doku.php?id=varia:test_include. Basically, all CRAN and Bioconductor .Rd files would be converted automatically to read-only Wiki pages (regularly updated) that are themselves included in fully
[R] Within-Subjects ANOVA comparisons of individual means
I am having problems with comparing individual means in a within-subjects ANOVA. From my understanding, TukeyHSD is not appropriate in this context. So I am trying to compute contrasts, as follows: seven subjects participated in each of 6 conditions (intervals). subject = factor(rep(c(1:7), each = 6)) interval = factor(rep(c(1:6), 7)) and here is the dependent variable: dv = c(3.3868, 3.1068, 1.7261, 1.5415, 1.7356, 0.7538, + 2.5957, 1.5666, 1.1984, 1.2761, 1.0022, 0.8597, + 3.9819, 3.1506, 1.5824, 1.7400, 1.4248, 0.6519, + 2.2521, 1.5248, 1.1209, 1.2193, 1.1994, 2.0910, + 2.4661, 1.3863, 1.3591, 0.9163, 1.3976, 1.7471, + 3.2486, 1.9492, 2.4228, 1.1276, 1.2836, 0.9814, + 1.7148, 1.7278, 2.7433, 1.4924, 1.0992, 0.7821) d = data.frame(subject, interval, dv) next I'm defining a contrast matrix: con = matrix(c(1, -1, 0, 0, 0, 0, 0, 1, -1, 0, 0, 0, 0, 0, 1, -1, 0, 0, 0, 0, 0, 1, -1, 0, 0, 0, 0, 0, 1, -1), nrow=6, ncol=5, byrow=F) con [,1] [,2] [,3] [,4] [,5] [1,]10000 [2,] -11000 [3,]0 -1100 [4,]00 -110 [5,]000 -11 [6,]0000 -1 contrasts(d$interval)=con and then I'm doing the ANOVA aovRes = aov(dv~interval+Error(subject/interval), data=d) summary(aovRes) Error: subject Df Sum Sq Mean Sq F value Pr(F) Residuals 6 2.48531 0.41422 Error: subject:interval Df Sum Sq Mean Sq F valuePr(F) interval 5 13.8174 2.7635 8.7178 3.417e-05 *** Residuals 30 9.5098 0.3170 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 but if I want to look at the contrasts, something has gone wrong: summary.aov(aovRes, split=list(interval = list(i1 vs i2 = 1, i2 vs i3 = 2, i3 vs i4 = 3, i4 vs i5 = 4, i5 vs i6 = 5))) Error in 1:object$rank : NA/NaN argument aovRes$contrasts NULL Can anybody help? Thank you very much, -Steffen __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Loading of namespace on load of .Rdata (was strange behaviour of load)
Last week Giovanni Parrinello posted a message asking why various packages were loaded when he loaded an .Rdata file. Brian Ripley replied saying he thought it was because the saved workspace contained a reference to the namespace of ipred. (Correspondence copied below). This begs the question: how did the reference to the namespace of ipred come to be in the .Rdata file? Brian did say it is likely to be because the workspace contained object(s) saved with environment the namespace of ipred - but how would this come about? In this case I think is because the .Rdata file contained an object whose *parent* environment was the namespace of ipred. Take the following example from ?bagging (having loaded ipred): data(BreastCancer) mod - bagging(Class ~ Cl.thickness + Cell.size + + Cell.shape + Marg.adhesion + + Epith.c.size + Bare.nuclei + + Bl.cromatin + Normal.nucleoli + + Mitoses, data=BreastCancer, coob=TRUE) environment(mod$mtrees[[1]]$btree$terms) environment: 024E8138 parent.env(environment(mod$mtrees[[1]]$btree$terms)) environment: namespace:ipred This occurs because the terms object is taken from the model frame which was evaluated within the environment of a function from the ipred package (here ipred:::irpart). Therefore I think the behaviour observed by Giovanni will only occur in unusual circumstances: when the workspace contains a formula object, a terms object, a function, or some other object with a non-NULL environment, which has been created in the environment of a packaged function. In particular, this would not always occur with a packaged model fitting function, e.g. (from ?loglm in MASS) library(MASS) minn38a - array(0, c(3,4,7,2), lapply(minn38[, -5], levels)) minn38a[data.matrix(minn38[,-5])] - minn38$f fm - loglm(~1 + 2 + 3 + 4, minn38a) environment(fm$terms) environment: R_GlobalEnv in this case because the terms component is obtained from the formula, whose environment is .GlobalEnv. So, I have two points on this (more for R-devel than R-help now) 1. There is a more general situation where it would be useful to load the namespace of a package after loading a saved workspace: when the workspace contains objects of a class for which special methods are required. E.g. if 'fm' from the example above were saved in a workspace, the namespace of MASS would not be loaded when the workspace was loaded into R. Thus unless MASS was loaded by the user, default methods would be used by summary(), print() etc rather than the specialised methods for objects of class loglm. Of course the user should quickly realise this, but there may be cases where the default method gives a convincing but incorrect or unhelpful result. An alternative would be to add an attribute to objects of class loglm (say), e.g. attr(loglmObject, .Environment) - environment(MASS) so that the namespace would automatically be loaded when it is required. [In fact alternatives such as environment(loglmObject) - environment(MASS) or loglmObject$anyoldthing - environment(MASS) would work just as well, but perhaps the first suggestion is neatest.]. What do others think of this idea? Should it (or an equivalent idea) be encouraged amongst package writers? 2. In the case highlighted by Giovanni, the namespace of ipred was loaded, but the package was not. This would be fine, except that the packages on which ipred depends *were* loaded. This seems inconsistent. I guess as long as there are packages without namespaces though, this is the only way to proceed. Perhaps in the meantime, package authors should be encouraged to use importFrom() rather than import()? Or perhaps where packages do have namespaces, only the namespace should be loaded in such a case. Heather From: Prof Brian Ripley [EMAIL PROTECTED] Date: 12 January 2006 08:21:35 GMT To: giovanni parrinello [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Subject: Re: [R] Strange behaviour of load On Wed, 11 Jan 2006, giovanni parrinello wrote: Dear All, simetimes when I load an Rdata I get this message ### Code: load('bladder1.RData') Carico il pacchetto richiesto: rpart ( Bad traslastion: Load required package-...) Carico il pacchetto richiesto: MASS Carico il pacchetto richiesto: mlbench Carico il pacchetto richiesto: survival Carico il pacchetto richiesto: splines Carico il pacchetto richiesto: 'survival' The following object(s) are masked from package:Hmisc : untangle.specials Carico il pacchetto richiesto: class Carico il pacchetto richiesto: nnet # So I have many unrequired packages loaded. Any idea? They are required! My guess is that you have object(s) saved with environment the namespace of some package, and loading that namespace is pulling these in. The only CRAN package which requires mlbench appears to be ipred, and that requires all of those except splines,
[R] Coercing a list to integer?
Dear group, I am nearly beside myself. After an entire night spent on a niggling little detail, I am no closer to to the truth. I loaded an Excel file in .csv form into R. It apparentely loads as a list, but not the kind of list you can use. Oh no, it converts into a list that cannot be converted into an integer, numeric, or vector, only a matrix, whihc is useless without integers. How can I get a list of the form [1] 1,2,3,4,5 into the form [1] 1 [2] 2 [3] 3 [4] 4 [5] 5? Depending on hwo you define a list, apparentely, it goes one way or the other. x - list(1:5) means you have [1] 1,2,3,4,5 y - list(1,2,3,4,5) means you have [1] 1 [2] 2 [3] 3 [4] 4 [5] 5 Can anyone help?# I woudl greatly appreciate it. Sincerely, Norman Goodacre - [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] se.fit in predict.nls
The option se.fit in predict.nls is currently ignored. Is there any other function available to calculate the error in the predictions? Thanks, Manuel __ LLama Gratis a cualquier PC del Mundo. Llamadas a fijos y móviles desde 1 céntimo por minuto. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Coercing a list to integer?
?unlist y - list(1,2,3,4,5) y [[1]] [1] 1 [[2]] [1] 2 [[3]] [1] 3 [[4]] [1] 4 [[5]] [1] 5 unlist(y) [1] 1 2 3 4 5 On 1/18/06, Norman Goodacre [EMAIL PROTECTED] wrote: Dear group, I am nearly beside myself. After an entire night spent on a niggling little detail, I am no closer to to the truth. I loaded an Excel file in .csv form into R. It apparentely loads as a list, but not the kind of list you can use. Oh no, it converts into a list that cannot be converted into an integer, numeric, or vector, only a matrix, whihc is useless without integers. How can I get a list of the form [1] 1,2,3,4,5 into the form [1] 1 [2] 2 [3] 3 [4] 4 [5] 5? Depending on hwo you define a list, apparentely, it goes one way or the other. x - list(1:5) means you have [1] 1,2,3,4,5 y - list(1,2,3,4,5) means you have [1] 1 [2] 2 [3] 3 [4] 4 [5] 5 Can anyone help?# I woudl greatly appreciate it. Sincerely, Norman Goodacre - [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Jim Holtman Cincinnati, OH +1 513 247 0281 What the problem you are trying to solve? [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Coercing a list to integer?
i am not sure i clearly understood... do you want to coerce a vector into a list of its elements ?! like this: split(1:5,1:5) Norman Goodacre a écrit : Dear group, I am nearly beside myself. After an entire night spent on a niggling little detail, I am no closer to to the truth. I loaded an Excel file in .csv form into R. It apparentely loads as a list, but not the kind of list you can use. Oh no, it converts into a list that cannot be converted into an integer, numeric, or vector, only a matrix, whihc is useless without integers. How can I get a list of the form [1] 1,2,3,4,5 into the form [1] 1 [2] 2 [3] 3 [4] 4 [5] 5? Depending on hwo you define a list, apparentely, it goes one way or the other. x - list(1:5) means you have [1] 1,2,3,4,5 y - list(1,2,3,4,5) means you have [1] 1 [2] 2 [3] 3 [4] 4 [5] 5 Can anyone help?# I woudl greatly appreciate it. Sincerely, Norman Goodacre - [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Possible improvement in lm
1. Try using lm(...whatever..., na.action = na.exclude) 2. Be sure to read the note on Using Time Series in ?lm 3. The dyn package will accept ts, irts, its and zoo class time series and output time series for the residuals. Just preface lm with dyn$. e.g. library(dyn) # test data set.seed(1) x - ts(1:10, start = 2000, freq = 4) x[5] - NA y - x + rnorm(10) # regress series y against series x y.lm - dyn$lm(y ~ x) resid(y.lm) # note that residuals are a time series On 1/18/06, Vivek Satsangi [EMAIL PROTECTED] wrote: Folks, I do a series of regressions (one for each quarter in the dataset) and then go and extract the residuals from each stored lm object that is returned as follows: vResiduals - as.vector(unlist(resid(lQuarterlyRegressions[[i]]))); Here lQuarterlyRegressions is a vector of objects returned by lm(). Next, I may go find outliers using identify() on a plot or do some other analysis which tells me which row of the quarterly data I need to take a closer look at. However, if I try to match some point in one of the quarters that I have with its residual, then I have to drop the points from my current Data which have NA's for either the explanatory variables or the explained, so that the vector or residuals and the data have the same indexes. This lead to some serious confusion/bugs for me, and I am wondering if it might not be better for lm to put an NA into those rows where the point was dropped because of NA's in the explanatory or explained variables (currently it just returns nothing at that index). Ofcourse, there might be some arguments against this idea, and I would be interested to hear them. Thank you for your time and attention, -- Vivek Satsangi Student, Rochester, NY USA __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] some EPS rotated in journal preview
The problem is a well-known one in viewers looking at whole pages, especially PS - PDF converters. R figures are particularly vulnerable as they have text running both horizontally and vertically (with normal axes). Please do follow exactly the advice on the postscript help page. The postscript produced for a single R plot is EPS (_Encapsulated PostScript_) compatible, and can be included into other documents, e.g., into LaTeX, using '\includegraphics{filename}'. For use in this way you will probably want to set 'horizontal = FALSE, onefile = FALSE, paper = special'. If you have done that, suggest to your publisher that they turn autorotation off. On Wed, 18 Jan 2006, Knut Krueger wrote: I am trying to send a manuscript to a journal. One of the figures build by R is in the right orientation and 4 are rotated clockwise 90 deg in the preview. I used the right click save to PS option and I used the command line postscript(c:/temp/fig04.eps,bg=transparent,onefile = TRUE ,pointsize=20,paper = letter,height=8,width=8,horizontal=FALSE,family = Helvetica, font = Helvetica) I treed Horizontal=TRUE Ghostsript show the rotated image but not the preview from the journal. :-( Is there anything to change that the - unknown system - of the journal will be forced to display the image in the right direction? -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] some EPS rotated in journal preview
Knut Krueger wrote: I am trying to send a manuscript to a journal. One of the figures build by R is in the right orientation and 4 are rotated clockwise 90 deg in the preview. I used the right click save to PS option and I used the command line postscript(c:/temp/fig04.eps,bg=transparent,onefile = TRUE ,pointsize=20,paper = letter,height=8,width=8,horizontal=FALSE,family = Helvetica, font = Helvetica) Please use paper=special for files to be included in documents. I treed Horizontal=TRUE Ghostsript show the rotated image but not the preview from the journal. :-( Is there anything to change that the - unknown system - of the journal will be forced to display the image in the right direction? Always hard to tell what an unknown system is doing. We do not know of any rotation problems with R graphics, hence please ask the poeple running that unknown system. Uwe Ligges Regards Knut __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Possible improvement in lm
It seems you are looking for na.action = na.exclude. This is described in all good books on S/R (e.g. MASS4, p. 141) and also in an answer on this list already this week. It is even on the help pages for residuals.lm and fitted. On Wed, 18 Jan 2006, Vivek Satsangi wrote: Folks, I do a series of regressions (one for each quarter in the dataset) and then go and extract the residuals from each stored lm object that is returned as follows: vResiduals - as.vector(unlist(resid(lQuarterlyRegressions[[i]]))); Here lQuarterlyRegressions is a vector of objects returned by lm(). Next, I may go find outliers using identify() on a plot or do some other analysis which tells me which row of the quarterly data I need to take a closer look at. However, if I try to match some point in one of the quarters that I have with its residual, then I have to drop the points from my current Data which have NA's for either the explanatory variables or the explained, so that the vector or residuals and the data have the same indexes. This lead to some serious confusion/bugs for me, and I am wondering if it might not be better for lm to put an NA into those rows where the point was dropped because of NA's in the explanatory or explained variables (currently it just returns nothing at that index). Ofcourse, there might be some arguments against this idea, and I would be interested to hear them. Thank you for your time and attention, -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] princomp() with missing values in panel data?
thank you. I am still not sure how to get the scores in princomp, though: ds= as.data.frame( cbind(rnorm(10),rnorm(10)) ) names(ds)=c(x1,x2) ds[5,]=c(NA,NA) pc= princomp( formula = ~ ds$x1 + ds$x2, na.action=na.omit) ds$pc1 = pc$scores[,1] #-- error, scores has 9 obs, ds has 10 obs is there an elegant method to do this, or do I need to learn how to operate with pc$loadings? (may I also humbly suggest that the default behavior or $scores should be to contain NA in row 5?) Incidentally, R is a lot cleverer than I understand. pc$loadings by itself gives me wonderfully intuitive output, with names, text, different components---but I can still use p$loadings[,2]. I presume that the array operator on the loadings object is overloaded. very nice. [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Documentation in R 2.2.1 for Windows of: rcmd build --docs=
I recommend that the allowable values for docs be included in the short help output generated by the command 'rcmd build --help'. The html help file for build does not contain this information (which is fine since the HTML help system supports multiple operating systems). I did find the options I needed in the html help for install. Given the voluminous (... and inherently somewhat disconnected...) documentation of R, I'm sure someone can point me to other locations where these options are documented. However, since the options seem to be OS specific and there is already help output available from the OS specific command, I recommend that the options be included in the short help listing generated in association with the command. Thanks for all of your hard work. Chuck __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] negative predicted values in poisson glm
Dear R helpers, running the following code of a glm model of the family poisson, gives predicted values 0. Why? library(MASS) library(stats) library(mvtnorm) library(pscl) data(bioChemists) poisson_glm - glm(art ~ fem + mar + kid5 + phd + ment, data = bioChemists, family = poisson) predicted.values = predict(poisson_glm) range(predicted.values) Thank you in advance for any hints. Best regards, P. Olsson [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] some EPS rotated in journal preview
On Wed, 2006-01-18 at 13:04 +0100, Knut Krueger wrote: I am trying to send a manuscript to a journal. One of the figures build by R is in the right orientation and 4 are rotated clockwise 90 deg in the preview. I used the right click save to PS option and I used the command line postscript(c:/temp/fig04.eps,bg=transparent,onefile = TRUE ,pointsize=20,paper = letter,height=8,width=8,horizontal=FALSE,family = Helvetica, font = Helvetica) I treed Horizontal=TRUE Ghostsript show the rotated image but not the preview from the journal. :-( Is there anything to change that the - unknown system - of the journal will be forced to display the image in the right direction? Regards Knut One of the first things to do is to use 'onefile = FALSE', 'horizontal = FALSE' and paper = special'. This is in the Details section of ?postscript, which provides guidance on the creation of EPS files for inclusion in publications. I would try that to see if that provides a more consistent formatting of the plots. You don't indicate what you are using to create the manuscript itself (ie. Word, LaTeX, or ?) to help us in considering other possibilities (such as autorotation). If the above does not help, please provide a reproducible example of the plot code and what you are using for the manuscript. HTH, Marc Schwartz __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Data frame index?
It's worth noting that there are quite a few for loops inside the code used by matrix indexing of data frames. I think a single for-loop over the columns is as good as any, something like DF - data.frame(x=1, y=rep(a, 4), z = 3) ind - c(1,3,3,1) # only numeric cols for(i in unique(ind)) DF[ind==i, i] - 0 DF x y z 1 0 a 3 2 1 a 0 3 1 a 0 4 0 a 3 On Wed, 18 Jan 2006, Duncan Murdoch wrote: On 1/18/2006 2:35 AM, Kenneth Cabrera wrote: Hi, R users: I have a data.frame (not a matrix), I got a vector with the same length as the number of records (rows) of the data frame, and each element of that vector is the column number (in a specific range of columns) of the corresponding record that I must set to zero. How can I do this without a for loop? It sounds as though you've found that you can use two-column matrix indexing on a data frame for reading but not assigning. You create a matrix where the first column is the row number, and the second column is the column number. Then indexing by that selects those particular elements in order. For instance, if you have named your vector of columns cols, you'd do my.data.frame[ cbind(1:rows, cols) ] - 0 Here's an example: df x y 1 1 a 2 1 a 3 1 a 4 1 a 5 1 a 6 1 a 7 1 a 8 1 a 9 1 a 10 1 a df[cbind(1:4,c(1,2,1,2))] [1] 1 a 1 a But df[cbind(1:4,c(1,2,1,2))] - 0 Error in [-.data.frame(`*tmp*`, cbind(1:4, c(1, 2, 1, 2)), value = 0) : only logical matrix subscripts are allowed in replacement To get around this, construct the logical matrix using this method, then use it as an index: mat - matrix(FALSE, 10, 2) mat[cbind(1:4,c(1,2,1,2))] - TRUE df[mat] - 0 Warning message: invalid factor level, NAs generated in: [-.factor(`*tmp*`, thisvar, value = 0) df xy 1 0a 2 1 NA 3 0a 4 1 NA 5 1a 6 1a 7 1a 8 1a 9 1a 10 1a If your columns are all numeric, you won't get the warning I got. Duncan Murdoch __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Coercing a list to integer?
On Wed, 2006-01-18 at 12:27 +, Norman Goodacre wrote: Dear group, I am nearly beside myself. After an entire night spent on a niggling little detail, I am no closer to to the truth. I loaded an Excel file in .csv form into R. It apparentely loads as a list, but not the kind of list you can use. Oh no, it converts into a list that cannot be converted into an integer, numeric, or vector, only a matrix, whihc is useless without integers. How can I get a list of the form [1] 1,2,3,4,5 into the form [1] 1 [2] 2 [3] 3 [4] 4 [5] 5? Depending on hwo you define a list, apparentely, it goes one way or the other. x - list(1:5) means you have [1] 1,2,3,4,5 y - list(1,2,3,4,5) means you have [1] 1 [2] 2 [3] 3 [4] 4 [5] 5 Can anyone help?# I woudl greatly appreciate it. Sincerely, Norman Goodacre Presuming that you used read.csv() or similar, the imported CSV object should be a simple data frame. It is not truly clear here what your problem is relative to how you want to use the imported data. The difference between: list(1:5) [[1]] [1] 1 2 3 4 5 and list(1,2,3,4,5) [[1]] [1] 1 [[2]] [1] 2 [[3]] [1] 3 [[4]] [1] 4 [[5]] [1] 5 is that in the first case, you have a list with one element, which is a vector containing 5 elements: str(list(1:5)) List of 1 $ : int [1:5] 1 2 3 4 5 whereas in the second case, you have a list with five elements, each of which is a vector with one element: str(list(1,2,3,4,5)) List of 5 $ : num 1 $ : num 2 $ : num 3 $ : num 4 $ : num 5 Note also the not so subtle difference where in the first case, the result of the sequence 1:5 yields integers and in the second case, the elements are doubles (numeric), which is the default data type in R. Please provide additional details on what it is you are trying to do here and we can attempt to offer more specific guidance. HTH, Marc Schwartz __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Coercing a list to integer?
of course unlist()/as.list() do that !!! (sorry, i am tired) Norman Goodacre a écrit : Dear group, I am nearly beside myself. After an entire night spent on a niggling little detail, I am no closer to to the truth. I loaded an Excel file in .csv form into R. It apparentely loads as a list, but not the kind of list you can use. Oh no, it converts into a list that cannot be converted into an integer, numeric, or vector, only a matrix, whihc is useless without integers. How can I get a list of the form [1] 1,2,3,4,5 into the form [1] 1 [2] 2 [3] 3 [4] 4 [5] 5? Depending on hwo you define a list, apparentely, it goes one way or the other. x - list(1:5) means you have [1] 1,2,3,4,5 y - list(1,2,3,4,5) means you have [1] 1 [2] 2 [3] 3 [4] 4 [5] 5 Can anyone help?# I woudl greatly appreciate it. Sincerely, Norman Goodacre - [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Windows package upates
Dear list Having just started to use the Windows version, I am very impressed with it's package handling as well as the gui. So I tried to see what was due for update and packages such as Hmisc, Matrix and others came up. But when I had updated them - which took a few goes as something hung between here and Bristol - I noticed that the default packages such as nmle, MASS had disappeared. I re-installed them but is this a glitch or a feature? Best wishes John John Logsdon Try to make things as simple Quantex Research Ltd, Manchester UK as possible but not simpler [EMAIL PROTECTED] [EMAIL PROTECTED] +44(0)161 445 4951/G:+44(0)7717758675 www.quantex-research.com __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] some EPS rotated in journal preview
Uwe Ligges schrieb: Always hard to tell what an unknown system is doing. We do not know of any rotation problems with R graphics, hence please ask the poeple running that unknown system. I tried to ask, they converted the files to rotated tiff (as a free service) but used old files from the first submission, and they did not give me any answer why one file is not rotated the other 3 are rotated ... Sorry that I could not give you any more hints, but I will ask them to switch off auto rotation. Regards Knut __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Own Color Palette
On Wed, 18 Jan 2006, Robert Michael Rausch wrote: I would like to generate a contour-plot according to a master plot. The problem is that the rainbow-palette included in R does not answer this purpose. I need a darker blue, no turquoise, relatively less green, more yellow and more red. Haw can I adjust the rainbow? Alternatively: How can I generate my own palette with at least 100 colors with smooth transitions?o:p/o:p colorRampPalette will interpolate within a list of colors. You might look at the RColorBrewer package to find a list of colors to interpolate within. -thomas __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] negative predicted values in poisson glm
On Wed, 18 Jan 2006, P. Olsson wrote: Dear R helpers, running the following code of a glm model of the family poisson, gives predicted values 0. Why? Look at the help page for predict.glm, particularly the type argument. Then exponentiate the results. -thomas library(MASS) library(stats) library(mvtnorm) library(pscl) data(bioChemists) poisson_glm - glm(art ~ fem + mar + kid5 + phd + ment, data = bioChemists, family = poisson) predicted.values = predict(poisson_glm) range(predicted.values) Thank you in advance for any hints. Best regards, P. Olsson [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Thomas Lumley Assoc. Professor, Biostatistics [EMAIL PROTECTED] University of Washington, Seattle __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] some EPS rotated in journal preview
Prof Brian Ripley schrieb: The problem is a well-known one in viewers looking at whole pages, especially PS - PDF converters. R figures are particularly vulnerable as they have text running both horizontally and vertically (with normal axes). Please do follow exactly the advice on the postscript help page. The postscript produced for a single R plot is EPS (_Encapsulated PostScript_) compatible, and can be included into other documents, e.g., into LaTeX, using '\includegraphics{filename}'. For use in this way you will probably want to set 'horizontal = FALSE, onefile = FALSE, paper = special'. postscript(c:/temp/fig04.eps,bg=transparent,onefile = TRUE ,pointsize=20,paper = letter,height=8,width=8,horizontal=FALSE,family = Helvetica, font = Helvetica) postscript(c:/temp/fig04.eps,bg=transparent,height=8,width=8,bg=transparent,pointsize=20,horizontal = FALSE,onefile = FALSE, paper = special,family = Helvetica, font = Helvetica) There is a error from boxplot without ,height=8,width=8 error in plot.new() : Grafikränder zu groß ( margins to big) but I am afraid that they have autorotation on. I will aks the journal to switch it off. Regards Knut __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] negative predicted values in poisson glm
P. Olsson [EMAIL PROTECTED] writes: Dear R helpers, running the following code of a glm model of the family poisson, gives predicted values 0. Why? library(MASS) library(stats) library(mvtnorm) library(pscl) data(bioChemists) poisson_glm - glm(art ~ fem + mar + kid5 + phd + ment, data = bioChemists, family = poisson) predicted.values = predict(poisson_glm) range(predicted.values) Thank you in advance for any hints. The prediction is on the link scale. -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] some EPS rotated in journal preview
Marc Schwartz schrieb: On Wed, 2006-01-18 at 13:04 +0100, Knut Krueger wrote: One of the first things to do is to use 'onefile = FALSE', 'horizontal = FALSE' and paper = special'. I am afraid the problem is on the journals side, because the wrong postscript line (with letter )is working and I changed the postscript options between both examples no change in the behaviour of the pdf creator of the journal. This is in the Details section of ?postscript, which provides guidance on the creation of EPS files for inclusion in publications. I would try that to see if that provides a more consistent formatting of the plots. You don't indicate what you are using to create the manuscript itself (ie. Word, LaTeX, or ?) to help us in considering other possibilities (such as auto rotation). The journal collects the figures once by once after the manuscript in the original file format. to solve the problem I could change it to tiff and submit it, but the tiff files are look ing not as good as eps files. If the above does not help, please provide a reproducible example of the plot code and what you are using for the manuscript. for your information: the same problem occurs if I plot the data to the graphic device and use the right button - create postscript file One is not rotated the others are rotated. Hope the code helps without data: colored.barlabels -function(x.barplot,x.barcolors,x.xrowlist,cexaxis=0.7,x.thick=F) { countmax-length(x.xrowlist) if (countmax 0) { x.labels-c(1:countmax) for (count in 1:countmax) { for(count2 in 1:countmax)x.labels[count2]- x.labels[count] - x.xrowlist[count] axis(1, at=x.barplot, tick=x.thick, labels=x.labels, col.axis=x.barcolors[count], cex.axis=cexaxis) } }#end if (countmax 0) else Message from colored.barlabel: Nothing to do } postscript(c:/temp//fig04.ps,bg=transparent,onefile = TRUE ,pointsize=20,paper = letter,height=8,width=8,horizontal=FALSE,family = Helvetica, font = Helvetica) barcolors - c(boxcolor1_3,boxcolor1_3,boxcolor1_3,boxcolor4,boxcolor5,boxcolor6,boxcolor7,boxcolor8,boxcolor9,boxcolor10) xrow -c(mean1,mean2,mean8,mean71,mean72,mean78,mean79,mean85,mean86,mean92) sdrow-c(sd1,sd2,sd8,sd71,sd72,sd78,sd79,sd85,sd86,sd92) ci.l-xrow-sdrow/2 ci.h-xrow+sdrow/2 xrowlist - c(Day 1,Day 2,Day 3,Day 4,Day 5,Day 6,Day 7,Day 8,Day 9,Day 10) t.barplot-barplot2(xrow,col = barcolors,cex.lab=1.2,plot.ci=TRUE,ci.l=ci.l,ci.u=ci.h,ylab=mean time till following in sec.) barcolors - c(black,black,black,black,black,black,black,black,black,black,) colored.barlabels(t.barplot,barcolors,xrowlist,0.7,F) dev.off() this is the not rotated figure ant the follwoing is rotated outl.text -function(boxdata,input = FALSE) { countmax-length(boxdata$out) if (countmax 0) { for (count in 1:countmax) { outl.name -boxdata$out[count] if (input==TRUE) outl.name - readline(paste(text for value boxplot Nr:,toString(boxdata$group[count]), - value: ,toString(boxdata$out[count]),: )) text(boxdata$group[count], boxdata$out[count],pos=4,cex=0.8,outl.name) } }#end if (countmax 0) else Message from outl.text: no outlier } postscript(c:/r/anschluss/plots/fig02.ps,height=8,width=8,bg=transparent,pointsize=20,horizontal = FALSE,onefile = FALSE, paper = special,family = Helvetica, font = Helvetica) day1-Day_1 day2-Day_2 day3-Day_8 boxplot.names -c(day 1,day 2,day 3) # Vektor für namen erzeugen boxdata - boxplot(day1,day2,day3, col = boxcolor,names=boxplot.names,ylab=time till following in sec.,cex.lab=size.x.lab,xlab=,cex.lab=size.y.lab ) outl.text(boxdata,input=FALSE) dev.off() Regards Knut __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Windows package upates
On 1/18/2006 9:51 AM, John Logsdon wrote: Dear list Having just started to use the Windows version, I am very impressed with it's package handling as well as the gui. So I tried to see what was due for update and packages such as Hmisc, Matrix and others came up. But when I had updated them - which took a few goes as something hung between here and Bristol - I noticed that the default packages such as nmle, MASS had disappeared. I re-installed them but is this a glitch or a feature? That's not supposed to happen; my guess is that it's related to those hangs you mention. Do you think it tried to update those, but failed because of problems with the mirror? Duncan Murdoch __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Possible improvement in lm
I'm afraid you're a day late and a dollar short: see ?na.exclude. lm() has been around longer than you have, maybe, and is thus pretty well optimized. Not perfect, mind you, but I think it unlikely that casual suggestions haven't already been considered. -- Bert Gunter Genentech Non-Clinical Statistics South San Francisco, CA The business of the statistician is to catalyze the scientific learning process. - George E. P. Box -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Vivek Satsangi Sent: Wednesday, January 18, 2006 4:08 AM To: r-help@stat.math.ethz.ch Subject: [R] Possible improvement in lm Folks, I do a series of regressions (one for each quarter in the dataset) and then go and extract the residuals from each stored lm object that is returned as follows: vResiduals - as.vector(unlist(resid(lQuarterlyRegressions[[i]]))); Here lQuarterlyRegressions is a vector of objects returned by lm(). Next, I may go find outliers using identify() on a plot or do some other analysis which tells me which row of the quarterly data I need to take a closer look at. However, if I try to match some point in one of the quarters that I have with its residual, then I have to drop the points from my current Data which have NA's for either the explanatory variables or the explained, so that the vector or residuals and the data have the same indexes. This lead to some serious confusion/bugs for me, and I am wondering if it might not be better for lm to put an NA into those rows where the point was dropped because of NA's in the explanatory or explained variables (currently it just returns nothing at that index). Ofcourse, there might be some arguments against this idea, and I would be interested to hear them. Thank you for your time and attention, -- Vivek Satsangi Student, Rochester, NY USA __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] negative predicted values in poisson glm
On Wed, 18 Jan 2006 14:58:26 +0100, P. Olsson wrote: PO Dear R helpers, PO running the following code of a glm model of the family poisson, PO gives predicted values 0. Why? PO PO library(MASS) PO library(stats) PO library(mvtnorm) PO library(pscl) PO data(bioChemists) PO poisson_glm - glm(art ~ fem + mar + kid5 + phd + ment, data = PO bioChemists, family = poisson) PO predicted.values = predict(poisson_glm) PO range(predicted.values) PO use predicted.values = predict(poisson_glm, type=response) best wishes, Adelchi Azzalini [EMAIL PROTECTED] Dipart.Scienze Statistiche, Università di Padova, Italia tel. +39 049 8274147, http://azzalini.stat.unipd.it/ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] negative predicted values in poisson glm
P. Olsson olsson1 at gmail.com writes: Dear R helpers, running the following code of a glm model of the family poisson, gives predicted values 0. Why? Because by default predict.glm() gives answers on the link (log) scale. predict(poisson_glm,type=response) is probably what you're looking for. See ?predict.glm. Ben Bolker __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] linear contrasts with anova
I have some doubts about the validity of my procedure to estimeate linear contrasts ina a factorial design. For sake of semplicity, let's imagine a one way ANOVA with three levels. I am interested to test the significance of the difference between the first and third level (called here contrast C1) and between the first and the seconda level (called here contrast C2). I used the following procedure: --- reading data from a text file --- ar -read.table(C:/Programmi/R/myworks/contrasti/cont1.txt,header=TRUE) ar CC GROUP 1 3.0 0 2 3.0 0 3 4.0 0 4 5.0 0 5 6.0 0 6 7.0 0 7 3.0 0 8 2.0 0 9 1.0 1 10 6.0 1 11 5.0 1 12 7.0 1 13 2.0 1 14 3.0 1 15 1.5 1 16 1.7 1 17 17.0 2 18 12.0 2 19 15.0 2 20 16.0 2 21 12.0 2 22 23.0 2 23 19.0 2 24 21.0 2 --- creating a new array of data--- ar-data.frame(GROUP=factor(ar$GROUP),DIP=ar$CC) ar GROUP DIP 1 0 3.0 2 0 3.0 3 0 4.0 4 0 5.0 5 0 6.0 6 0 7.0 7 0 3.0 8 0 2.0 9 1 1.0 10 1 6.0 11 1 5.0 12 1 7.0 13 1 2.0 14 1 3.0 15 1 1.5 16 1 1.7 17 2 17.0 18 2 12.0 19 2 15.0 20 2 16.0 21 2 12.0 22 2 23.0 23 2 19.0 24 2 21.0 --- creating two dummy variables (C1 and C2) for linear contrasts--- ar-data.frame(GROUP=factor(ar$GROUP),C1=factor(ar$GROUP),C2=factor(ar$GROUP),DIP=ar$DIP) ar GROUP C1 C2 DIP 1 0 0 0 3.0 2 0 0 0 3.0 3 0 0 0 4.0 4 0 0 0 5.0 5 0 0 0 6.0 6 0 0 0 7.0 7 0 0 0 3.0 8 0 0 0 2.0 9 1 1 1 1.0 10 1 1 1 6.0 11 1 1 1 5.0 12 1 1 1 7.0 13 1 1 1 2.0 14 1 1 1 3.0 15 1 1 1 1.5 16 1 1 1 1.7 17 2 2 2 17.0 18 2 2 2 12.0 19 2 2 2 15.0 20 2 2 2 16.0 21 2 2 2 12.0 22 2 2 2 23.0 23 2 2 2 19.0 24 2 2 2 21.0 --- selecting the contrast levels--- ar$C1 - C(ar$C1, c(1,0,-1), how.many = 1) ar$C2 - C(ar$C2, c(1,-1,0), how.many = 1) --- contrast analysis of C2 --- r.aov8 -aov(DIP ~ C2 + GROUP , data = ar) anova(r.aov8) Analysis of Variance Table Response: DIP Df Sum Sq Mean Sq F valuePr(F) C2 1 2.102.10 0.2622 0.614 GROUP 1 917.00 917.00 114.3460 5.915e-10 *** Residuals 21 168.418.02 --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 --- contrast analysis of C1 --- r.aov9 -aov(DIP ~ C1 + GROUP , data = ar) anova(r.aov9) Analysis of Variance Table Response: DIP Df Sum Sq Mean Sq F valuePr(F) C1 1 650.25 650.25 81.083 1.175e-08 *** GROUP 1 268.85 268.85 33.525 9.532e-06 *** Residuals 21 168.418.02 --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 --- anova of the global design --- r.aov10 -aov(DIP ~ GROUP , data = ar) anova(r.aov10) Analysis of Variance Table Response: DIP Df Sum Sq Mean Sq F valuePr(F) GROUP 2 919.10 459.55 57.304 3.121e-09 *** Residuals 21 168.418.02 --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 I would like to know if there is a more economic procedure with R to do linear contrasts. Every comments will be well accepted. Thank you very much and best regards Marco Tommasi [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Influence measure + lme ?
Hi all, Does lme has function to compute the cook's distance or influence measure like lm? I can't find them. Thanks. Yen Lin [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Windows package upates
On Wed, 18 Jan 2006, John Logsdon wrote: Dear list Having just started to use the Windows version, I am very impressed with it's package handling as well as the gui. So I tried to see what was due for update and packages such as Hmisc, Matrix and others came up. But when I had updated them - which took a few goes as something hung between here and Bristol - I noticed that the default packages such as nmle, MASS had disappeared. I re-installed them but is this a glitch or a feature? Were they amongst the updates (I suspect so)? R tries hard to recover from failures to download, but it cannot do so if you kill it at a crucial time. (For source-code installs, 2.3.0 will be even better at recovering as it will catch attempts to kill and clean up.) -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Logistftest to select diagnostic genes
Hi, all, Anyone has experience on Logistf package? I am using logistftest in Logistf package to selelct diagnosis genes. The result seems not the same as I expected. I have 10 gene expression data for 27 tumor 1 and 11 tumor 0. I want to select the best one using Maximum likelihood ratio test in logistic regression model. This is the way my code works: 1. Read in 10 genes as independent variables and tumor type (1 or 0) as dependent variable. 2 Fit in a 10 variable logistic model and calculate it's likelihood. 3 Loop through the 10 genes, each time take one gene out of the model and calculate the likelihood without the gene. And compare each new likelihood with the likelihood of the original original 10 variable model to get likelihood ratio test for each gene. I guess gene 8 should have best discrimination power because it can totally seperate the two group tumors. But the test result shows the likelihood ratio for Ratio 8 is not the biggest. Where am I wrong? The package is not good for this problem? But the package document says this package takes care of separation and small sample size by finite parameter estimates and profile penalized log likelihood. Or my assumption is wrong? I am totally lost. Thank you very much for reading this email! Your suggestion or comments will be appreciated. Lingsheng The fear of the LORD is the beginning of wisdom, and knowledge of the Holy One is understanding. --Proverbs 10:10 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] linear contrasts with anova
group=factor(rep(c(0:2), each = 8)) ar = data.frame(group, dip) con = matrix(c(1, -1, 0, 1, 0, -1), nrow=3, ncol=2, byrow=F) contrasts(ar$group)=con aovRes = aov(dip~group, ar) summary.aov(aovRes, split=list(group = list(0 vs 1 = 1, 0 vs 3 = 2))) Df Sum Sq Mean Sq F valuePr(F) group2 919.10 459.55 57.3041 3.121e-09 *** group: 0 vs 1 1 2.102.10 0.2622 0.614 group: 0 vs 3 1 917.00 917.00 114.3460 5.915e-10 *** Residuals 21 168.418.02 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 I only don't know why it does not work with within-subject designs. I have posted this question before, does anybody know? -steffen I have some doubts about the validity of my procedure to estimeate linear contrasts ina a factorial design. For sake of semplicity, let's imagine a one way ANOVA with three levels. I am interested to test the significance of the difference between the first and third level (called here contrast C1) and between the first and the seconda level (called here contrast C2). I used the following procedure: --- reading data from a text file --- ar -read.table(C:/Programmi/R/myworks/contrasti/cont1.txt,header=TRUE) ar CC GROUP 1 3.0 0 2 3.0 0 3 4.0 0 4 5.0 0 5 6.0 0 6 7.0 0 7 3.0 0 8 2.0 0 9 1.0 1 10 6.0 1 11 5.0 1 12 7.0 1 13 2.0 1 14 3.0 1 15 1.5 1 16 1.7 1 17 17.0 2 18 12.0 2 19 15.0 2 20 16.0 2 21 12.0 2 22 23.0 2 23 19.0 2 24 21.0 2 --- creating a new array of data--- ar-data.frame(GROUP=factor(ar$GROUP),DIP=ar$CC) ar GROUP DIP 1 0 3.0 2 0 3.0 3 0 4.0 4 0 5.0 5 0 6.0 6 0 7.0 7 0 3.0 8 0 2.0 9 1 1.0 10 1 6.0 11 1 5.0 12 1 7.0 13 1 2.0 14 1 3.0 15 1 1.5 16 1 1.7 17 2 17.0 18 2 12.0 19 2 15.0 20 2 16.0 21 2 12.0 22 2 23.0 23 2 19.0 24 2 21.0 --- creating two dummy variables (C1 and C2) for linear contrasts--- ar-data.frame(GROUP=factor(ar$GROUP),C1=factor(ar$GROUP),C2=factor(ar$GROUP),DIP=ar$DIP) ar GROUP C1 C2 DIP 1 0 0 0 3.0 2 0 0 0 3.0 3 0 0 0 4.0 4 0 0 0 5.0 5 0 0 0 6.0 6 0 0 0 7.0 7 0 0 0 3.0 8 0 0 0 2.0 9 1 1 1 1.0 10 1 1 1 6.0 11 1 1 1 5.0 12 1 1 1 7.0 13 1 1 1 2.0 14 1 1 1 3.0 15 1 1 1 1.5 16 1 1 1 1.7 17 2 2 2 17.0 18 2 2 2 12.0 19 2 2 2 15.0 20 2 2 2 16.0 21 2 2 2 12.0 22 2 2 2 23.0 23 2 2 2 19.0 24 2 2 2 21.0 --- selecting the contrast levels--- ar$C1 - C(ar$C1, c(1,0,-1), how.many = 1) ar$C2 - C(ar$C2, c(1,-1,0), how.many = 1) --- contrast analysis of C2 --- r.aov8 -aov(DIP ~ C2 + GROUP , data = ar) anova(r.aov8) Analysis of Variance Table Response: DIP Df Sum Sq Mean Sq F valuePr(F) C2 1 2.102.10 0.2622 0.614 GROUP 1 917.00 917.00 114.3460 5.915e-10 *** Residuals 21 168.418.02 --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 --- contrast analysis of C1 --- r.aov9 -aov(DIP ~ C1 + GROUP , data = ar) anova(r.aov9) Analysis of Variance Table Response: DIP Df Sum Sq Mean Sq F valuePr(F) C1 1 650.25 650.25 81.083 1.175e-08 *** GROUP 1 268.85 268.85 33.525 9.532e-06 *** Residuals 21 168.418.02 --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 --- anova of the global design --- r.aov10 -aov(DIP ~ GROUP , data = ar) anova(r.aov10) Analysis of Variance Table Response: DIP Df Sum Sq Mean Sq F valuePr(F) GROUP 2 919.10 459.55 57.304 3.121e-09 *** Residuals 21 168.418.02 --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 I would like to know if there is a more economic procedure with R to do linear contrasts. Every comments will be well accepted. Thank you very much and best regards Marco Tommasi [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list
[R] Help with plot.svm from e1071
Hi. I'm trying to plot a pair of intertwined spirals and an svm that separates them. I'm having some trouble. Here's what I tried. library(mlbench) library(e1071) Loading required package: class raw - mlbench.spirals(200,2) spiral - data.frame(class=as.factor(raw$classes), x=raw$x[,1], y=raw$x[,2]) m - svm(class~., data=spiral) plot(m, spiral) Error in -x$index : invalid argument to unary operator So we delve into e1071:::plot.svm. When I run the code in plot.svm everything is fine up until points(formula, data = data[-x$index, ], pch = dataSymbol, col = symbolPalette[colind[-x$index]]) That gives me the same error message, Error in -x$index : invalid argument to unary operator. The weird thing is that I can run either of the those statements in isolation data[-x$index, ] symbolPalette[colind[-x$index]] and neither gives me an error. I looked in the two points functions I can see (points.default and points.formula) but neither calls x$index. I was following along the documentation for plot.svm, which has a simple example (that works) ## a simple example library(MASS) data(cats) m - svm(Sex~., data = cats) plot(m, cats) I don't see what the difference between their example and mine. Can anyone help me? Thank you, Josh. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] ICC for Binary data
Hello R users: I am fairly new to R and am trying to figure out how to compute an intraclass correlation (ICC) and/or design effect for binary data? More specifically, I am trying to determine the amount of clustering in a data set - that is, whether certain treatment programs tend to work with more or less severe clients. The outcome variable is dichotomous (low severity / high severity) and the grouping variable is the treatment program. Any suggestions for some R code or the appropriate R package would be greatly appreciated. Thanks, Brian [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] some EPS rotated in journal preview
On Wed, 2006-01-18 at 17:09 +0100, Knut Krueger wrote: Marc Schwartz schrieb: On Wed, 2006-01-18 at 13:04 +0100, Knut Krueger wrote: One of the first things to do is to use 'onefile = FALSE', 'horizontal = FALSE' and paper = special'. I am afraid the problem is on the journals side, because the wrong postscript line (with letter )is working and I changed the postscript options between both examples no change in the behaviour of the pdf creator of the journal. This is in the Details section of ?postscript, which provides guidance on the creation of EPS files for inclusion in publications. I would try that to see if that provides a more consistent formatting of the plots. You don't indicate what you are using to create the manuscript itself (ie. Word, LaTeX, or ?) to help us in considering other possibilities (such as auto rotation). The journal collects the figures once by once after the manuscript in the original file format. to solve the problem I could change it to tiff and submit it, but the tiff files are look ing not as good as eps files. If the above does not help, please provide a reproducible example of the plot code and what you are using for the manuscript. for your information: the same problem occurs if I plot the data to the graphic device and use the right button - create postscript file One is not rotated the others are rotated. Hope the code helps without data: SNIP Unfortunately, it may have precluded my being able to replicate exactly what you are seeing, despite being intimately familiar with one of the functions (barplot2 ;-) that you are using. I did create some data that would generally fit what you are doing in both cases, however, probably because I am on Linux and you are on Windows, where both device and GS differences may be problematic, I could not get the second to be rotated relative to the page. If you want, e-mail me (offlist) the rotated EPS file that is created in the second case and I can use that to try to replicate the problem as well as review the EPS code to see if something sticks out. Under Linuxen, one can set: GS_OPTIONS=-dAutoRotatePages=/None export GS_OPTIONS which disables autorotation on the system processing the EPS graphics. I would envision that there is a similar situation on Windows, using either an environment variable or via the command line, depending upon what your publisher is using for tools. A review of the Windows documentation for GS would be helpful here. Another possibility is that the rotation is being caused by some confounding in the text in your second plot. GS will base the rotation properties on the dominant text in the file. Given your plots, it may be possible that some of the text label components in the second plot are causing GS to perform the rotation. If this is the case, disabling the autorotation as above should help. The use of 'paper = special' is very helpful in general, as it enables tight bounding boxes around the EPS graphic for inclusion in other documents. I would highly recommend using that as your default. One other quick comment, which is relative to the readability of your code. Strategically placed spaces ( ) and line breaks would help tremendously. Especially since most e-mail clients will create line breaks at 72 chars (or similar), which can make the flow of the code difficult to review. HTH, Marc Schwartz __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Converting a Perl Array of Arrayrefs to an R array or matrix using RS Perl
Dear R/RS-Perl users, I have a perl script in which I parse a large number of files and construct an array of arrayrefs from the data in the files. I then pass that construct to R using the RS Perl interface. I want to be able to use the construct as an R array or matrix so that I can use the R function colSums. So far, I've tried constructing an R matrix with dummy values, and then populating each row of the new matrix with the nested Arrays of the passed construct. But that just converts the R matrix to type list, so I get an error when I try the function colSums. I then tried converting each nested Array using as.numeric() before putting the values into the matrix. But for some reason as.numeric doesn't convert the nested arrays to numeric, so I still get an error when I try colSums. If I check the values (say matrix[1][1]) the correct values are in the correct locations, but it still won't perform colSums. I've been hammering away on this for longer than I would like to admit. Any help you could provide would be greatly appreciated. Best regards, Aaron __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Bootstrapping help
Andrew, Thanks for the suggestion! This seems to have fixed things. I was wondering if you could explain why this works and what was wrong. If I issue the command my.boot-boot(dataframe,cs,R=999) and in order to what effect the command you told me use has I then do something like dataframe[my.boot$weights,] my.boot$weights looks to be a vector where element is 1/sample size.R reports that [1] var1var2 var3 var4 var5 0 rows (or 0-length row.names) Which indicates to me that I then have a dataframe with no data in it! (Am I wrong about this?) What is going on here? Why did this work? Sorry for the basic questions. Ben --- Benjamin Ridenhour School of Biological Sciences Washigton State University P.O. Box 644236 Pullman, WA 99164-4236 Phone (509)335-7218 Nothing in biology makes sense except in the light of evolution. -T. Dobzhansky - Original Message From: Andrew Robinson [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Sent: Tue Jan 17 22:20:55 2006 Subject: Re: [R] Bootstrapping help The first thing you are doing wrong is that you are not including a copy of cs for us to see ;). Based on what you have written, I speculate that cs does not use the index correctly. if so then a simple, although inefficient, workaround is to rewrite cs: cs - function(dataframe, index) { dataframe - dataframe[index,] ... } Good luck, Andrew On Tue, Jan 17, 2006 at 10:02:49PM -0800, Ben Ridenhour wrote: Hello, I am new to using R and I am having problems get boot() to work properly. Here is what I am trying to do: I have statistic called cs. cs takes a data matrix (154 x 5) and calculates 12 different scores for me. cs outputs the data as a vector (12 x 1). cs doesn't really use weights, per se, however I have included this as one of the 2 arguments cs can take. I try performing a bootstrap by issuing: myout-boot(data, cs,R=999) I have tried other versions where I specify stype=w, etc... The problem I get is that the dataset does not seem to be resampled. I end up with 999 replicates that have the exact same value of the output of cs. In the end I have something like Bootstrap Statistics : originalbiasstd. error t1* 0.865122275 1.698641e-14 0 t2* -0.005248414 -9.627715e-17 0 t3* -0.052833740 -8.812395e-16 0 t4* 0.807040121 1.287859e-14 0 t5* 0.542082588 -9.103829e-15 0 t6* -0.018617838 -7.285839e-17 0 t7* 0.006409704 1.422473e-16 0 t8* 0.529874453 8.104628e-15 0 t9* 0.074804390 2.359224e-16 0 t10* -0.007153634 1.301043e-16 0 t11* -0.018241243 -2.359224e-16 0 t12* 0.049409513 -1.200429e-15 0 Clearly the bootstrap is not working. What am I doing wrong? Thanks, Ben --- Benjamin Ridenhour School of Biological Sciences Washigton State University P.O. Box 644236 Pullman, WA 99164-4236 Phone (509)335-7218 Nothing in biology makes sense except in the light of evolution. -T. Dobzhansky [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Andrew Robinson Department of Mathematics and StatisticsTel: +61-3-8344-9763 University of Melbourne, VIC 3010 Australia Fax: +61-3-8344-4599 Email: [EMAIL PROTECTED] http://www.ms.unimelb.edu.au [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Help with plot.svm from e1071
Hi. I'm trying to plot a pair of intertwined spirals and an svm that separates them. I'm having some trouble. Here's what I tried. library(mlbench) library(e1071) Loading required package: class raw - mlbench.spirals(200,2) spiral - data.frame(class=as.factor(raw$classes), x=raw$x[,1], y=raw$x[,2]) m - svm(class~., data=spiral) plot(m, spiral) Error in -x$index : invalid argument to unary operator So we delve into e1071:::plot.svm. When I run the code in plot.svm everything is fine up until points(formula, data = data[-x$index, ], pch = dataSymbol, col = symbolPalette[colind[-x$index]]) That gives me the same error message, Error in -x$index : invalid argument to unary operator. The weird thing is that I can run either of the those statements in isolation data[-x$index, ] symbolPalette[colind[-x$index]] and neither gives me an error. I looked in the two points functions I can see (points.default and points.formula) but neither calls x$index. I was following along the documentation for plot.svm, which has a simple example (that works) ## a simple example library(MASS) data(cats) m - svm(Sex~., data = cats) plot(m, cats) I don't see what the difference between their example and mine. Can anyone help me? Thank you, Josh. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] some EPS rotated in journal preview
Marc Schwartz (via MN) schrieb: Unfortunately, it may have precluded my being able to replicate exactly what you are seeing, despite being intimately familiar with one of the functions (barplot2 ;-) that you are using. I did create some data that would generally fit what you are doing in both cases, however, probably because I am on Linux and you are on Windows, where both device and GS differences may be problematic, I could not get the second to be rotated relative to the page. I think I described the problem not clearly or misunderstood your text. There was no Figure rotated on my system - only in the PDF file from the jurnal. If you want, e-mail me (offlist) the rotated EPS file that is created in the second case and I can use that to try to replicate the problem as well as review the EPS code to see if something sticks out. Under Linuxen, one can set: GS_OPTIONS=-dAutoRotatePages=/None export GS_OPTIONS ... I will try to findout this on windows system, but maybe I need some time. Fist I must send the paper out to the journal. I will attach both Tiff for the reviewers and eps for printing . One other quick comment, which is relative to the readability of your code. Strategically placed spaces ( ) and line breaks would help tremendously. Especially since most e-mail clients will create line breaks at 72 chars (or similar), which can make the flow of the code difficult to review. sure I only used copy and paste from the sriptfile and did noth thought about the autowrap fo the e-mailclient. thanks very much Knut __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] some EPS rotated in journal preview
On Wed, 2006-01-18 at 20:48 +0100, Knut Krueger wrote: Marc Schwartz (via MN) schrieb: Unfortunately, it may have precluded my being able to replicate exactly what you are seeing, despite being intimately familiar with one of the functions (barplot2 ;-) that you are using. I did create some data that would generally fit what you are doing in both cases, however, probably because I am on Linux and you are on Windows, where both device and GS differences may be problematic, I could not get the second to be rotated relative to the page. I think I described the problem not clearly or misunderstood your text. There was no Figure rotated on my system - only in the PDF file from the jurnal. That is correct. The rotation occurs during the post processing of the EPS file in the full document (or single page) when using latex+dvips +ps2pdf (or similar tool chain). Not typically the standalone EPS file itself. Regards, Marc __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Breakpoints for multiple variables using Segmented
Hi all, I am using the package Segmented to estimate logistic regression models with unknown breakpoints (see Muggeo 2003 Statistics in Medicine 22:3055-3071). In the documentation it suggests that it might be possible to include several variables with breakpoints in the same model: Z = a vector or a matrix meaning the (continuous) explanatory variable(s) having segmented relationships with the response. However, the syntax for including multiple Z and psi (starting values for the break-point(s)) is not stated. Does anyone have any suggestions? Here is an example of correct code for detecting single breakpoint: model.seg-segmented.glm(obj = model.glm, Z = predictor_variable, psi = 2 , it.max = 50) Thanks very much for your help. Matthew G. Betts, Ph.D. NB Cooperative Fish and Wildlife Research Unit Faculty of Forestry and Environmental Management University of New Brunswick UNB Tweedale Centre Hugh John Flemming Forestry Complex 1350 Regent St., Fredericton, N.B. E3C 2G6 (506) 447-3408 http://www.unb.ca/web/acwern/people/mbetts/mbetts.htm __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Plotting an lme( ) object
R for Windows, version 2.2. I am trying the following code: BLAH.lme4 - lme(fixed = , data = , random = ) plot(BLAH.lme4, resid(.,type = p) ~ fitted(.) | GROUP, id = 0.05, adj = -0.3, idLabels = BLAH$VALUE, main = Pearson Residuals vs. Fitted Values, by Group) When I use plot( ), it generates a nice plot. However, no matter what I use for adj = ... I cannot get the observation labels to shift up, down, left or right. How can I do this? I am trying to follow along with the example in Pinheiro and Bates, page 176. Thanks much, Greg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] r-help, how can i use my own distance matrix without using dist()
Dear R-helpers, i am a beginner of R and i am using cluster package to do hierarchical clustering i am wondering if i can use my own distance matrix to do the hierarchical clustering without using dist() function. if i have my own distance matrix, how can i ask hclust() function to recongnize it( as the output of dist() function). thank you very much and i looking forward to hearing from you. Marshall __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] r-help, how can i use my own distance matrix without usin g dist()
Use something like hclust(as.dist(mydist), ...) ought to work. Andy -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of [EMAIL PROTECTED] Sent: Wednesday, January 18, 2006 4:47 PM To: r-help@stat.math.ethz.ch Subject: [R] r-help, how can i use my own distance matrix without using dist() Dear R-helpers, i am a beginner of R and i am using cluster package to do hierarchical clustering i am wondering if i can use my own distance matrix to do the hierarchical clustering without using dist() function. if i have my own distance matrix, how can i ask hclust() function to recongnize it( as the output of dist() function). thank you very much and i looking forward to hearing from you. Marshall __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] do.call, missing arguments and [
x - array(1:30, c(4,5,3)) x[1,,] How can I do the same thing with do.call? do.call([, list(x, 1)) == x[1] do.call([, list(x, 1, NULL, NULL)) == x[1, NULL, NULL] I guess you can't, because of the special way that [ deals with arguments. How can I index an array programmatically for arrays and indices of varying dimensionality? Thanks, Hadley __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Plotting an lme( ) object
I apologize for the second posting, my other email address died today. I am using R for Windows, version 2.2. Here is my code: plot(FOO.lme4, resid(.,type = p) ~ fitted(.) | LOT, id = 0.05, adj = -0.3, idLabels = RO54F$VALUE, main = Pearson Residuals vs. Fitted Values, by Lot, between = list(x = .5, y = .5)) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] r-help, how can i use my own distance matrix without using dist()
see ?dist there's an example x - matrix(rnorm(100), nrow=5) m - as.matrix(dist(x)) d - as.dist(m) 'as.dist' is what you're probably looking for regards, Marco [EMAIL PROTECTED] wrote: Dear R-helpers, i am a beginner of R and i am using cluster package to do hierarchical clustering i am wondering if i can use my own distance matrix to do the hierarchical clustering without using dist() function. if i have my own distance matrix, how can i ask hclust() function to recongnize it( as the output of dist() function). thank you very much and i looking forward to hearing from you. Marshall __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html - Photo Books. You design it and well bind it! [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Plotting an lme( ) object
I apologize for the second posting, my other email address died today. I am using R for Windows, version 2.2. Here is my code: plot(FOO.lme4, resid(.,type = p) ~ fitted(.) | GROUP, id = 0.05, adj = -0.3, idLabels = FOO$value, main = Pearson Residuals vs. Fitted Values, by Group, between = list(x = .5, y = .5)) The plot looks fine, but the adj = -0.3 option seems to have no effect on the labels that are added to identify potential outliers. I would like to offset the FOO$value text that is currently being displayed right on top of several pearson residuals. How can I do this? Thanks much, Greg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] do.call, missing arguments and [
On 1/18/06, hadley wickham [EMAIL PROTECTED] wrote: x - array(1:30, c(4,5,3)) x[1,,] How can I do the same thing with do.call? do.call([, list(x, 1, TRUE, TRUE)) seems to work for me. Deepayan do.call([, list(x, 1)) == x[1] do.call([, list(x, 1, NULL, NULL)) == x[1, NULL, NULL] I guess you can't, because of the special way that [ deals with arguments. How can I index an array programmatically for arrays and indices of varying dimensionality? __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Canonical Variance Analysis by any other name?
I've been asked about Canonical Variance Analysis (CVA). I don't see any reference to it searching the R site. Does it go by other names? Genstat describes it thus: Canonical variates analysis operates on a within-group sums of squares and products matrix, calculated from a set of variates and factor that specifies the grouping of units. It finds linear combinations of the variates that maximize the ratio of between-group to within-group variation, thereby giving functions that can be used to discriminate between the groups. It's probably not particularly difficult to do, so I suspect someone has a package for doing it. What other name might I search for? Thnx -- 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 ~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Bootstrapping help
Ben, although I appended a smiley to my first note, the message was serious. If you don't show us what you're doing, we can't help you. Please provide an example in which you: 1) generate a small dataframe similar in structure to yours 2) provide cs 3) show the boot statement that applies cs to the example dataframe. Also, it seems that you are unfamiliar with the use of indexing and datframes. Please read the Introduction to R, carefully, it is freely available on CRAN. You have asked R to provide you with all the rows that are numbered 1/sample size.R; since the row numbers are integers there aren't any. And, please say hello to Andrew Storfer and Melanie Murphy from me. Andrew On Wed, Jan 18, 2006 at 11:39:36AM -0800, Ben Ridenhour wrote: Andrew, Thanks for the suggestion! This seems to have fixed things. I was wondering if you could explain why this works and what was wrong. If I issue the command my.boot-boot(dataframe,cs,R=999) and in order to what effect the command you told me use has I then do something like dataframe[my.boot$weights,] my.boot$weights looks to be a vector where element is 1/sample size.R reports that [1] var1var2 var3 var4 var5 0 rows (or 0-length row.names) Which indicates to me that I then have a dataframe with no data in it! (Am I wrong about this?) What is going on here? Why did this work? Sorry for the basic questions. Ben -- Andrew Robinson Department of Mathematics and StatisticsTel: +61-3-8344-9763 University of Melbourne, VIC 3010 Australia Fax: +61-3-8344-4599 Email: [EMAIL PROTECTED] http://www.ms.unimelb.edu.au __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Bootstrapping help
Thanks for responding :) Again... I understand how indexing works (basically as in any other programming language), that is why I am so confused as to why that statement made my bootstraping work! It seems like, if anything, it would completely screw up everything. Here is the (now working) cs function after I amended it to what you said to do (i.e. I added the _very confusing_ statement data-data[w,]): cs-function(data, w){ data-data[w,]; y-data[1]; x1-data[2]; x2-data[3]; x3-data[4]; z-data[5]; c1-x1*z; c2-x2*z; c3-x3*z; X-cbind(x1,x2,x3,z,c1,c2,c3); regcoef-lsfit(X,y)$coefficients; bx1-regcoef[[2]]; bx2-regcoef[[3]]; bx3-regcoef[[4]]; bz-regcoef[[5]]; bc1-regcoef[[6]]; bc2-regcoef[[7]]; bc3-regcoef[[8]]; fx-bx1*x1+bx2*x2+bx3*x3; gy-bz*z; hxy-bc1*c1+bc2*c2+bc3*c3; sfx1-cov(fx,x1); sfx2-cov(fx,x2); sfx3-cov(fx,x3); sgx1-cov(gy,x1); sgx2-cov(gy,x2); sgx3-cov(gy,x3); shx1-cov(hxy,x1); shx2-cov(hxy,x2); shx3-cov(hxy,x3); sTx1-cov(y,x1); sTx2-cov(y,x2); sTx3-cov(y,x3); dataout-c(sfx1,sgx1,shx1,sTx1,sfx2,sgx2,shx2,sTx2,sfx3,sgx3,shx3,sTx3); dataout } An example data frame would be mydata-data.frame(Y=rnorm(20,0,1),X1=rnorm(20,0,1),X2=rnorm(20,0,1),X3=rnorm(20,0,1),Z=rnorm(20,0,1)) The boot statement is boot(mydata,cs,R=999) Why does this rather mysterious indexing statement data-data[w,] make the bootstrap work when it didn't beforehand? Thanks, Ben ps. I'll tell Melanie and Storfer hello. --- Benjamin Ridenhour School of Biological Sciences Washigton State University P.O. Box 644236 Pullman, WA 99164-4236 Phone (509)335-7218 Nothing in biology makes sense except in the light of evolution. -T. Dobzhansky - Original Message From: Andrew Robinson [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Sent: Wednesday, January 18, 2006 15:45:15 Subject: Re: [R] Bootstrapping help Ben, although I appended a smiley to my first note, the message was serious. If you don't show us what you're doing, we can't help you. Please provide an example in which you: 1) generate a small dataframe similar in structure to yours 2) provide cs 3) show the boot statement that applies cs to the example dataframe. Also, it seems that you are unfamiliar with the use of indexing and datframes. Please read the Introduction to R, carefully, it is freely available on CRAN. You have asked R to provide you with all the rows that are numbered 1/sample size.R; since the row numbers are integers there aren't any. And, please say hello to Andrew Storfer and Melanie Murphy from me. Andrew On Wed, Jan 18, 2006 at 11:39:36AM -0800, Ben Ridenhour wrote: Andrew, Thanks for the suggestion! This seems to have fixed things. I was wondering if you could explain why this works and what was wrong. If I issue the command my.boot-boot(dataframe,cs,R=999) and in order to what effect the command you told me use has I then do something like dataframe[my.boot$weights,] my.boot$weights looks to be a vector where element is 1/sample size.R reports that [1] var1var2 var3 var4 var5 0 rows (or 0-length row.names) Which indicates to me that I then have a dataframe with no data in it! (Am I wrong about this?) What is going on here? Why did this work? Sorry for the basic questions. Ben -- Andrew Robinson Department of Mathematics and StatisticsTel: +61-3-8344-9763 University of Melbourne, VIC 3010 Australia Fax: +61-3-8344-4599 Email: [EMAIL PROTECTED] http://www.ms.unimelb.edu.au [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Bootstrapping help
Ben, Ok, it's clear now, thanks. Note that your boot call boot(mydata,cs,R=999) does not specify an stype argument. The boot help file notes that the default value for stype is i, which means that boot will pass an index to the function, not a weight, regardless of whether you call it w, i, or whatever. The index that boot sends to the function is then used to index the dataframe, thus selecting rows randomly with replacement. Previously you passed the dataframe to the function, which did not alter it, so it passed through undisturbed. In this incarnation the data-data[w,] command provides you with the (pseudo-)random sample with replacement of the data. I hope that this clears up the confusion. Cheers, Andrew ps it's always good to provide a brief bit of sample code when you ask a question. Also, let me recommend that you omit semi-colons and space the code to make it easier to read. Thus cs - function(data, w) { data-data[w, ] ... On Wed, Jan 18, 2006 at 04:35:47PM -0800, Ben Ridenhour wrote: Thanks for responding :) Again... I understand how indexing works (basically as in any other programming language), that is why I am so confused as to why that statement made my bootstraping work! It seems like, if anything, it would completely screw up everything. Here is the (now working) cs function after I amended it to what you said to do (i.e. I added the _very confusing_ statement data-data[w,]): cs-function(data, w){ data-data[w,]; y-data[1]; x1-data[2]; x2-data[3]; x3-data[4]; z-data[5]; c1-x1*z; c2-x2*z; c3-x3*z; X-cbind(x1,x2,x3,z,c1,c2,c3); regcoef-lsfit(X,y)$coefficients; bx1-regcoef[[2]]; bx2-regcoef[[3]]; bx3-regcoef[[4]]; bz-regcoef[[5]]; bc1-regcoef[[6]]; bc2-regcoef[[7]]; bc3-regcoef[[8]]; fx-bx1*x1+bx2*x2+bx3*x3; gy-bz*z; hxy-bc1*c1+bc2*c2+bc3*c3; sfx1-cov(fx,x1); sfx2-cov(fx,x2); sfx3-cov(fx,x3); sgx1-cov(gy,x1); sgx2-cov(gy,x2); sgx3-cov(gy,x3); shx1-cov(hxy,x1); shx2-cov(hxy,x2); shx3-cov(hxy,x3); sTx1-cov(y,x1); sTx2-cov(y,x2); sTx3-cov(y,x3); dataout-c(sfx1,sgx1,shx1,sTx1,sfx2,sgx2,shx2,sTx2,sfx3,sgx3,shx3,sTx3 ); dataout } An example data frame would be mydata-data.frame(Y=rnorm(20,0,1),X1=rnorm(20,0,1),X2=rnorm(20,0,1) ,X3=rnorm(20,0,1),Z=rnorm(20,0,1)) The boot statement is boot(mydata,cs,R=999) Why does this rather mysterious indexing statement data-data[w,] make the bootstrap work when it didn't beforehand? Thanks, Ben ps. I'll tell Melanie and Storfer hello. --- Benjamin Ridenhour School of Biological Sciences Washigton State University P.O. Box 644236 Pullman, WA 99164-4236 Phone (509)335-7218 Nothing in biology makes sense except in the light of evolution. -T. Dobzhansky http://www.ms.unimelb.edu.au -- Andrew Robinson Department of Mathematics and StatisticsTel: +61-3-8344-9763 University of Melbourne, VIC 3010 Australia Fax: +61-3-8344-4599 Email: [EMAIL PROTECTED] http://www.ms.unimelb.edu.au __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] do.call, missing arguments and [
do.call([, list(x, 1, TRUE, TRUE)) seems to work for me. Oh! Of course. I knew someone obvious was staring me in the face. Thanks, Hadley __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Powell's unconstrained derivative-free nonlinear least squares routine, VA05AD
G'day David, DK == David Kinniburgh [EMAIL PROTECTED] writes: DK I have used Mike Powell's optimization routine (VA05AD) from DK the Harwell Subroutine Library (HSL) for more than 20 years. DK [...] DK Now that I have converted to R, I will miss my trusted DK friend. I have started using nls() but have not accumulated DK enough experience to compare the two. It would be great if DK VA05AD could be an option in there (algorithm=Powell). To DK this end, I recently enquired of the custodians of the HSL DK whether it would be possible to make it freely available to DK the R community. The answer was basically 'yes'. [...] You better ask again. :-) I presume VA05AD is in what is now called the HSL archive? If it is part of HSL 2004, then I don't see a way how it could be incorporated into R given the information on their web site. Even for the HSL archive, it is stated at http://hsl.rl.ac.uk/archive/hslarchive.html that [...] HSL Archive codes are now available without charge to anyone, so long as they are not then incorporated into a commercial product; the latter is, of course, still permitted subject to a commercial licence. [...] which sounds as if it would be possible to interface (parts of) HSL archive to R. But then that URL goes on ans states: Access to the Archive is by means of a short-lived individual password-controlled account. Potential users are asked for brief details of the use they intend to make of the package(s) they aim to download. Users are also asked to accept a conditions-of-use form, and are not permitted to divulge their userid and password to anyone else, nor to distribute any codes they obtain to a third party. IANAL, but as I understand the GPL, the last part of that sentence will make it pretty impossible to incorporate into R (or distribute as an R package) (parts of) the HSL archive routines. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Mathematics and Statistics+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: [EMAIL PROTECTED] Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] mcmcsamp() in lmer
I am working with lmer() in the latest release of Matrix, doing various things including writing a function called mcsamp() that acts as a wrapper for mcmcsamp() and automatically runs multiple chains, diagnoses convergence, and stores the result as a bugs object so it can be plotted. I recognize that at this point, mcmcsamp() is somewhat of a placeholder (since it doesn't work on a lot of models) but I'm sure it will continue to be improved so I'd like to be able to work with it, as a starting point if necessary. Anyway, I couldn't get mcmcsamp() to work with the saveb=TRUE option. Here's a simple example: y - 1:10 group - rep (c(1,2), c(5,5)) M1 - lmer (y ~ 1 + (1 | group)) # works fine mcmcsamp (M1) # works fine mcmcsamp (M1, saveb=TRUE) This last gives an error message: Error in colnames-(`*tmp*`, value = c((Intercept), log(sigma^2), : length of 'dimnames' [2] not equal to array extent Thanks for your help. Andrew -- Andrew Gelman Professor, Department of Statistics Professor, Department of Political Science [EMAIL PROTECTED] www.stat.columbia.edu/~gelman Tues, Wed, Thurs: Social Work Bldg (Amsterdam Ave at 122 St), Room 1016 212-851-2142 Mon, Fri: International Affairs Bldg (Amsterdam Ave at 118 St), Room 711 212-854-7075 Mailing address: 1255 Amsterdam Ave, Room 1016 Columbia University New York, NY 10027-5904 212-851-2142 (fax) 212-851-2164 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] glmmPQL: Na/NaN/Inf in foreign function call
I have not seen a reply to this post, and unfortunately, I'm not smart enough to help you. If you would still like help from this listserve, please provide a terse, reproducible example that someone can in seconds copy from your email into R and presumably get the same error message you report. If you have that, then anyone with interest in glmmPQL can try to help you. Without that, your email must reach someone who has gotten that specific error message from glmmPQL with sufficient regularity that the error message has become solidly recorded in their brains. Moreover, this person must also have the time and interest to reply. Don't gamble on a long shot. Play the odds. PLEASE do read the posting guide! www.R-project.org/posting-guide.html. spencer graves David Reitter wrote: I'm using glmmPQL, and I still have a few problems with it. In addition to the issue reported earlier, I'm getting the following error and I was wondering if there's something I can do about it. Error in logLik.reStruct(object, conLin) : Na/NaN/Inf in foreign function call (arg 3) ... Warnings: 1: Singular precistion matrix in level -1, block 4 (...) 4: The interaction terms are primed ~ log(dist) * role random = ~ dist | target.utt / prime.utt The family is binomial (logit). I ensured that log(dist) is always [0; 1]. Role is a factor (binary). target.utt and prime.utt are categories as well. This is version 2.1.1 - is the latest version more reliable? Thanks for any help you can give me. Dave __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] numericDeriv() giving a vector when multiple variables input
R Help List -- I have defined two time-series-vector-valued-functions, let them be f and g, and want to find the numeric derivative of f with respect to the variable x where f depends on x through g: (d/dx)(f (g(x) ) Moreover, x is a vector I tried this out the long way (naming every element of the x vector and then making the 'theta' argument in numericDeriv() the character vector of all these names) and the result is just one time series vector; I was hoping for a matrix. Also weirdly, if I instead make theta equal to just one of the named elements of x, I get the same time series vector; the same happens for any subset of the named elements My call to numericDeriv looks like this (vphi acts as x, swz.kalman.vectoracts as g, decision.ts.vector acts as f): *** numericDeriv( expr = decision.ts.vector( a = swz.kalman.vector( zeta1=zeta[1], u = Phi$u, phi = Phi$Phi, vphi = c(varphi1,varphi2,varphi3,varphi4,varphi5,varphi6,varphi7,varphi8,varphi9,varphi10, varphi11,varphi12,varphi13,varphi14,varphi15,varphi16,varphi17,varphi18,varphi19, varphi20,varphi21,varphi22,varphi23,varphi24,varphi25,varphi26,varphi27,varphi28, varphi29,varphi30,varphi31,varphi32,varphi33,varphi34,varphi35,varphi36,varphi37, varphi38,varphi39,varphi40,varphi41,varphi42), alpha.prior = specs$alpha.prior ), phi = Phi$Phi, lambda = specs$lambda, delta = specs$delta, pi.star = specs$pi.star, u.2star = specs$u.2star ), theta = c(varphi1,varphi2,varphi3,varphi4,varphi5,varphi6,varphi7,varphi8,varphi9,varphi10, varphi11,varphi12,varphi13,varphi14,varphi15,varphi16,varphi17,varphi18,varphi19, varphi20,varphi21,varphi22,varphi23,varphi24,varphi25,varphi26,varphi27,varphi28, varphi29,varphi30,varphi31,varphi32,varphi33,varphi34,varphi35,varphi36,varphi37, varphi38,varphi39,varphi40,varphi41,varphi42) ); *** As you can see, it includes some calls to other objects that are lying around. Maybe from this email someone can tell me my mistake with how I've called things; otherwise, I'm happy to send along the .R files to whomever is so kind as offer guidance. I appreciate your time, Seth -- Seth Pruitt Department of Economics University of California, San Diego [EMAIL PROTECTED] http://dss.ucsd.edu/~sjpruitt [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Legend Outside Plot Dimension
Dear All, I'm trying to attach a legend outside the plot (Inside plot OK), but failed. Any help is very much appreciated. Thanks. Abd. Rahman Kassim, PhD Forest Management Ecology Program Forestry Conservation Division Forest Research Institute Malaysia Kepong 52109 Selangor MALAYSIA * Checked by TrendMicro Interscan Messaging Security. For any enquiries, please contact FRIM IT Department. * [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Homogenous groups building - Randomisation
I'm sorry, but I still do not understand, but I will make a guess: You have data on n characteristics from each of many patients, and you want to form N groups, with roughly k patients per group from the next roughly k*N patients who walk into your clinic. If this is your problem, and if you have a training set of data you can use to define group selection criteria, I suggest you try various versions of cluster analysis, including hclust in the stats package and other functions described with See also in the hclust documentation. If you don't have this and would still like help from this listserve, please think very carefully about what it would take for someone to produce a useful response to your question. Je voudrais vous aidez, mais ce demand plus que j'ai. Hope this helps, spencer graves [EMAIL PROTECTED] wrote: Dear R-users, We expect to form N homogeneous groups of n features from an experimentation including N*n data. The aim is to prevent group effects. How to do that with R functionalitites ? Does anyone know any way enabling this ? Example : 100 patients are observed. 3 biochemical parameters are mesured for each one (Red and white globules ans glycemia). Patient RG RW Gly 1 3.4 1.38 1.62 2 1.8 1.19 1.55 3 1.9 1.26 1.77 4 3.0 1.29 1.72 5 1.9 1.09 1.72 6 3.3 1.31 1.63 ... ... ... ... These are freakish data. How to compute 10 equivalent groups ? For example with closed parameters means between groups : Group RGmean RWmean Glymean 1 1.5 1.22 1.68 2 1.3 1.29 1.75 3 1.6 1.25 1.63 4 1.2 1.23 1.70 ... ... ... ... Best regards. Alexandre MENICACCI Bioinformatics - FOURNIER PHARMA 50, rue de Dijon - 21121 Daix - FRANCE [EMAIL PROTECTED] tél : 03.80.44.76.17 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] creating objects with a slot of class formula, using new
Hi, This works fine. setClass(a,representation(b=list,c=list)) new(a,b=list(7)) An object of class a Slot b: [[1]] [1] 7 Slot c: list() But, now suppose you want a slot to accept an object of class formula... setClass(a,representation(b=list,c=formula)) new(a,b=list(7)) Error in validObject(.Object) : invalid class a object: invalid object for slot c in class a: got class NULL, should be or extend class formula Why can't new handle this? Why must the slot be defined? If I call new without any named arguments, it works fine new(a) An object of class a Slot b: list() Slot c: NULL If I call new with only a formula, it works fine. new(a,c=formula(x~y)) An object of class a Slot b: list() Slot c: x ~ y How can I get R to do this? setClass(a,representation(b=list,c=formula)) new(a,b=list(7)) An object of class a Slot b: [[1]] [1] 7 Slot c: NULL Thanks, Steve platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major2 minor1.1 year 2005 month06 day 20 language R [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] singular convergence(7)?
In general, singular convergence means you were trying to estimate k parameters when the data would support estimation k-1 or fewer. Beyond this, RSiteSearch(singlular convergence) produce 53 hits for me just now, and RSiteSearch(lme singular convergence) produced 7. If none of those answer your question and you would still like help from this listserve, please submit a terse, reproducible example consistent with www.R-project.org/posting-guide.html. You can increase your chances of getting a useful reply by making it easier for potential respondents to comment. With a simple, terse, reproducible example, anyone interested in lme can copy your code into R and see your error message in seconds. With a slightly greater effort, they might be able to tell you exactly what you want to know. Without such a simple, reproducible example, you are throwing darts blindly, hoping one will find the bulls-eye of a target you only hope is there. hope this helps. spencer graves Chia, Yen Lin wrote: Hi all, I just wonder what singular convergence means. Thanks. Yen Lin Error in lme.formula(Data ~ 1, random = ~1 | Wafer/fie/loc, subset = Wafer == : singular convergence (7) [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Legend Outside Plot Dimension
use xpd argument in par(), as follows: ?par par(xpd=T, mar=par()$mar+c(0,0,0,4)) plot(1,1) legend(1.5,1,point,pch=1) Abd Rahman Kassim a écrit : Dear All, I'm trying to attach a legend outside the plot (Inside plot OK), but failed. Any help is very much appreciated. Thanks. Abd. Rahman Kassim, PhD Forest Management Ecology Program Forestry Conservation Division Forest Research Institute Malaysia Kepong 52109 Selangor MALAYSIA * Checked by TrendMicro Interscan Messaging Security. For any enquiries, please contact FRIM IT Department. * [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] creating objects with a slot of class formula, using new
Create a virtual class that can either be a formula or be NULL: setClassUnion(formulaOrNULL, c(formula, NULL)) setClass(a,representation(b=list,c=formulaOrNULL)) new(a, b = list(7)) On 1/18/06, Steven Lacey [EMAIL PROTECTED] wrote: Hi, This works fine. setClass(a,representation(b=list,c=list)) new(a,b=list(7)) An object of class a Slot b: [[1]] [1] 7 Slot c: list() But, now suppose you want a slot to accept an object of class formula... setClass(a,representation(b=list,c=formula)) new(a,b=list(7)) Error in validObject(.Object) : invalid class a object: invalid object for slot c in class a: got class NULL, should be or extend class formula Why can't new handle this? Why must the slot be defined? If I call new without any named arguments, it works fine new(a) An object of class a Slot b: list() Slot c: NULL If I call new with only a formula, it works fine. new(a,c=formula(x~y)) An object of class a Slot b: list() Slot c: x ~ y How can I get R to do this? setClass(a,representation(b=list,c=formula)) new(a,b=list(7)) An object of class a Slot b: [[1]] [1] 7 Slot c: NULL Thanks, Steve platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major2 minor1.1 year 2005 month06 day 20 language R [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] sapply on data frames - $X vs.[X]
I have a program where the following code works fine: df.opt$Delta - apply(df.opt, 1, function(x) EuropeanOption(x[CallPutString], x[UnderlyingPrice], x[StrikePrice], 0, interest.rate, x[FractionalExpiration], x[Vol])[[2]]) However, the code below fails: #df.opt$Delta - apply(df.opt, 1, function(x) EuropeanOption(x$CallPutString, #x[UnderlyingPrice], x[StrikePrice], 0, interest.rate, x[FractionalExpiration], #x[Vol])[[2]]) Here df.opt is a data frame and European Option is a function in the package RQuantLib. The only difference between the two pieces of code is that a column of the data frame is accessed using $ in the first case and brackets in the second case. Sorry if this is a trivial question, but could somebody tell me the reason for this behavior? FS __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] sapply on data frames - $X vs.[X]
with apply loop, x is a vector you are indexing in the first sample code with its names. so it failed when you tried indexing it with $ which is dedicated to data frames. Fernando Saldanha a écrit : I have a program where the following code works fine: df.opt$Delta - apply(df.opt, 1, function(x) EuropeanOption(x[CallPutString], x[UnderlyingPrice], x[StrikePrice], 0, interest.rate, x[FractionalExpiration], x[Vol])[[2]]) However, the code below fails: #df.opt$Delta - apply(df.opt, 1, function(x) EuropeanOption(x$CallPutString, #x[UnderlyingPrice], x[StrikePrice], 0, interest.rate, x[FractionalExpiration], #x[Vol])[[2]]) Here df.opt is a data frame and European Option is a function in the package RQuantLib. The only difference between the two pieces of code is that a column of the data frame is accessed using $ in the first case and brackets in the second case. Sorry if this is a trivial question, but could somebody tell me the reason for this behavior? FS __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] sapply on data frames - $X vs.[X]
$ can be used with data frames but not arrays and apply is not sending a data frame to your function. For example, try this: apply(df.opt, 1, is.data.frame) By the way, in addition the difference you mention there are a bunch of # marks in your second statement. On 1/19/06, Fernando Saldanha [EMAIL PROTECTED] wrote: I have a program where the following code works fine: df.opt$Delta - apply(df.opt, 1, function(x) EuropeanOption(x[CallPutString], x[UnderlyingPrice], x[StrikePrice], 0, interest.rate, x[FractionalExpiration], x[Vol])[[2]]) However, the code below fails: #df.opt$Delta - apply(df.opt, 1, function(x) EuropeanOption(x$CallPutString, #x[UnderlyingPrice], x[StrikePrice], 0, interest.rate, x[FractionalExpiration], #x[Vol])[[2]]) Here df.opt is a data frame and European Option is a function in the package RQuantLib. The only difference between the two pieces of code is that a column of the data frame is accessed using $ in the first case and brackets in the second case. Sorry if this is a trivial question, but could somebody tell me the reason for this behavior? FS __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] creating objects with a slot of class formula, using new
On 18 Jan 2006, [EMAIL PROTECTED] wrote: But, now suppose you want a slot to accept an object of class formula... setClass(a,representation(b=list,c=formula)) new(a,b=list(7)) Error in validObject(.Object) : invalid class a object: invalid object for slot c in class a: got class NULL, should be or extend class formula Why can't new handle this? Why must the slot be defined? Because 'formula' is not a formal S4 class. new(formula) doesn't work. And because isVirtualClass(getClass(formula)) returns TRUE, which I find surprising since one can have an instance of a formula. One workaround is to define a prototype: setClass(a, representation(b=list, c=formula), prototype=prototype(c=formula(~1))) + seth __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Legend Outside Plot Dimension
Dear Jacques, Thanks for the promt response. Abd. Rahman - Original Message - From: Jacques VESLOT [EMAIL PROTECTED] To: Abd Rahman Kassim [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Sent: Wednesday, January 18, 2006 9:15 PM Subject: Re: [R] Legend Outside Plot Dimension use xpd argument in par(), as follows: ?par par(xpd=T, mar=par()$mar+c(0,0,0,4)) plot(1,1) legend(1.5,1,point,pch=1) Abd Rahman Kassim a écrit : Dear All, I'm trying to attach a legend outside the plot (Inside plot OK), but failed. Any help is very much appreciated. Thanks. Abd. Rahman Kassim, PhD Forest Management Ecology Program Forestry Conservation Division Forest Research Institute Malaysia Kepong 52109 Selangor MALAYSIA * Checked by TrendMicro Interscan Messaging Security. For any enquiries, please contact FRIM IT Department. * [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html * Checked by TrendMicro Interscan Messaging Security. For any enquiries, please contact FRIM IT Department. * * Checked by TrendMicro Interscan Messaging Security. For any enquiries, please contact FRIM IT Department. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html