Re: [R] Analog to matlab repmat function
On Mon, 29 May 2006, Daniil Ivanov wrote: Oops, mistake, should be data = data - repmat(mean(data,1),[130 1 1 1]); Sorry. I am afraid I do not know what that actually does in matlab, except that it looks like subtracting (not extracting) a mean. Functions which I would expect you to be using are sweep(), rowMeans or colMeans(). Something like arr - sweep(arr, 2:4, colMeans(arr, 1)) perhaps? (Using 'data', the name of a system R object, is a good way to sidetrack your audience and so best avoided.) On 5/29/06, Daniil Ivanov [EMAIL PROTECTED] wrote: Hello, I'm trying to switch from Matlab to R-project, and having some difficulties. I make a use of multidimensional matrices. For example, I need to extract mean from one of the dimensions: % we have matrix data of size: 130 x 11 x 350 x 2 data = data - repmat(mean(data,3),[130 1 1 1]); In R project I managed to do that in a very pervarsive way: # mean(data,3) in R base - apply(data,c(2,3,4),mean) Probably rowMeans(data, 1) is meant. # repmat(...,[130 1 1 1]) base - rep(base,130) dim(base) - c(11,2,350,130) base - aperm(base,c(4,1,2,3)) # data = data - repmat data - data - base; Could you please show me a more elegant way of doing the same in R :) Thanks, Daniil. __ R-help@stat.math.ethz.ch 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] combinatorial programming problem
SpG == Spencer Graves [EMAIL PROTECTED] on Sun, 28 May 2006 16:21:53 -0700 writes: SpG I'm not sure I understand your question, but SpG are you asking how to index choose(k, r) objects? SpG Almost 3 years ago, I asked a question like this. Andy SpG Liaw referred me to nchoosek(vsn) SpG (http://finzi.psych.upenn.edu/R/Rhelp02a/archive/12518.html). SpG This produces a matrix of dimension (r, choose(k, r)). SpG With this matrix, you could convert an integer between SpG 1 and choose(k, r) into an r-vector by table look-up. SpG Reading the code for nchoosek might help you further if SpG this does not seem appropriate for you. Note that *if* the above is the answer, I'd rather recommend to use combn() from package combinat, since a (slightly improved) version of combn() has been part of R-devel (to become R 2.4.0 in October) for a while. Also, combn() from combinat precedes nchoosek() historically and is also faster. Martin Maechler, ETH Zurich SpG I found this just now using 'RSiteSearch(all SpG subsets of a size)', which produced 102 hits. Another SpG one that looked like it might help you is SpG http://finzi.psych.upenn.edu/R/Rhelp02a/archive/1717.html;. SpG Hope this helps, Spencer Graves SpG Kjetil Brinchmann Halvorsen wrote: Hola! I am programming a class (S3) symarray for storing the results of functions symmetric in its k arguments. Intended use is for association indices for more than two variables, for instance coresistivity against antibiotics. There is one programming problem I haven't solved, making an inverse of the index function indx() --- se code below. It could for instance return the original k indexes in strictly increasing order, to make indx() formally invertible. Any ideas? Kjetil Code: # Implementing an S3 class for symarrays with array rank r for dimension # [k, k, ..., k] with k=r repeated r times. We do not store the diagonal. # Storage requirement is given by {r, k}= choose(k, r) # where r=array rank, k=maximum index symarray - function(data=NA, dims=c(1,1)){ r - dims[1] k - dims[2] if(r k) stop(symarray needs dimension larger than array rank) len - choose(k, r) out - data[1:len] attr(out, dims) - dims class(out) - symarray out } # Index calculation: indx - function(inds, k){ r - length(inds) if(r==1) return(inds) else { if(inds[1]==1) { return( indx(inds[-1]-1, k-1 ) ) } else { return( indx(c(inds[1]-1, seq(from=k-r+2, by=1, to=k)), k) + indx( inds[-1]-inds[1], k-inds[1] )) } } } # end indx # Methods for assignment and indexing: [.symarray - function(x, inds, drop=FALSE){ dims - attr(x, dims) k - dims[2] inds - indx(inds, k) res - NextMethod([, x) res } [-.symarray - function(x, inds, value){ dims - attr(x, dims) k - dims[2] inds - indx(inds, k) res - NextMethod([-, x) res } __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html SpG __ SpG R-help@stat.math.ethz.ch mailing list SpG https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do SpG read the posting guide! SpG http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Question regarding reading arrayvision files in limma
Dear Debayan, I only use Genepix files but I will try to help you with the experience I have concerning layout troubleshooting and limma. First : RG$R, RG$Rb, RG$G and RG$Gb seem strange to me. As I understand these data, it has to be a vector of the same length than the number of spots. It seems to me like if it was the layout of your chip, and as 3*6*10*14=2520 which is very different from 7195 that you have in RG$genes, maybe the problem is here (that there is more spots that announced by the layout ? Second : I had problems when there is blocks with incomplete line(s). It seems to me (but I did not push so far so maybe Gordon Smyth will correct me) that the function getLayout() assumes that the block on the chip are entire, so if it is not the case it is waiting for more spots that it effectively finds in the data. Finally, what do you have in RG$ngrid.c, RG$ngrid.r etc .. ? HTH Florence. On 5/25/06, Debayan Datta [EMAIL PROTECTED] wrote: Hi Everyone, I have been trying to read some Arrayvision files( 2 channel cDNA) and am having some problem. My code is : setwd('C:/work/data/limma/ndd1'); files - c('ndd1_1.txt','ndd1_2.txt','ndd1_3.txt'); RG=read.maimages(files,arrayvision,sep=\t); #Normalisation MA=normalizeWithinArrays(RG); #plotPrintTipLoess(MA); #Fit Linear model and Empirical Bayes method fit=lmFit(MA); fit - eBayes(fit) ; res=topTable(fit,sort.by=P,number=7200); hist(res$P.Value,breaks=50); The error message I get is : Read ndd1_1.txt Read ndd1_2.txt Read ndd1_3.txt Error in switch(method, loess = { : Layout argument not specified I then created a layout file and added in the following code, but I still got the same error message. RG$genes=readGAL(); gal=readGAL(); layout=getLayout(gal); Also, when I type RG , I get, RG An object of class RGList $R [1] 3 $Rb [1] 6 $G [1] 10 $Gb [1] 14 $targets FileName ndd1_1 ndd1_1.txt ndd1_2 ndd1_2.txt ndd1_3 ndd1_3.txt $source [1] arrayvision $genes Block Row ColumnID Name 1 1 1 1 my_ID my_Description 2 1 1 2 my_ID my_Description 3 1 1 3 my_ID my_Description 4 1 1 4 my_ID my_Description 5 1 1 5 my_ID my_Description 7195 more rows ... Thus, its not reading the data in the 2 forground and background channels, just the column indices. Could you please suggest what I am doing wrong? Thank you for your time! Debayan __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] query: lme
Good R-users, I have difficulties accessing the variance components for an lme fit when the variance covariance matrix of the random effects is not positive definite. For example, i fit the following model: ggg - lme (ST~ -1 + as.factor(endp):Z.sas + as.factor(endp), data=dat2a, random=~-1 + as.factor(endp) + as.factor(endp):Z.sas|as.factor(trials), correlation = corSymm(form=~1|as.factor(trials)/as.factor(id)), weights=varIdent(form=~1|endp)) intervals(ggg, which=var-cov) when i try to access the variance components using the 'intervals' function i get the following error message: Error in intervals.lme(ggg, which = var-cov) : Cannot get confidence intervals on var-cov components: Non-positive definite approximate variance-covariance Is there a way out of this? or better still Is there another function through which i can access these variance components other than the intervals function? Kind regards Pryseley - Sneak preview the all-new Yahoo.com. It's not radically different. Just radically better. [[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] Load a program from internet
Hello, I know how to load a program from my computer with the command source. example : source(home/david/prog.r) But is it possible to load a program from the internet ? example : source(http://davidblog/prog.r;) -- David [[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] variation on vioplot?
I have been using the package vioplot to make boxplot/densityplots. Now I am looking for a variation of this, and I was wondering if someone could give me any tips on how to do what I want. I have nine groups of values, each group containing two sets of numbers which show the difference between two values, such as this: Start dev Stop dev -1 10 5 -2 0 -4 and so on. This has uptil now been shown as vioplots with start deviation above the other, and I have had three groups within each vioplot. I would now like to have these two next to each other, each with their own zero point. Does anybody have any tips on how to achieve this? Any help much appreciated. Karin -- Karin Lagesen, PhD student [EMAIL PROTECTED] http://www.cmbn.no/rognes/ __ R-help@stat.math.ethz.ch 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] Load a program from internet
David Hajage [EMAIL PROTECTED] writes: Hello, I know how to load a program from my computer with the command source. example : source(home/david/prog.r) But is it possible to load a program from the internet ? example : source(http://davidblog/prog.r;) Yes ( e.g., Bioconductor uses this mechanism, see http://www.bioconductor.org/docs/install-howto.html, and - not to put too fine a point on it: Did you consider reading the part of ?source that says file: a connection or a character string giving the pathname of the file or URL to read from. ) -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] troubles with kzft
I need to compute the Fourier transformation of the autocovariance function of a panel composed by 196 time series. With autocovarance function i mean a matrix Gamma(u) which contains the covariance of all the variables at time (t) with all the variables evaluated at time (t-u). It can be represented by a tri-dimensional matrix whose dimensions are (196,196,K), where K is the maximum lag length I am considering (in my case 18). Once computed the Fourier transformation I need to compute its eigenvalues at any frequencies. My R-scripts do not work, and I am trying to use the function 'kzft' iin the package 'kzft'. The problem is that i do not understand well how to use it. In particular, i tryed with this script (though the meaning of some of the coefficients is non clear to me): aa-coeff.kzft(100,1) frq-seq(0, 2*pi, by=0.05) tf1-transfun.kzft(100,1,frq,0.0025) X-kzft(Gamma[(1:196),(1:196),], 100, 1, 1, 0.005) thanks in advance. Lorenzo -- credo nella ragione umana, e nella liberta' e nella giustizia che dalla ragione scaturiscono. (sciascia) [[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] Load a program from internet
ok... I've juste understood that my computer (at work) doesn't let R access to internet... That's why I was wondering if it works ! I've made a function which adds exposed to risk on your survival curve... If you want to try : source(http://pepereblog.free.fr/blog/share/tigeR%202.0.r;) tigeR.graph(aml$time, aml$status, aml$x, tps = seq(0,72, 24), IC = 24, decalIC = 0.5, legend = names(table(aml$x)), main = Survival) tigeR.results(aml$time, aml$status, aml$x, latex = F) It's in french, but you can modify it if you want... Thank you all for your answers, and sorry for my english. 2006/5/29, David Hajage [EMAIL PROTECTED]: and sorry for my english !... 2006/5/29, David Hajage [EMAIL PROTECTED]: thank you. Same result... It doesn't work ! But i will try at home... I'm sure il will work ! Or you can try : I've made a function which adds exposed to risk on your survival curve... source(http://pepereblog.free.fr/blog/share/tigeR%202.0.r ) tigeR.graph(aml$time, aml$status, aml$x, tps = seq(0,72, 24), IC = 24, decalIC = 0.5, legend = names(table(aml$x)), main = Survival) tigeR.results(aml$time, aml$status, aml$x, latex = F) It's in french, but you can modify it if you want... Thank for your response. 2006/5/29, Christian Ritz [EMAIL PROTECTED]: Hi David, try: source(url( http://davidblog/prog.r )) Christian -- David -- David -- David [[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] mtext in trellis-loop colorkey
Hello, I would like to create a sequence of plots (using a for loop). I read in the FAQ that print() has to be used in order to obtain any output. This works perfectly fine as long as I only consider one function call in the loop, but I would like to add mtext() to the each plot in the loop. Unfortunately, this did not work. Any suggestions? As you can see from the provided example, there is another problem with such animations: You do not see any difference in the plots (because only the height changes). Is there any possibility to keep a fixed scale for the colorkey (fixed labels and also fixed colors) and to see the different heights of the function from the colors (e.g. the first plot should be mainly gray (as it is the lowest), the last one mainly white (as it is the highest)). Thanks very much! marius Here is a complete minimal example: remove(list=objects()) library(lattice) for(i in 1:4){ output_file_path-paste(~/Desktop/test_,i,.png,sep=) x-rep(seq(-3,3,length=50),50) y-rep(seq(-3,3,length=50),each=50) z-x*y+10*i trellis.device(png,color=F,file=output_file_path) print(wireframe (z~x*y,drape=T,distance=0,zoom=0.84,cuts=100,col.regions=gray (100:400/400),colorkey=list(tick.number=6))) #print(mtext(paste(Parameter= ,1,sep=),side=3,line=0)) #This does not work! dev.off() } __ R-help@stat.math.ethz.ch 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] Analog to matlab repmat function
Thanks a lot to all of you! Now I see, that use of R (and S) is very different from use of Matlab. There is no one-to-one correspondence. I'm shy to ask so stupid questions, but name of language R makes it hard to find relevant links trough google. Thanks for pointing to S Poetry, I would never find it alone, since I was looking for R-related info. Thanks, Daniil. On 5/29/06, Patrick Burns [EMAIL PROTECTED] wrote: S Poetry may be of use to you. Patrick Burns [EMAIL PROTECTED] +44 (0)20 8525 0696 http://www.burns-stat.com (home of S Poetry and A Guide for the Unwilling S User) Daniil Ivanov wrote: Hello, I'm trying to switch from Matlab to R-project, and having some difficulties. I make a use of multidimensional matrices. For example, I need to extract mean from one of the dimensions: % we have matrix data of size: 130 x 11 x 350 x 2 data = data - repmat(mean(data,3),[130 1 1 1]); In R project I managed to do that in a very pervarsive way: # mean(data,3) in R base - apply(data,c(2,3,4),mean) # repmat(...,[130 1 1 1]) base - rep(base,130) dim(base) - c(11,2,350,130) base - aperm(base,c(4,1,2,3)) # data = data - repmat data - data - base; Could you please show me a more elegant way of doing the same in R :) Thanks, Daniil. __ R-help@stat.math.ethz.ch 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] Analog to matlab repmat function
Actually google does understand R. For example, 1. google for: R and you will get the R home page.as first hit 2. google for R repmat and you will find Robin Hankin's lexicon that translates between R and octave/maxima which includes a repmat function as 7th hit. On 5/29/06, Daniil Ivanov [EMAIL PROTECTED] wrote: Thanks a lot to all of you! Now I see, that use of R (and S) is very different from use of Matlab. There is no one-to-one correspondence. I'm shy to ask so stupid questions, but name of language R makes it hard to find relevant links trough google. Thanks for pointing to S Poetry, I would never find it alone, since I was looking for R-related info. Thanks, Daniil. On 5/29/06, Patrick Burns [EMAIL PROTECTED] wrote: S Poetry may be of use to you. Patrick Burns [EMAIL PROTECTED] +44 (0)20 8525 0696 http://www.burns-stat.com (home of S Poetry and A Guide for the Unwilling S User) Daniil Ivanov wrote: Hello, I'm trying to switch from Matlab to R-project, and having some difficulties. I make a use of multidimensional matrices. For example, I need to extract mean from one of the dimensions: % we have matrix data of size: 130 x 11 x 350 x 2 data = data - repmat(mean(data,3),[130 1 1 1]); In R project I managed to do that in a very pervarsive way: # mean(data,3) in R base - apply(data,c(2,3,4),mean) # repmat(...,[130 1 1 1]) base - rep(base,130) dim(base) - c(11,2,350,130) base - aperm(base,c(4,1,2,3)) # data = data - repmat data - data - base; Could you please show me a more elegant way of doing the same in R :) Thanks, Daniil. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] newbie question: ROW average
Hi, I am new to R programming. I have a 992 x 74 matrix. I would like to form a new matrix by averging each 4 rows from the original one. How can I use 'apply' instead of usual mean inside the nested for loop? Thanks in advance. regards, ezhil __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] newbie question: ROW average
look at ?rowMeans; you can also use apply(mat, 1, mean) but rowMeans() is better. 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: A Ezhil [EMAIL PROTECTED] To: r-help r-help@stat.math.ethz.ch Sent: Monday, May 29, 2006 1:24 PM Subject: [R] newbie question: ROW average Hi, I am new to R programming. I have a 992 x 74 matrix. I would like to form a new matrix by averging each 4 rows from the original one. How can I use 'apply' instead of usual mean inside the nested for loop? Thanks in advance. regards, ezhil __ R-help@stat.math.ethz.ch 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
Re: [R] mtext in trellis-loop colorkey
mtext() is not a lattice function. Could you not use the 'main' argument? Peter Ehlers Marius Hofert wrote: Hello, I would like to create a sequence of plots (using a for loop). I read in the FAQ that print() has to be used in order to obtain any output. This works perfectly fine as long as I only consider one function call in the loop, but I would like to add mtext() to the each plot in the loop. Unfortunately, this did not work. Any suggestions? As you can see from the provided example, there is another problem with such animations: You do not see any difference in the plots (because only the height changes). Is there any possibility to keep a fixed scale for the colorkey (fixed labels and also fixed colors) and to see the different heights of the function from the colors (e.g. the first plot should be mainly gray (as it is the lowest), the last one mainly white (as it is the highest)). Thanks very much! marius Here is a complete minimal example: remove(list=objects()) library(lattice) for(i in 1:4){ output_file_path-paste(~/Desktop/test_,i,.png,sep=) x-rep(seq(-3,3,length=50),50) y-rep(seq(-3,3,length=50),each=50) z-x*y+10*i trellis.device(png,color=F,file=output_file_path) print(wireframe (z~x*y,drape=T,distance=0,zoom=0.84,cuts=100,col.regions=gray (100:400/400),colorkey=list(tick.number=6))) #print(mtext(paste(Parameter= ,1,sep=),side=3,line=0)) #This does not work! dev.off() } __ R-help@stat.math.ethz.ch 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] newbie question: ROW average
A Ezhil [EMAIL PROTECTED] writes: Hi, I am new to R programming. I have a 992 x 74 matrix. I would like to form a new matrix by averging each 4 rows from the original one. How can I use 'apply' instead of usual mean inside the nested for loop? How about dim(M) - c(4,248,74) mn - apply(M, c(2,3), mean) (Yes, this is a bit confusing and easy to get wrong...) -- 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] query: lme
I'm not exactly certain, but it seems to me that you're including a factor in the LHS of the random term. You might write down the model that you're fitting, and reflect upon it. Of course, there may be a perfectly good reason, but that does make for an exceptionally complicated model. You might try to reinvent your question so that it is matched by a simpler model. Of course, I'm only speculating, because you haven't told us anything about the data, or your intentions. In general, you might find that the variance covariance matrix of the random effects not being positive definite could be a sign of a model that is not a good match with the available data. Cheers Andrew On Mon, May 29, 2006 at 02:15:54AM -0700, Pryseley Assam wrote: Good R-users, I have difficulties accessing the variance components for an lme fit when the variance covariance matrix of the random effects is not positive definite. For example, i fit the following model: ggg - lme (ST~ -1 + as.factor(endp):Z.sas + as.factor(endp), data=dat2a, random=~-1 + as.factor(endp) + as.factor(endp):Z.sas|as.factor(trials), correlation = corSymm(form=~1|as.factor(trials)/as.factor(id)), weights=varIdent(form=~1|endp)) intervals(ggg, which=var-cov) when i try to access the variance components using the 'intervals' function i get the following error message: Error in intervals.lme(ggg, which = var-cov) : Cannot get confidence intervals on var-cov components: Non-positive definite approximate variance-covariance Is there a way out of this? or better still Is there another function through which i can access these variance components other than the intervals function? Kind regards Pryseley - Sneak preview the all-new Yahoo.com. It's not radically different. Just radically better. [[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 __ R-help@stat.math.ethz.ch 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] mtext in trellis-loop colorkey
Marius Hofert wrote: Hello, I would like to create a sequence of plots (using a for loop). I read in the FAQ that print() has to be used in order to obtain any output. This works perfectly fine as long as I only consider one function call in the loop, but I would like to add mtext() to the each plot in the loop. Unfortunately, this did not work. Any suggestions? As you can see from the provided example, there is another problem with such animations: You do not see any difference in the plots (because only the height changes). Is there any possibility to keep a fixed scale for the colorkey (fixed labels and also fixed colors) and to see the different heights of the function from the colors (e.g. the first plot should be mainly gray (as it is the lowest), the last one mainly white (as it is the highest)). Thanks very much! marius Here is a complete minimal example: remove(list=objects()) library(lattice) for(i in 1:4){ output_file_path-paste(~/Desktop/test_,i,.png,sep=) x-rep(seq(-3,3,length=50),50) y-rep(seq(-3,3,length=50),each=50) z-x*y+10*i trellis.device(png,color=F,file=output_file_path) print(wireframe (z~x*y,drape=T,distance=0,zoom=0.84,cuts=100,col.regions=gray (100:400/400),colorkey=list(tick.number=6))) #print(mtext(paste(Parameter= ,1,sep=),side=3,line=0)) #This does not work! dev.off() } I think the page argument will help you the title (or simply use main; see ?xyplot). For the colorkey, set the zlim and at arguments. library(lattice) for(i in 1:4) { output_file_path - paste(./test_, i, .png, sep = ) x - rep(seq(-3, 3, length = 50), 50) y - rep(seq(-3, 3, length = 50), each = 50) z - x * y + 10 * i trellis.device(png, color = FALSE, file = output_file_path) w - wireframe(z ~ x * y, drape = TRUE, distance = 0, zoom = 0.84, cuts = 100, col.regions = gray(100:400/400), page = function(n) { label - paste(Parameter = , i, sep = ) ltext(0.5, 1, label) }, zlim = c(0, 50), at = seq(0, 50, length = 10), colorkey = list(tick.number = 6)) print(w) dev.off() } HTH, --sundar __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] newbie question: ROW average
Peter Dalgaard wrote: How about dim(M) - c(4,248,74) mn - apply(M, c(2,3), mean) Hey! That's sexy! Much better than my kludgy suggestion! cheers, Rolf __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] newbie question: ROW average
Dimitris Rizopoulos wrote: look at ?rowMeans; you can also use apply(mat, 1, mean) but rowMeans() is better. By my reading of the question, this is not what Ezhil wants. He said: ``I have a 992 x 74 matrix. I would like to form a new matrix by averaging each 4 rows from the original one.'' I.e. he wants (I think) the first row of the new matrix to be the mean of the first 4 rows of the old one, the second row of the new matrix to be the mean of rows 5 through 8 of the old one, and so on. One way this could be done is via m.new - t(apply(array(t(m.old),dim=c(74,4,992/4)),c(1,3),mean)) cheers, Rolf Turner [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] newbie question: ROW average
yes you're right; it was my mistake. 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://med.kuleuven.be/biostat/ http://www.student.kuleuven.be/~m0390867/dimitris.htm - Original Message - From: Rolf Turner [EMAIL PROTECTED] To: [EMAIL PROTECTED]; [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Sent: Monday, May 29, 2006 1:55 PM Subject: Re: [R] newbie question: ROW average Dimitris Rizopoulos wrote: look at ?rowMeans; you can also use apply(mat, 1, mean) but rowMeans() is better. By my reading of the question, this is not what Ezhil wants. He said: ``I have a 992 x 74 matrix. I would like to form a new matrix by averaging each 4 rows from the original one.'' I.e. he wants (I think) the first row of the new matrix to be the mean of the first 4 rows of the old one, the second row of the new matrix to be the mean of rows 5 through 8 of the old one, and so on. One way this could be done is via m.new - t(apply(array(t(m.old),dim=c(74,4,992/4)),c(1,3),mean)) cheers, Rolf Turner [EMAIL PROTECTED] Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] newbie question: ROW average
Try this: rowsum(mat, gl(nrow(mat)/4, 4)) / 4 On 5/29/06, A Ezhil [EMAIL PROTECTED] wrote: Hi, I am new to R programming. I have a 992 x 74 matrix. I would like to form a new matrix by averging each 4 rows from the original one. How can I use 'apply' instead of usual mean inside the nested for loop? Thanks in advance. regards, ezhil __ R-help@stat.math.ethz.ch 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] query: lme
As noted below, it is hard to diagnose the specific issue. But, as a recommendation, I would suggest using lmer and then investigating the parameters of the model using the MCMCsamp() function. You can then do all diagnostics using the various functions in the coda package as MCMCsamp() returns an object of class mcmc. -Original Message- From: [EMAIL PROTECTED] on behalf of Andrew Robinson Sent: Mon 5/29/2006 7:22 AM To: Pryseley Assam Cc: R-Help Discussion Subject:Re: [R] query: lme I'm not exactly certain, but it seems to me that you're including a factor in the LHS of the random term. You might write down the model that you're fitting, and reflect upon it. Of course, there may be a perfectly good reason, but that does make for an exceptionally complicated model. You might try to reinvent your question so that it is matched by a simpler model. Of course, I'm only speculating, because you haven't told us anything about the data, or your intentions. In general, you might find that the variance covariance matrix of the random effects not being positive definite could be a sign of a model that is not a good match with the available data. Cheers Andrew On Mon, May 29, 2006 at 02:15:54AM -0700, Pryseley Assam wrote: Good R-users, I have difficulties accessing the variance components for an lme fit when the variance covariance matrix of the random effects is not positive definite. For example, i fit the following model: ggg - lme (ST~ -1 + as.factor(endp):Z.sas + as.factor(endp), data=dat2a, random=~-1 + as.factor(endp) + as.factor(endp):Z.sas|as.factor(trials), correlation = corSymm(form=~1|as.factor(trials)/as.factor(id)), weights=varIdent(form=~1|endp)) intervals(ggg, which=var-cov) when i try to access the variance components using the 'intervals' function i get the following error message: Error in intervals.lme(ggg, which = var-cov) : Cannot get confidence intervals on var-cov components: Non-positive definite approximate variance-covariance Is there a way out of this? or better still Is there another function through which i can access these variance components other than the intervals function? Kind regards Pryseley - Sneak preview the all-new Yahoo.com. It's not radically different. Just radically better. [[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 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] TsayData
Hi, I'm trying to work with TsayData in fSeries package. How can i fetch any time series data of this package. Please advice. Thanks, Sumanta Basak. Send instant messages to your online friends http://in.messenger.yahoo.com __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] 29/05/06
LUCKYDAY LOTTO 1102 VAN DER TRAAT AMSTERDAM THE NETHERLANDS www.luckyday.nl Dear Winner, We are please to announce you as one of the 10 lucky winners in the luckyday Lotto draw held yesterday the 2nd of May 2006.All 10 winning Addresses were randomly selected from a batch of 50,000,000 international Emails. Your email address emerged alongside 9 others as category 2 Winners in this year's Annual luckyday Lotto Draw. Consequently, you have therefore been approved for a total pay out of (one million EUROS) only. The following particulars Are attached to your lotto payment order: (i) Winning numbers: 37-13-43-85-67-11 (ii) Email ticket number: FL754/22/76 (iii) Lotto code number: FL23112 (iv) The file Ref number: FL/04/5237900/NL Please contact the under listed claims officer via email and/or phone as soon as possible for the immediate release of your winnings and also send him the information requested below: Mr.Pedro Lopez Complete Solutions Ltd Tel-0031 644751015 Email- [EMAIL PROTECTED] 1. IN FULL: 2. ADDRESS:. 3. NATIONALITY:. 4. AGE:. 5. OCCUPATION: 6. PHONE/FAX:... 7.COUNTRY OF RESIDENCE. N.B: Please quote your Reference number in all correspondence with the claims officer. Winners are advice to keep this award confidential until prize are claimed to avert incidence of impersonation by unscrupulous elements. Once again on behalf of all our staff, CONGRATULATIONS!!! Sincerely, Paul Green Promotions Manager The Luckyday Lotto Company The Luckyday Lotto Awards is proudly sponsored by the Microsoft Corporation, the Intel Group, Toshiba, Dell Computers,Mackintosh and a Conglomeration of other international IT companies. The Luckyday lotto Internet Draw is held once in a year and is so organized to encourage the use of The Internet and computers worldwide. We are proud to say that over 200 Million Euros are won annually in more than 150 countries worldwide. We also encourage you to visit our website at www.luckyday.nl and Take your chance to play and become part of our annual winners. Millions are won on a yearly basis on one of the world's most famous lotto Service. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] newbie question: ROW average
Hi, Thank you all (Dimitris, Peter, Rolf, Gabor) for your suggestions. I tried with all your suggestions. I am getting different answers when I use: rowsum(mat, gl(nrow(mat)/4, 4)) / 4 and m.new - t(apply(array(t(m.old),dim=c(74,4,992/4)),c(1,3),mean)) When I tried with (assuming 'M' is my old matrix): dim(M) - c(4,248,74) mn - apply(M, c(2,3), mean) the following error occured: Error: dim- : dims [product 73408] do not match the length of object [74] When I manually checked the answers, it seems that rowsum(mat, gl(nrow(mat)/4, 4)) / 4 gives me the correct answer. Thanks again for your time suggestions. Regards, Ezhil --- Dimitris Rizopoulos [EMAIL PROTECTED] wrote: yes you're right; it was my mistake. 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://med.kuleuven.be/biostat/ http://www.student.kuleuven.be/~m0390867/dimitris.htm - Original Message - From: Rolf Turner [EMAIL PROTECTED] To: [EMAIL PROTECTED]; [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch Sent: Monday, May 29, 2006 1:55 PM Subject: Re: [R] newbie question: ROW average Dimitris Rizopoulos wrote: look at ?rowMeans; you can also use apply(mat, 1, mean) but rowMeans() is better. By my reading of the question, this is not what Ezhil wants. He said: ``I have a 992 x 74 matrix. I would like to form a new matrix by averaging each 4 rows from the original one.'' I.e. he wants (I think) the first row of the new matrix to be the mean of the first 4 rows of the old one, the second row of the new matrix to be the mean of rows 5 through 8 of the old one, and so on. One way this could be done is via m.new - t(apply(array(t(m.old),dim=c(74,4,992/4)),c(1,3),mean)) cheers, Rolf Turner [EMAIL PROTECTED] Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] newbie question: ROW average
A Ezhil [EMAIL PROTECTED] writes: When I tried with (assuming 'M' is my old matrix): dim(M) - c(4,248,74) mn - apply(M, c(2,3), mean) the following error occured: Error: dim- : dims [product 73408] do not match the length of object [74] In that case, M clearly wasn't a 992x74 matrix! -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] parameter-restrictions in OPTIM
Hallo, I'm trying to optimize a function A which calls another function B internally. This function A has got two input parameters, say a1 and a2. Using OPTIM there ist no problem restricting these parameters: 0 = a1 = 1 and 0 = a2 = 1. My problem is an additional restriction 0 = a1+a2 = 1. I'd be grateful if you could show me an easy way to find a solution using OPTIM? If not, which other function is applicable? Thanks for your help, Florian Staab _ Dipl.-Kfm. Florian Staab Universität Osnabrück Fachbereich Wirtschaftswissenschaften Statistik / Empirische Wirtschaftsforschung Rolandstr. 8 (Raum 208) D-49069 Osnabrück Tel: 0541-969-2756 Fax: 0541-969-12756 Email: [EMAIL PROTECTED] Homepage: http://nts4.oec.uni-osnabrueck.de/stat1/stat1.htm _ [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] OT: Monograph on Statistical Programming
Dear all, my question might be a bit off-topic. Is there anything like a standard textbook on statistical programming? With that I don't mean anything like MASS, S Programming, Programming with Data, ... (no offense meant, they are fantastic books and each of those three helped me a great deal). Rather in the direction of the Numerical Recipes Series addressing how to implement various functions. There are bits and pieces all over the place (The Art of Computer Programming, Numerical Recipes, Algorithms in ...(Sedgewick), etc.)... But is there one book like The Art of Statistical Programming? Maybe writing such a book could be impossible because the choice of implementation language could be crucial and also the choice of algorithms could be problematic. Nevertheless I hope that someone dared to tackle such a project. Thanks, Roland -- This mail has been sent through the MPI for Demographic Rese...{{dropped}} __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] mtext in trellis-loop colorkey
On 5/29/06, Marius Hofert [EMAIL PROTECTED] wrote: Hello, I would like to create a sequence of plots (using a for loop). I read in the FAQ that print() has to be used in order to obtain any output. This works perfectly fine as long as I only consider one function call in the loop, but I would like to add mtext() to the each plot in the loop. Unfortunately, this did not work. Any suggestions? Others have already responded to this. As you can see from the provided example, there is another problem with such animations: You do not see any difference in the plots (because only the height changes). Is there any possibility to keep a fixed scale for the colorkey (fixed labels and also fixed colors) and to see the different heights of the function from the colors (e.g. the first plot should be mainly gray (as it is the lowest), the last one mainly white (as it is the highest)). It seems to me that you are missing the primary point of Trellis graphics, which is not having to manually manage such details in a multipanel plot. Consider the following (which should produce 4 files with a png device): library(lattice) x-rep(seq(-3,3,length=50),50) y-rep(seq(-3,3,length=50),each=50) z1 - x * y + 10 * 1 z2 - x * y + 10 * 2 z3 - x * y + 10 * 3 z4 - x * y + 10 * 4 wireframe(z1 + z2 + z3 + z4 ~ x * y, outer = TRUE, drape = TRUE, zlab = z, layout = c(1, 1), distance=0,zoom=0.84, cuts=100, col.regions=gray(100:400/400), colorkey=list(tick.number=6)) Caveat: this won't work if you want the z-axis completely filled up in each panel. -Deepayan __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] parameter-restrictions in OPTIM
You could reparametrize or use constrOptim, Hth, Ingmar From: Florian Staab [EMAIL PROTECTED] Date: Mon, 29 May 2006 12:08:59 +0200 To: r-help@stat.math.ethz.ch Subject: [R] parameter-restrictions in OPTIM Hallo, I'm trying to optimize a function A which calls another function B internally. This function A has got two input parameters, say a1 and a2. Using OPTIM there ist no problem restricting these parameters: 0 = a1 = 1 and 0 = a2 = 1. My problem is an additional restriction 0 = a1+a2 = 1. I'd be grateful if you could show me an easy way to find a solution using OPTIM? If not, which other function is applicable? Thanks for your help, Florian Staab _ Dipl.-Kfm. Florian Staab Universität Osnabrück Fachbereich Wirtschaftswissenschaften Statistik / Empirische Wirtschaftsforschung Rolandstr. 8 (Raum 208) D-49069 Osnabrück Tel: 0541-969-2756 Fax: 0541-969-12756 Email: [EMAIL PROTECTED] Homepage: http://nts4.oec.uni-osnabrueck.de/stat1/stat1.htm _ [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] how to multiply a constant to a matrix?
comtech == comtech usa [EMAIL PROTECTED] on Tue, 23 May 2006 23:14:21 -0700 writes: comtech imagine when you have complicated matrix algebra computation using R, comtech you cannot prevent some middle-terms become quadratic and absorbs into one comtech scalar, right? comtech if R cannot intelligently determine this, and you have to manually add comtech drop everywhere, comtech do you think it is reasonable? I don't think you can find any reasonable example where `` you have to manually add drop everywhere ''. To the contrary: In my experience, the S behavior (which won't be changed -- forget about that if you still consider it!!), is even rather useful, pointing to ``infelicities'' in your code most of the time, when matrix dimensions do not match. Further note, that c( 1x1--matrix ) also produces a scalar and gives slightly more readable formulae --- however in the present case, c() is only advisable if you *know* you don't have a ``proper'' (n x k, with n, k 1) matrix. Martin Maechler, ETH Zurich comtech On 5/23/06, Patrick Burns [EMAIL PROTECTED] wrote: I think drop(B/D) * solve(A) would be a more transparent approach. It isn't that R can not do what you want, it is that it is saving you from shooting yourself in the foot in your attempt. What you are doing is not really a matrix computation. Patrick Burns [EMAIL PROTECTED] +44 (0)20 8525 0696 http://www.burns-stat.com (home of S Poetry and A Guide for the Unwilling S User) Michael wrote: This is very strange: I want compute the following in R: g = B/D * solve(A) where B and D are quadratics so they are just a scalar number, e.g. B=t(a) %*% F %*% a; I want to multiply B/D to A^(-1), but R just does not allow me to do that and it keeps complaining that nonconformable array, etc. I tried the following two tricks and they worked: as.numeric(B/D) * solve(A) diag(as.numeric(B/D), 5, 5) %*% solve (A) But if R cannot intelligently do scalar and matrix multiplication, it is really problemetic. It basically cannot be used to do computations, since in complicated matrix algebras, you have to distinguish where is scalar, and scalars obtained from quadratics cannot be directly used to multiply another matrix, etc. It is going to a huge mess... Any thoughts? [[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 comtech [[alternative HTML version deleted]] comtech __ comtech R-help@stat.math.ethz.ch mailing list comtech https://stat.ethz.ch/mailman/listinfo/r-help comtech PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Is there a way to draw 3d plot?
Hi, I am not aware what this function does in matlab, but I was wondering why cloud or wireframe in the lattice package has not entered the discussion. Sorry if I am being naive. Ritwik. http://darwin.cwru.edu/~rsinha [[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] parameter-restrictions in OPTIM
One trick is to add a penalty to your function for violations of the linear constraint. However, if the constraint is binding, then gradient methods may not perform well since the objective won't be differentiable there. If the constraint looks to be binding, then you can reparameterize to a single parameter. Patrick Burns [EMAIL PROTECTED] +44 (0)20 8525 0696 http://www.burns-stat.com (home of S Poetry and A Guide for the Unwilling S User) Florian Staab wrote: Hallo, I'm trying to optimize a function A which calls another function B internally. This function A has got two input parameters, say a1 and a2. Using OPTIM there ist no problem restricting these parameters: 0 = a1 = 1 and 0 = a2 = 1. My problem is an additional restriction 0 = a1+a2 = 1. I'd be grateful if you could show me an easy way to find a solution using OPTIM? If not, which other function is applicable? Thanks for your help, Florian Staab _ Dipl.-Kfm. Florian Staab Universität Osnabrück Fachbereich Wirtschaftswissenschaften Statistik / Empirische Wirtschaftsforschung Rolandstr. 8 (Raum 208) D-49069 Osnabrück Tel: 0541-969-2756 Fax: 0541-969-12756 Email: [EMAIL PROTECTED] Homepage: http://nts4.oec.uni-osnabrueck.de/stat1/stat1.htm _ [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] lme, best model without convergence
I'm sorry. My post was in error on that point. I did rm(fm1.), then tried again and got the same error message you got. Conclusion: Either the 'returnObject' argument doesn't work in the version(s) of lme that you and I are using or we've flunked another literacy exam while RTFM. With luck, someone else will read this post and answer this question for both of us. If I had time to work on this, I'd make local copies of the code and walk through them line by line using 'debug' until I figured out how to solve the problem. I've outlined how to do this in previous posts, and RSiteSearch(graves debug lme) led me to http://finzi.psych.upenn.edu/R/Rhelp02a/archive/72784.html;, in case you want to try that. Hope this helps, Spencer Graves sessionInfo() Version 2.3.0 (2006-04-24) i386-pc-mingw32 attached base packages: [1] methods stats graphics grDevices utils datasets [7] base other attached packages: nlme 3.1-73 # Dear Spencer Graves, thank you for your reply and the hint to the 'control' argument of the lme function. Thank you for your self-contained example. Running this example I get the verbose output. So I am able to trace the iterations and maybe I am able to detect ridges. However, I just updated the nlme-package at the Mainz CRAN mirror (on R-windows-2.2.1 version), but I do not get the result object from the not-converging fit: lmeCtlList - + lmeControl(maxIter=2, + msMaxIter=3, tolerance=1e-4, niter=4, + msTol=1e-5, nlmStepMax=5, + ,msVerbose=TRUE + ,returnObject=TRUE + ) fm1. - lme(distance ~ age, data = Orthodont, + control=lmeCtlList) # random is ~ age 0 318.310: -0.202572 0.0501133 2.17193 1 318.297: -0.202070 0.0662548 2.16740 2 318.286: -0.207284 0.0593806 2.15302 3 318.257: -0.215916 0.0783917 2.12677 3 318.257: -0.215916 0.0783917 2.12677 Error in lme.formula(distance ~ age, data = Orthodont, control = lmeCtlList) : iteration limit reached without convergence (9) fm1. Error: object fm1. not found ## lmeControl package:nlme R Documentation Control Values for lme Fit Description: The values supplied in the function call replace the defaults and a list with all possible arguments is returned. The returned list is used as the 'control' argument to the 'lme' function. ## Note the last sentence. Now consider the following example, combining modifications of your example with the example from ?lme: lmeCtlList - + lmeControl(maxIter=2, +msMaxIter=3, tolerance=1e-4, niter=4, + msTol=1e-5, nlmStepMax=5, + ,msVerbose=TRUE + ,returnObject=TRUE + ) fm1. - lme(distance ~ age, data = Orthodont, + control=lmeCtlList) # random is ~ age 0 318.310: -0.202572 0.0501133 2.17193 1 318.297: -0.202071 0.0662548 2.16740 2 318.286: -0.207284 0.0593805 2.15302 3 318.257: -0.215916 0.0783916 2.12677 3 318.257: -0.215916 0.0783916 2.12677 Error in lme.formula(distance ~ age, data = Orthodont, control = lmeCtlList) : iteration limit reached without convergence (9) fm1. Linear mixed-effects model fit by REML Data: Orthodont Log-restricted-likelihood: -221.3183 Fixed: distance ~ age (Intercept) age 16.761 0.6601852 Random effects: Formula: ~age | Subject Structure: General positive-definite StdDevCorr (Intercept) 2.3270360 (Intr) age 0.2264279 -0.609 Residual1.3100396 Number of Observations: 108 Number of Groups: 27 Hope this helps. Spencer Graves p.s. First thanks for this example. I didn't know until I read your question that one could actually get lme to return something when it didn't converge. Second, you might have gotten a quicker reply if you had included a simple, self-contained example, as suggested in the posting guide! www.R-project.org/posting-guide.html. Instead, since I didn't know the answer, I had to invent one. If your example had been self contained, you might have gotten an earlier reply from someone who didn't have the time or inclination to invent an example like I did but who might have otherwise been able to solve the problem. Thomas Wutzler wrote: Dear R-help list readers, I am fitting mixed models with the lme function of the nlme package. If I get convergence depends on how the method (ML/REM) and which (and how much) parameters will depend randomly on the cluster-variable. How get the bist fit without convergence? I set the parameters msVerbose and returnObject to TRUE: lmeControl(maxIter=5, msMaxIter=200, tolerance=1e-4, niter=50, msTol=1e-5, nlmStepMax=500, ,msVerbose=TRUE ,returnObject=TRUE ) However, the
[R] line over multiplot
Dear Listmembers, I do have a multiplot with 1 row and 4 columns (4 scatter plots with same x- and y-axe plotted in one row). I would like to draw a horizontal line across the whole multiplot at a certain y-value. Unfortunately the normal abline command stops between the plot. I hope my problem is clear and there is somebody who can help me Thanks in advance Thomas H. __ R-help@stat.math.ethz.ch 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] line over multiplot
Thomas Hoffmann wrote: Dear Listmembers, I do have a multiplot with 1 row and 4 columns (4 scatter plots with same x- and y-axe plotted in one row). I would like to draw a horizontal line across the whole multiplot at a certain y-value. Unfortunately the normal abline command stops between the plot. I hope my problem is clear and there is somebody who can help me par(mfrow=c(1,4)) replicate(4, plot(1:10)) par(xpd=NA) abline(h=5) Uwe Ligges Thanks in advance Thomas H. __ R-help@stat.math.ethz.ch 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] multiple comparisons of time series data
Kyle, You might try the Wilcoxon Rank Sum test (and there is also the paired rank sum test) that may be useful. Both are found in R. There is an application of the test in the textbook by Loucks, D.P., Stedinger J.R., and Haith, D., 1981. Water Resources Systems Planning and Analysis, Prentice-Hall, Eaglewood Cliffs, New Jersey. I hope this helps… Go HOKIES! Tom Kyle Hall wrote: I am interested in a statistical comparison of multiple (5) time series' generated from modeling software (Hydrologic Simulation Program Fortran). The model output simulates daily bacteria concentration in a stream. The multiple time series' are a result of varying our representation of the stream within the model. Our main question is: Do the different methods used to represent a stream produce different results at a statistically significant level? We want to compare each otput time series to determine if there is a difference before looking into the cause within the model. In a previous study, the Kolmogorov-Smirnov k-sample test was used to compare multiple time series'. I am unsure about the strength of the Kolmogorov-Smirnov test and I have set out to determine if there are any other tests to compare multiple time series'. I know htat R has the ks.test but I am unsure how this test handles multiple comparisons. Is there something similar to a pairwise.t.test with a bonferroni corection, only with time series data? Does R currently (v 2.3.0) have a comparison test that takes into account the strong serial correlation of time series data? Kyle Hall Graduate Research Assistant Biological Systems Engineering Virginia Tech __ R-help@stat.math.ethz.ch 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 E Adams National Weather Service Ohio River Forecast Center 1901 South State Route 134 Wilmington, OH 45177 EMAIL: [EMAIL PROTECTED] VOICE: 937-383-0528 FAX:937-383-0033 __ R-help@stat.math.ethz.ch 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] useR! 2006: final program
Dear useRs, the final scientific program for useR! 2006, the second R user conference taking place in Vienna June 15-17, is now available online at http://www.R-project.org/useR-2006/ The registration is open until May 31, so you still have two days to register for this exciting event where you can meet 400 other useRs, attend 160 presentations and participate in some of the pre-conference tutorials given by prominent members of the R community. We look forward to meeting you in Vienna! The organizing team Torsten, Achim, David, Bettina, Fritz and Kurt. ___ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-announce __ R-help@stat.math.ethz.ch 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] combinatorial programming problem
On Mon, 29 May 2006, Martin Maechler wrote: SpG == Spencer Graves [EMAIL PROTECTED] on Sun, 28 May 2006 16:21:53 -0700 writes: SpG I'm not sure I understand your question, but SpG are you asking how to index choose(k, r) objects? SpG Almost 3 years ago, I asked a question like this. Andy SpG Liaw referred me to nchoosek(vsn) SpG (http://finzi.psych.upenn.edu/R/Rhelp02a/archive/12518.html). SpG This produces a matrix of dimension (r, choose(k, r)). SpG With this matrix, you could convert an integer between SpG 1 and choose(k, r) into an r-vector by table look-up. SpG Reading the code for nchoosek might help you further if SpG this does not seem appropriate for you. Note that *if* the above is the answer, I'd rather recommend to use combn() from package combinat, since a (slightly improved) version of combn() has been part of R-devel (to become R 2.4.0 in October) for a while. Also, combn() from combinat precedes nchoosek() historically and is also faster. And if the problem is large enough to bog down combn() (e.g. choose(25,12) combinations) look at this approach: http://article.gmane.org/gmane.comp.lang.r.general/62065 Martin Maechler, ETH Zurich SpG I found this just now using 'RSiteSearch(all SpG subsets of a size)', which produced 102 hits. Another SpG one that looked like it might help you is SpG http://finzi.psych.upenn.edu/R/Rhelp02a/archive/1717.html;. SpG Hope this helps, Spencer Graves SpG Kjetil Brinchmann Halvorsen wrote: Hola! I am programming a class (S3) symarray for storing the results of functions symmetric in its k arguments. Intended use is for association indices for more than two variables, for instance coresistivity against antibiotics. There is one programming problem I haven't solved, making an inverse of the index function indx() --- se code below. It could for instance return the original k indexes in strictly increasing order, to make indx() formally invertible. Any ideas? Kjetil Code: # Implementing an S3 class for symarrays with array rank r for dimension # [k, k, ..., k] with k=r repeated r times. We do not store the diagonal. # Storage requirement is given by {r, k}= choose(k, r) # where r=array rank, k=maximum index symarray - function(data=NA, dims=c(1,1)){ r - dims[1] k - dims[2] if(r k) stop(symarray needs dimension larger than array rank) len - choose(k, r) out - data[1:len] attr(out, dims) - dims class(out) - symarray out } # Index calculation: indx - function(inds, k){ r - length(inds) if(r==1) return(inds) else { if(inds[1]==1) { return( indx(inds[-1]-1, k-1 ) ) } else { return( indx(c(inds[1]-1, seq(from=k-r+2, by=1, to=k)), k) + indx( inds[-1]-inds[1], k-inds[1] )) } } } # end indx # Methods for assignment and indexing: [.symarray - function(x, inds, drop=FALSE){ dims - attr(x, dims) k - dims[2] inds - indx(inds, k) res - NextMethod([, x) res } [-.symarray - function(x, inds, value){ dims - attr(x, dims) k - dims[2] inds - indx(inds, k) res - NextMethod([-, x) res } __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html SpG __ SpG R-help@stat.math.ethz.ch mailing list SpG https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do SpG read the posting guide! SpG http://www.R-project.org/posting-guide.html [ Part 3.15: Included Message ] Charles C. Berry(858) 534-2098 Dept of Family/Preventive Medicine E mailto:[EMAIL PROTECTED] UC San Diego http://biostat.ucsd.edu/~cberry/ La Jolla, San Diego 92093-0717 __ R-help@stat.math.ethz.ch 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 and ViM
dear all, to close the open question i asked more than a month ago, i would like to tell you my conclusions: exchanging emails with larry clapp who wrote funnel.pl, i found out that the tasks i wanted to be done by funnel.pl could actually, and much more, also be done by gnuscreen (http://www.gnu.org/software/screen/), which is a 'a full-screen window manager that multiplexes a physical terminal between several processes'. quite interesting thing. anyway. what i'm doing now to use vim with r is, to open a screen session with r and vim running within it: %screen R starts screen with R in it. then i open another window with 'ctrl-a c'. 'ctrl-a' start screen internal commands. in this new shell i can start vim normally to edit my .r-file. to change windows i use 'ctrl-a space'. now i can split the window with 'ctrl-a S'. to jump between the splitted halfs you can use 'ctrl-a tab'. i normally resize the r window with 'ctrl-a :resize 8', to make 8 lines high. to send text from one window to another i put into my .screenrc - file: bind y at R# paste . this makes it possible to use 'ctrl-a y' to send your buffer to the R window. to get something into your buffer you make use of the screen copy mode 'ctrl-a [' starts the copy mode. now you can start copying by pressing 'space', move the cursor and then pressing 'space' again. if you want to cut off the line numbers you can press 'v' before pressing 'space' the second time. now you can send the buffer as i told above. this works pretty well so far! what i would like to improve further is, to be able to copy whole lines and blocks with simple commands, maybe with sending it to r directly. any ideas? best, michael On Apr 19, 2006, at 8:42 AM, Martin Maechler wrote: PD == Peter Dalgaard [EMAIL PROTECTED] on 19 Apr 2006 01:06:02 +0200 writes: PD Jose Quesada [EMAIL PROTECTED] writes: Hmm, how timely. I posted yesterday my solution to integrate R and vim. The message is in my sent box but I don't think it showed up in the list... Here it is again: Hi All, If you use vim to edit R code, you may be interested in this. I have put together a personalized syntax file, some code templates, and a way to send code from Vim to R using autoHotKeys (windows). http://www.andrew.cmu.edu/user/jquesada/RvimSuite/instructions.html Actually, the little autoHotKeys can be useful even if you don't use vim just to send the example R code from the help pages to the console. Best wishes, -Jose PS: @list moderators Any idea why my message (from a @gmail account) appeared in the sent box but never on the list? PD I saw it the first time... Could it be that your mail PD reader is set up to ignore mail from yourself? PD Rather than resending stuff, it is preferable to check PD the archives PD https://stat.ethz.ch/pipermail/r-help/2006-April/date.html PD https://stat.ethz.ch/pipermail/r-help/2006-April/092457.html PD which are easily reachable via R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help Indeed. Please do check the archives. Now back to the subject: Jose, I think your main contribution is based on autoHotKeys and that only works on Windoze, right? Michael explicitly mentioned he's working in Mac OS X. Martin __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] contour help
Hi I'd like to draw some level plots but I just cant for the live of me figure out how to control the axis here. The axis for both the contour plots below appear to be the same, although one is a subset of the other. I was wondering how I could explicitly specify the xand y axis locations here. al dat-matrix(runif(50)*20, ncol = 5) contour(dat) contour(dat[1:2, 1:2]) __ R-help@stat.math.ethz.ch 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] contour help
Hi I'd like to draw some level plots but I just cant for the live of me figure out how to control the axis here. The axis for both the contour plots below appear to be the same, although one is a subset of the other. I was wondering how I could explicitly specify the xand y axis locations here. al dat-matrix(runif(50)*20, ncol = 5) contour(dat) contour(dat[1:2, 1:2]) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Is there a way to draw 3d plot?
It may be that since the original poster is coming from matlab others considered it best to start with the persp function that does not require additional packages or the other complications that come with lattice/trellis graphics. Actually my response indirectly brings cloud and wireframe into the discussion, if you look at the help page for rotate.persp in the TeachingDemos package you will also see help for rotate.cloud and rotate.wireframe functions. -Original Message- From: Ritwik Sinha [mailto:[EMAIL PROTECTED] Sent: Mon 5/29/2006 9:19 AM To: Greg Snow Cc: R-help@stat.math.ethz.ch Subject: Re: [R] Is there a way to draw 3d plot? Hi, I am not aware what this function does in matlab, but I was wondering why cloud or wireframe in the lattice package has not entered the discussion. Sorry if I am being naive. Ritwik. http://darwin.cwru.edu/~rsinha [[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] line over multiplot
It may be as simple as setting par(xpd=NA) see ?par for what that does. If you want more control (prevent the line from going through titles, etc) then look at the examples for the cnvrt.coords function in the TeachingDemos package. Hope this helps, -Original Message- From: [EMAIL PROTECTED] on behalf of Thomas Hoffmann Sent: Mon 5/29/2006 10:10 AM To: r-help@stat.math.ethz.ch Subject: [R] line over multiplot Dear Listmembers, I do have a multiplot with 1 row and 4 columns (4 scatter plots with same x- and y-axe plotted in one row). I would like to draw a horizontal line across the whole multiplot at a certain y-value. Unfortunately the normal abline command stops between the plot. I hope my problem is clear and there is somebody who can help me Thanks in advance Thomas H. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] OT: Monograph on Statistical Programming
I'd say most books on statistical computing would sort of fit your description, but few will have code examples. The ones I have (and liked) are Thisted, Monahan (which has Fortran code available), and the Handbook of Statistical Computing (edited volume). Prof. Gentle has at least two volumes in a series published so far as well (recently, not the one from 1980). I'm sure there are a few others. HTH, Andy _ From: [EMAIL PROTECTED] on behalf of Rau, Roland Sent: Mon 5/29/2006 9:37 AM To: r-help@stat.math.ethz.ch Subject: [R] OT: Monograph on Statistical Programming [Broadcast] Dear all, my question might be a bit off-topic. Is there anything like a standard textbook on statistical programming? With that I don't mean anything like MASS, S Programming, Programming with Data, ... (no offense meant, they are fantastic books and each of those three helped me a great deal). Rather in the direction of the Numerical Recipes Series addressing how to implement various functions. There are bits and pieces all over the place (The Art of Computer Programming, Numerical Recipes, Algorithms in ...(Sedgewick), etc.)... But is there one book like The Art of Statistical Programming? Maybe writing such a book could be impossible because the choice of implementation language could be crucial and also the choice of algorithms could be problematic. Nevertheless I hope that someone dared to tackle such a project. Thanks, Roland -- This mail has been sent through the MPI for Demographic Rese...{{dropped}} __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html 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] Two problems with write.foreign (SPSS)
Hello, Working on a little primer to ease the transition from SPSS to R, I have encountered two problems with write.foreign. One is cosmetic (but still annoying): the text data files will contain NA for missing values. Each time SPSS encounters an NA, it will print a warning. This could be easily avoided by calling write.table (which I suppose write.foreign uses) with na=,,. That's right, two commas. The second is more serious: character variables are not read at all because they are not specified correctly in the SPSS syntax file. A string variable has to be specified with a preceding asterisk and a type and length indicator. Suppose mydf contained a numerical variable a and a string variable b with max(nchar(mydf$b))==10. What write.foreign(mydf,...,package=SPSS) writes is: DATA FILE filename.dat free / a b. What it should write, however, is: DATA FILE filename.dat free / a * b (A10). __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html