Re: [R] comparing lm(), survreg( ... , dist="gaussian") and survreg( ... , dist="lognormal")
On Mon, 2 May 2005, Charles Annis, P.E. wrote: I have tried everything I can think of and hope not to appear too foolish when my error is pointed out to me. I have some real data (18 points) that look linear on a log-log plot so I used them for a comparison of lm() and survreg. There are no suspensions. survreg.df <- data.frame(Cycles=c(2009000, 577000, 145000, 376000, 37000, 979000, 1742, 71065000, 46397000, 70168000, 6912, 68798000, 72615000, 133051000, 38384000, 15204000, 1558000, 14181000), stress=c(90, 100, 110, 90, 100, 80, 70, 60, 56, 62, 62, 59, 56, 53, 59, 70, 90, 70), event=rep(1, 18)) sN.lm<- lm(log(Cycles) ~ log10(stress), data=survreg.df) and vvv gaussian.survreg<- survreg(formula=Surv(time=log(Cycles), event) ~ log10(stress), dist="gaussian", data=survreg.df) produce identical parameter estimates and differ slightly in the residual standard error and scale, which is accounted for by scale being the MLE and thus biased. Correcting by sqrt(18/16) produces agreement. Using predict() for the lm, and predict.survreg() for the survreg model and correcting for the differences in stdev, produces identical plots of the fit and the upper and lower confidence intervals. All of this is as it should be. I trust you called predict() on both and let R choose the method. And, vv lognormal.survreg<- survreg(formula=Surv(time=(Cycles), event) ~ log10(stress), dist="lognormal", data=survreg.df) produces summary() results that are identical to the earlier call to survreg(), except for the call, of course. The parameter estimates and SE are identical. Again this is as I would expect it. But since the call uses Cycles, rather than log(Cycles) predict.survreg() returns $fit in Cycles units, rather than logs, and of course the fits are identical when plotted on a log-log grid and also agree with lm() Here is the fly in the ointment: The upper and lower confidence intervals, based on the $se.fit for the dist="lognormal" are quite obviously different from the other two methods, and although I have tried everything I could imagine I cannot reconcile the differences. How did you do this? (BTW, I assume you mean upper and lower confidence limits< for the predicted means.) For the predictions and standard errors are (or should be) on the response scale, a non-linear function of the parameters. In that case it is normal to form confidence limits on the linear predictor scale and transform. I believe that the confidence bounds for both models should agree. After all, both calls to survreg() produce identical parameter estimates. They will, if computed on the same basis. On log-scale (to avoid large numbers) pr1 <- predict(lognormal.survreg, se.fit=T) log(cbind(pr1$fit - 1.96*pr1$se.fit, pr1$fit + 1.96*pr1$se.fit)) pr2 <- predict(gaussian.survreg, se.fit=T) cbind(pr2$fit - 1.96*pr2$se.fit, pr2$fit + 1.96*pr2$se.fit) are really pretty close. The main difference is a slight shift, which comes about because the mean of a log(X) is not log(mean(X)). Note that the second set at the preferred ones. Transforming to log scale before making the confidence limits: cbind(log(pr1$fit) - 1.96*pr1$se.fit/pr1$fit, log(pr1$fit) + 1.96*pr1$se.fit/pr1$fit) does give identical answers. Consider care is needed in interpreting what predict() is actually predicting in non-linear models. For both glm() and survreg() it is closer to the median of the uncertainty in the predictions than to the mean. -- 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] "Special" characters in URI
On Tue, 3 May 2005, Gorjanc Gregor wrote: I am crossposting this to R-help and BioC, since it is relevant to both groups. I don't see the relevance to R-help. But the answer to your subject is unambiguous: valid URLs do not contain `special' characters -- they must be encoded. See RFC1738 at e.g. ftp://ftp.funet.fi/pub/doc/rfc/rfc1738.txt At some point (probably 2.2.0) I intend to ensure that the mapping to file:// URLs that is done is a few places is encoded as necessary. This will likely result in a utility function filePathToURL or some such. -- 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] Combining numeric vs numeric & numeric vs factor graphs into one ps/pdf file
Dear R community, xyplot (lattice) has been great in displaying tons of data for my research. I have used the following two xyplot commands (with example dataframe) to create two separate postscript/pdf files with respect to the variable "acft" and subset "status": test.df [[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: Fwd: Re: [R] eigenvalues of a circulant matrix
I am not sure where your problem is, but I get these eigenvalues (no NaNs): > eigen(X)$values [1] 33.6609067 33.6609067 33.3657337 33.3657337 33.1688489 33.1688489 [7] 32.9120603 32.9120603 32.6971529 32.6971529 32.4486757 32.4486757 [13] 32.2266513 32.2266513 31.9811880 31.9811880 31.7549774 31.7549774 [19] 31.5106802 31.5106802 31.2814636 31.2814636 31.0374623 31.0374623 [25] 30.8058168 30.8058168 30.5616180 30.5616180 30.3278637 30.3278637 [31] 30.0831665 30.0831665 29.8474910 29.8474910 29.6020821 29.6020821 [37] 29.3645937 29.3645937 29.1183207 29.1183207 28.8790883 28.8790883 [43] 28.6318199 28.6318199 28.3908814 28.3908814 28.1425171 28.1425171 [49] 27.8998862 27.8998862 27.6503394 27.6503394 27.4060152 27.4060152 [55] 27.1552045 27.1552045 26.9091770 26.9091770 26.6570228 26.6570228 [61] 26.4092748 26.4092748 26.1557083 26.1557083 25.9062128 25.9062128 [67] 25.6511554 25.6511554 25.3998797 25.3998797 25.1432660 25.1432660 [73] 24.8901666 24.8901666 24.6319253 24.6319253 24.3769595 24.3769595 [79] 24.1170167 24.1170167 23.8601266 23.8601266 23.5984086 23.5984086 [85] 23.3395425 23.3395425 23.0759756 23.0759756 22.8150622 22.8150622 [91] 22.5495580 22.5495580 22.2865330 22.2865330 22.0190111 22.0190111 [97] 21.7537982 21.7537982 21.4841613 21.4841613 21.2166787 21.2166787 [103] 20.9448339 20.9448339 20.6749916 20.6749916 20.4008361 20.4008361 [109] 20.1285318 20.1285318 19.8519563 19.8519563 19.5770843 19.5770843 [115] 19.2979657 19.2979657 19.0204113 19.0204113 18.7386197 18.7386197 [121] 18.4582550 18.4582550 18.1736507 18.1736507 17.8903299 17.8903299 [127] 17.6027615 17.6027615 17.3163357 17.3163357 17.0256338 17.0256338 [133] 16.7359239 16.7359239 16.4419046 16.4419046 16.1487284 16.1487284 [139] 15.8511833 15.8511833 15.5543204 15.5543204 15.2530227 15.2530227 [145] 14.9522432 14.9522432 14.6469402 14.6469402 14.3419740 14.3419740 [151] 14.0323742 14.0323742 13.7229182 13.7229182 13.4087011 13.4087011 [157] 13.0944062 13.0944062 12.7751956 12.7751956 12.4556716 12.4556716 [163] 12.1310382 12.1310382 11.8058191 11.8058191 11.4752626 11.4752626 [169] 11.1438118 11.1438118 10.8067417 10.8067417 10.4684198 10.4684198 [175] 10.1241331 10.1241331 9.7781642 9.7781642 9.4258124 9.4258124 [181] 9.0712521 9.0712521 8.7097773 8.7097773 8.3454594 8.3454594 [187] 7.9735359 7.9735359 7.5979459 7.5979459 7.2138560 7.2138560 [193] 6.8250204 6.8250204 6.4264559 6.4264559 6.0216592 6.0216592 [199] 5.6053886 5.6053886 5.1807131 5.1807131 4.7419171 4.7419171 [205] 4.2912801 4.2912801 3.8221094 3.8221094 3.3349406 3.3349406 [211] 2.8206308 2.8206308 2.2746353 2.2746353 1.6787279 1.6787279 [217] 0.9988213 Again, what is the issue here? Janusz. On Mon, 2 May 2005, Globe Trotter wrote: > Looks like the files did not go through again. In any case, here is the kinv: > please cut and paste and save to a file: > > > >-1.16801E-03 -2.24310E-03 -1.16864E-03 -2.24634E-03 -1.17143E-03 >-2.25358E-03 -1.17589E-03 -2.26484E-03 -1.18271E-03 -2.27983E-03 >-1.19124E-03 -2.29896E-03 -1.20164E-03 -2.32206E-03 -1.21442E-03 >-2.34911E-03 -1.22939E-03 -2.38073E-03 -1.24626E-03 -2.41702E-03 >-1.26596E-03 -2.45828E-03 -1.28801E-03 -2.50458E-03 -1.31296E-03 >-2.55646E-03 -1.34048E-03 -2.61444E-03 -1.37127E-03 -2.67887E-03 >-1.40531E-03 -2.75026E-03 -1.44311E-03 -2.82930E-03 -1.48481E-03 >-2.91652E-03 -1.53081E-03 -3.01281E-03 -1.58131E-03 -3.11930E-03 >-1.63727E-03 -3.23708E-03 -1.69907E-03 -3.36712E-03 -1.76720E-03 >-3.51113E-03 -1.84251E-03 -3.67073E-03 -1.92580E-03 -3.84787E-03 >-2.01834E-03 -4.04507E-03 -2.12087E-03 -4.26509E-03 -2.23531E-03 >-4.51127E-03 -2.36357E-03 -4.78743E-03 -2.50664E-03 -5.09847E-03 >-2.66813E-03 -5.45027E-03 -2.85019E-03 -5.84987E-03 -3.05664E-03 >-6.30596E-03 -3.29224E-03 -6.82972E-03 -3.56187E-03 -7.43448E-03 >-3.87322E-03 -8.13766E-03 -4.23449E-03 -8.96182E-03 -4.65684E-03 >-9.93567E-03 -5.15519E-03 -1.10980E-02 -5.74887E-03 -1.25006E-02 >-6.46346E-03 -1.42143E-02 -7.33466E-03 -1.63391E-02 -8.41211E-03 >-1.90180E-02 -9.76709E-03 -2.24632E-02 -1.15055E-02 -2.70006E-02 >-1.37894E-02 -3.31497E-02 -1.68780E-02 -4.17865E-02 -2.12092E-02 >-5.44795E-02 -2.75722E-02 -7.42814E-02 -3.75180E-02 -0.107820 > -5.44778E-02 > -0.171908 -8.74660E-02 -0.320830 -0.167792 -0.826079 -0.486405 > -6.28085 > 19.4495 -6.28085 -0.486405 -0.826079 -0.167792 -0.320830 > -8.74660E-02 > -0.171908 -5.44778E-02 -0.107820 -3.75180E-02 -7.42816E-02 > -2.75722E-02 >-5.44795E-02 -2.12092E-02 -4.17865E-02 -1.68780E-02 -3.31495E-02 >-1.37894E-02 -2.70004E-02 -1.15055E-02 -2.24632E-02 -9.76708E-03 >-1.90179E-02 -8.41210E-03 -1.63391E-02 -7.33465E-03 -1.42142E-02 >-6.46346E-03 -1.
Fwd: Re: [R] eigenvalues of a circulant matrix
Looks like the files did not go through again. In any case, here is the kinv: please cut and paste and save to a file: -1.16801E-03 -2.24310E-03 -1.16864E-03 -2.24634E-03 -1.17143E-03 -2.25358E-03 -1.17589E-03 -2.26484E-03 -1.18271E-03 -2.27983E-03 -1.19124E-03 -2.29896E-03 -1.20164E-03 -2.32206E-03 -1.21442E-03 -2.34911E-03 -1.22939E-03 -2.38073E-03 -1.24626E-03 -2.41702E-03 -1.26596E-03 -2.45828E-03 -1.28801E-03 -2.50458E-03 -1.31296E-03 -2.55646E-03 -1.34048E-03 -2.61444E-03 -1.37127E-03 -2.67887E-03 -1.40531E-03 -2.75026E-03 -1.44311E-03 -2.82930E-03 -1.48481E-03 -2.91652E-03 -1.53081E-03 -3.01281E-03 -1.58131E-03 -3.11930E-03 -1.63727E-03 -3.23708E-03 -1.69907E-03 -3.36712E-03 -1.76720E-03 -3.51113E-03 -1.84251E-03 -3.67073E-03 -1.92580E-03 -3.84787E-03 -2.01834E-03 -4.04507E-03 -2.12087E-03 -4.26509E-03 -2.23531E-03 -4.51127E-03 -2.36357E-03 -4.78743E-03 -2.50664E-03 -5.09847E-03 -2.66813E-03 -5.45027E-03 -2.85019E-03 -5.84987E-03 -3.05664E-03 -6.30596E-03 -3.29224E-03 -6.82972E-03 -3.56187E-03 -7.43448E-03 -3.87322E-03 -8.13766E-03 -4.23449E-03 -8.96182E-03 -4.65684E-03 -9.93567E-03 -5.15519E-03 -1.10980E-02 -5.74887E-03 -1.25006E-02 -6.46346E-03 -1.42143E-02 -7.33466E-03 -1.63391E-02 -8.41211E-03 -1.90180E-02 -9.76709E-03 -2.24632E-02 -1.15055E-02 -2.70006E-02 -1.37894E-02 -3.31497E-02 -1.68780E-02 -4.17865E-02 -2.12092E-02 -5.44795E-02 -2.75722E-02 -7.42814E-02 -3.75180E-02 -0.107820 -5.44778E-02 -0.171908 -8.74660E-02 -0.320830 -0.167792 -0.826079 -0.486405 -6.28085 19.4495 -6.28085 -0.486405 -0.826079 -0.167792 -0.320830 -8.74660E-02 -0.171908 -5.44778E-02 -0.107820 -3.75180E-02 -7.42816E-02 -2.75722E-02 -5.44795E-02 -2.12092E-02 -4.17865E-02 -1.68780E-02 -3.31495E-02 -1.37894E-02 -2.70004E-02 -1.15055E-02 -2.24632E-02 -9.76708E-03 -1.90179E-02 -8.41210E-03 -1.63391E-02 -7.33465E-03 -1.42142E-02 -6.46346E-03 -1.25005E-02 -5.74887E-03 -1.10979E-02 -5.15519E-03 -9.93565E-03 -4.65686E-03 -8.96166E-03 -4.23449E-03 -8.13766E-03 -3.87320E-03 -7.43444E-03 -3.56187E-03 -6.82968E-03 -3.29223E-03 -6.30590E-03 -3.05665E-03 -5.84988E-03 -2.85020E-03 -5.45026E-03 -2.66813E-03 -5.09848E-03 -2.50664E-03 -4.78742E-03 -2.36358E-03 -4.51160E-03 -2.23531E-03 -4.26516E-03 -2.12088E-03 -4.04506E-03 -2.01834E-03 -3.84791E-03 -1.92580E-03 -3.67066E-03 -1.84251E-03 -3.51094E-03 -1.76720E-03 -3.36700E-03 -1.69907E-03 -3.23702E-03 -1.63727E-03 -3.11926E-03 -1.58130E-03 -3.01269E-03 -1.53081E-03 -2.91633E-03 -1.48480E-03 -2.82912E-03 -1.44311E-03 -2.75018E-03 -1.40531E-03 -2.67891E-03 -1.37126E-03 -2.61459E-03 -1.34048E-03 -2.55652E-03 -1.31296E-03 -2.50461E-03 -1.28801E-03 -2.45833E-03 -1.26595E-03 -2.41717E-03 -1.24626E-03 -2.38063E-03 -1.22938E-03 -2.34904E-03 -1.21443E-03 -2.32207E-03 -1.20164E-03 -2.29903E-03 -1.19124E-03 -2.27985E-03 -1.18271E-03 -2.26497E-03 -1.17589E-03 -2.25364E-03 -1.17145E-03 -2.24634E-03 -1.16863E-03 -2.24284E-03 Thanks! --- Globe Trotter <[EMAIL PROTECTED]> wrote: > Date: Mon, 2 May 2005 19:51:24 -0700 (PDT) > From: Globe Trotter <[EMAIL PROTECTED]> > Subject: Re: [R] eigenvalues of a circulant matrix > To: r-help@stat.math.ethz.ch > > OK, here we go: > > I am submitting two attachments. The first is the datafile called kinv used > to > create my circulant matrix, using the following commands: > > > x<-scan("kinv") > y<-x[c(109:1,0:108)] > X=toeplitz(y) > eigen(X) > write(X,ncol=216,file="test.dat") > > reports the following columns full of NaN's: 18, 58, 194, 200. (Note that > eigen(X,symmetric=T) makes no difference and I get the same as above). > > The second attachment contains only the eigenvectors obtained on calling a > LAPACK routine directly (from C). The eigenvalues are essentially the same as > that obtained using R. Here, I use the LAPACK-recommended double precision > routine dspevd() routine for symmetric matrices in packed storage format. > Note > the absence of the NaN'sI would be happy to send my C programs to whoever > is interested. > > I am using > > :~> uname -a > Linux 2.6.11-1.14_FC3 #1 Thu Apr 7 19:23:49 EDT 2005 i686 i686 i386 GNU/Linux > > and R.2.0.1. > > Many thanks and best wishes! > > __ > R-help@stat.math.ethz.ch 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/
[R] comparing lm(), survreg( ... , dist="gaussian") and survreg( ... , dist="lognormal")
Dear R-Helpers: I have tried everything I can think of and hope not to appear too foolish when my error is pointed out to me. I have some real data (18 points) that look linear on a log-log plot so I used them for a comparison of lm() and survreg. There are no suspensions. survreg.df <- data.frame(Cycles=c(2009000, 577000, 145000, 376000, 37000, 979000, 1742, 71065000, 46397000, 70168000, 6912, 68798000, 72615000, 133051000, 38384000, 15204000, 1558000, 14181000), stress=c(90, 100, 110, 90, 100, 80, 70, 60, 56, 62, 62, 59, 56, 53, 59, 70, 90, 70), event=rep(1, 18)) sN.lm<- lm(log(Cycles) ~ log10(stress), data=survreg.df) and vvv gaussian.survreg<- survreg(formula=Surv(time=log(Cycles), event) ~ log10(stress), dist="gaussian", data=survreg.df) produce identical parameter estimates and differ slightly in the residual standard error and scale, which is accounted for by scale being the MLE and thus biased. Correcting by sqrt(18/16) produces agreement. Using predict() for the lm, and predict.survreg() for the survreg model and correcting for the differences in stdev, produces identical plots of the fit and the upper and lower confidence intervals. All of this is as it should be. And, vv lognormal.survreg<- survreg(formula=Surv(time=(Cycles), event) ~ log10(stress), dist="lognormal", data=survreg.df) produces summary() results that are identical to the earlier call to survreg(), except for the call, of course. The parameter estimates and SE are identical. Again this is as I would expect it. But since the call uses Cycles, rather than log(Cycles) predict.survreg() returns $fit in Cycles units, rather than logs, and of course the fits are identical when plotted on a log-log grid and also agree with lm() Here is the fly in the ointment: The upper and lower confidence intervals, based on the $se.fit for the dist="lognormal" are quite obviously different from the other two methods, and although I have tried everything I could imagine I cannot reconcile the differences. I believe that the confidence bounds for both models should agree. After all, both calls to survreg() produce identical parameter estimates. So I have missed something. Would some kind soul please point out my error? Thanks. Charles Annis, P.E. [EMAIL PROTECTED] phone: 561-352-9699 eFax: 614-455-3265 http://www.StatisticalEngineering.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] maximization help :
Given a vector : pvec=(p1,p2, p J) with sum(pvec)=1, all the elements are non-negative, that is, they are probabilities a matrix A ( N* J ), with the elements alpha(ij) are 0 or 1 I want to MAXIMIZE THE RESULT RESULT= product( i=1, to N [ sum ( alpha(ij)* pj , j =1,to J ) ] ) thus, I need to get pvec. how should I do ? for example say, A= 0 1 0 0 11 0 0 1 0 0 0 0 0 1 0 1 00 1 0 00 1 that is A is a matrix 6* 4thus pvec=(p1,p2,p3,p4) I want to get values of pvec such that , they can maximize p2 * ( p1 + p2 ) * p1 * p3 * (p1+p4) * p4 thanks __ R-help@stat.math.ethz.ch 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] eigenvalues of a circulant matrix
OK, here we go: I am submitting two attachments. The first is the datafile called kinv used to create my circulant matrix, using the following commands: x<-scan("kinv") y<-x[c(109:1,0:108)] X=toeplitz(y) eigen(X) write(X,ncol=216,file="test.dat") reports the following columns full of NaN's: 18, 58, 194, 200. (Note that eigen(X,symmetric=T) makes no difference and I get the same as above). The second attachment contains only the eigenvectors obtained on calling a LAPACK routine directly (from C). The eigenvalues are essentially the same as that obtained using R. Here, I use the LAPACK-recommended double precision routine dspevd() routine for symmetric matrices in packed storage format. Note the absence of the NaN'sI would be happy to send my C programs to whoever is interested. I am using :~> uname -a Linux 2.6.11-1.14_FC3 #1 Thu Apr 7 19:23:49 EDT 2005 i686 i686 i386 GNU/Linux and R.2.0.1. Many thanks and best wishes! __ R-help@stat.math.ethz.ch 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] xinch/yinch equivalent for log axis
Hi Steven J. Murdoch wrote: I would like to draw a vertical line from a given point, in user coordinates, to x inches before from another point, also in user coordinates. This is easy enough to do for linear scales, using code based on xinch/yinch, but I do not know how to do this for logarithmic scales. This code shows an example of what I mean[1]: split.screen(c(2,1)) screen(1) # Linear scale, works fine plot(1:100, cex=0.5, pch=19) plotheight <- diff(par("usr")[3:4]) igap = 0.5 # intended gap, in inches gap <- igap/par("pin")[2]*plotheight # gap in user units lines(c(20,20),c(100,60+gap)) lines(c(20,20),c(60-gap,1)) points(20,60, pch=18) screen(2) # Logarithmic scale, point no longer centered plot(1:100, cex=0.5, pch=19, log="y") plotheight <- diff(par("usr")[3:4]) igap = 10 gap <- igap/par("pin")[2]*plotheight lines(c(20,20),c(100,10+gap)) lines(c(20,20),c(10-gap,1)) points(20,10, pch=18) close.screen(all=TRUE) The top graph is as I want. The diamond is centered in the gap, and the gap is 1 inch high (2*igap). In bottom graph, using a log scale, the diamond is no longer centered in the gap and there is a non-linear relationship between the gap height in inches and igap. I understand why this is happening, and this is why the xinch and yinch functions raise a warning, but is there a way to handle it? For example, is there a function which will return the user coordinate of a point, x inches above a given point p (in user coordinates), for a logarithmic y axis? No function, but does this do what you want? The basic idea is to work in logged values (which is what par("usr") is in) then convert back to plot (broken into many steps to hopefully aid clarity) ... plot(1:100, cex=0.5, pch=19, log="y") plotheight <- diff(par("usr")[3:4]) igap <- 0.5 # intended gap, in inches pigap <- igap/par("pin")[2] # intended gap as propn of plot height dy <- 10 # y location for diamond ldy <- log10(dy) # logged y location for diamond upper <- ldy + pigap/2*plotheight # logged y location for top of gap lower <- ldy - pigap/2*plotheight # logged y location for bottom of gap lines(c(20,20),c(100, 10^upper)) lines(c(20,20),c(10^lower, 1)) points(20,10, pch=18) Paul [1] I actually want to use this in a larger script, where I leave a gap in the axis where the median in. This needs to be small but legible to the eye, so that is why I am defining it in inches. The source code is at: http://www.cl.cam.ac.uk/users/sjm217/projects/graphics/fancyaxis.R -- Dr Paul Murrell Department of Statistics The University of Auckland Private Bag 92019 Auckland New Zealand 64 9 3737599 x85392 [EMAIL PROTECTED] http://www.stat.auckland.ac.nz/~paul/ __ R-help@stat.math.ethz.ch 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] bivariate log-normal distribution
Dimitri Joe wrote on 5/2/2005 5:29 PM: Hi, I am looking for the density function of the bivariate log-normal distribution. Would anyone have it written? Thanks a lot, Dimitri [[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 ^ This is in the archive. http://finzi.psych.upenn.edu/R/Rhelp02a/archive/48580.html --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
[R] bivariate log-normal distribution
Hi, I am looking for the density function of the bivariate log-normal distribution. Would anyone have it written? Thanks a lot, Dimitri [[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] Nonparametric Tukey-type multiple comparisons "Nemenyi" test
help.search("Tukey"). The fifth hit is TukeyHSD(). Also take a look at tmd(lattice) Cheers Francisco From: "Rikki Dunsmore" <[EMAIL PROTECTED]> To: r-help@stat.math.ethz.ch Subject: [R] Nonparametric Tukey-type multiple comparisons "Nemenyi" test Date: Mon, 02 May 2005 16:45:20 -0700 I am trying to do a Nonparametric Tukey-type multiple comparison post-hoc test to determine which groups are significantly different. I have read the dialogue on this topic from the R-help, and am still not clear why no statistical packages include this test as an option? Is it not an appropriate test to conduct on non-normally distributed data? Is the only option to calculate it by hand using the (Zar 1996) formula? Thank you in advance for your help. -- Rikki Grober- Dunsmore National Marine Fisheries Service National Marine Protected Areas Center 110 Shaffer Rd. Santa Cruz, CA 95060 831-420-3991 Unless someone like you, Cares a whole awful lot, Nothing is going to get better, It's not. - The Lorax, by Dr. Seuss, 1971 __ R-help@stat.math.ethz.ch 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] "Special" characters in URI
Hello! I am crossposting this to R-help and BioC, since it is relevant to both groups. I wrote a wrapper for Entrez search utility (link for this is provided bellow), which can add some new search functionality to existing code in Bioconductor's package 'annotate'*. http://eutils.ncbi.nlm.nih.gov/entrez/query/static/esearch_help.html Entrez search utuility returns a XML document but I have a problem to use URI to retrieve that file, since URI can also contain characters, which should not be there according to http://www.faqs.org/rfcs/rfc2396.html I encountered problems with "[" and "]" as well as with space characters. However there might also be a problem with others i.e. reserved characters in URI syntax. My R example is: R> library("annotate") Loading required package: Biobase Loading required package: tools Welcome to Bioconductor Vignettes contain introductory material. To view, simply type: openVignette() For details on reading vignettes, see the openVignette help page. R> library(XML) R> tmp$term <- "gorjanc g[au]" R> tmp$URL <- "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc g[au]" R> tmp $term [1] "gorjanc g[au]" $URL [1] "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc g[au]" R> xmlTreeParse(tmp$URL, isURL=TRUE, handlers=NULL, asTree=TRUE) Error in xmlTreeParse(tmp$URL, isURL = TRUE, handlers = NULL, asTree = TRUE) : error in creating parser for http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc g[au] # so I have a problem with space and [ and ] # let's reduce a problem to just space or [] to be sure R> tmp$URL <- "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc g" R> xmlTreeParse(tmp$URL, isURL=TRUE, handlers=NULL, asTree=TRUE) Error in xmlTreeParse(tmp$URL, isURL = TRUE, handlers = NULL, asTree = TRUE) : error in creating parser for http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc g R> tmp$URL <- "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc[au]"; R> xmlTreeParse(tmp$URL, isURL=TRUE, handlers=NULL, asTree=TRUE) Error in xmlTreeParse(tmp$URL, isURL = TRUE, handlers = NULL, asTree = TRUE) : error in creating parser for http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc[au] # now show that it works fine without special chars R> tmp$URL <- "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc"; R> xmlTreeParse(tmp$URL, isURL=TRUE, handlers=NULL, asTree=TRUE) $doc $file [1] "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc"; $version [1] "1.0" $children ... # now show a workaround for space tmp$URL <- "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc%20g"; xmlTreeParse(tmp$URL, isURL=TRUE, handlers=NULL, asTree=TRUE) R> tmp$URL <- "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc%20g"; R> xmlTreeParse(tmp$URL, isURL=TRUE, handlers=NULL, asTree=TRUE) $doc $file [1] "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/esearch.fcgi?term=gorjanc%20g"; $version [1] "1.0" $children ... As can be seen from above there is a possibility to handle this special characters and I wonder if this has already been done somewhere? If not I thought on a function fixURLchar, which would replace reserved characters with ther escaped sequences. Any comments, pointers, ... ? from = c(" ", "\"", ",", "#"), to = c("%20", "%22", "%2c", "%23")) *When I'll solve problem I will send my code to 'annotate' maintainer and he can include it at his will in a package. Lep pozdrav / With regards, Gregor Gorjanc -- University of Ljubljana Biotechnical FacultyURI: http://www.bfro.uni-lj.si/MR/ggorjan Zootechnical Department mail: gregor.gorjanc bfro.uni-lj.si Groblje 3 tel: +386 (0)1 72 17 861 SI-1230 Domzale fax: +386 (0)1 72 17 888 Slovenia, Europe -- "One must learn by doing the thing; for though you think you know it, you have no certainty until you try." Sophocles ~ 450 B.C. __ R-help@stat.math.ethz.ch 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] Nonparametric Tukey-type multiple comparisons "Nemenyi" test
I am trying to do a Nonparametric Tukey-type multiple comparison post-hoc test to determine which groups are significantly different. I have read the dialogue on this topic from the R-help, and am still not clear why no statistical packages include this test as an option? Is it not an appropriate test to conduct on non-normally distributed data? Is the only option to calculate it by hand using the (Zar 1996) formula? Thank you in advance for your help. -- Rikki Grober- Dunsmore National Marine Fisheries Service National Marine Protected Areas Center 110 Shaffer Rd. Santa Cruz, CA 95060 831-420-3991 Unless someone like you, Cares a whole awful lot, Nothing is going to get better, It's not. - The Lorax, by Dr. Seuss, 1971 __ R-help@stat.math.ethz.ch 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] plot
Thanks a lot! Hui Berton Gunter wrote: ?points or ?xyplot in lattice package for a possibly better alternative -- 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 Hui Han Sent: Monday, May 02, 2005 4:08 PM To: R-help@stat.math.ethz.ch Subject: [R] plot Hi, I want to plot two sets of (x, y) in the same window, that is, plot of (x1,y1) and plot of (x2, y2) appear in the same window, so that I can compare the two plots. Can I get any help from you on how to do that? Thanks! Hui __ R-help@stat.math.ethz.ch 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] plot
?points or ?xyplot in lattice package for a possibly better alternative -- 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 Hui Han > Sent: Monday, May 02, 2005 4:08 PM > To: R-help@stat.math.ethz.ch > Subject: [R] plot > > Hi, > > I want to plot two sets of (x, y) in the same window, that > is, plot of > (x1,y1) and plot of (x2, y2) appear in the same window, so that I can > compare the two plots. > Can I get any help from you on how to do that? > > Thanks! > Hui > > __ > R-help@stat.math.ethz.ch 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] plot
Hi, I want to plot two sets of (x, y) in the same window, that is, plot of (x1,y1) and plot of (x2, y2) appear in the same window, so that I can compare the two plots. Can I get any help from you on how to do that? Thanks! Hui __ R-help@stat.math.ethz.ch 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] Reading in a dataset with uneven variable lengths
Is it OK to have NA's? If so you can use read.table(,fill=TRUE) or read.delim(). Notice that the default in the second option is fill=TRUE. Cheers Francisco From: "Liaw, Andy" <[EMAIL PROTECTED]> To: "'Owen, Jason'" <[EMAIL PROTECTED]>, r-help@stat.math.ethz.ch Subject: RE: [R] Reading in a dataset with uneven variable lengths Date: Mon, 2 May 2005 16:24:29 -0400 One possibility is to use read.fwf(). Andy > From: Owen, Jason > > Suppose I have a text file that I want to read into R like the > following: > > X Y > 649699 > 657891 > 714632 > 849727 > 721597 > 791868 > 874652 > 405978 >733 >549 >790 > > This is a simple example -- I could have a huge file with many > columns of unequal lengths. > > What is the best way to do it? I can't see how a data frame > can be used. I checked the FAQ and did a web search on the topic > but I came up empty. > > Jason > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! > http://www.R-project.org/posting-guide.html > > > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] summary(as.factor(x) - force to not sort the result according factor levels
Christoph Lehmann wrote on 5/2/2005 3:29 PM: Hi The result of a summary(as.factor(x)) (see example below) call is sorted according to the factor level. How can I get the result not sorted but in the original order of the levels in x? > test <- c(120402, 120402, 120402, 1323, 1323,200393, 200393, 200393, 200393, 200393) > summary(as.factor(test)) 1323 120402 200393 2 3 5 I need: 120402 1323 200393 32 5 thanks for a hint christoph Use the levels argument of ?factor. summary(factor(test, levels = unique(test))) --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] eigenvalues of a circulant matrix
The construction was y<-x[c(109:1,2:108)] so y is symmetric in the sense of the usual way of writing a function on integers mod n as a vector with 1-based indexing. I.e., y[i+1] = y[n-(i+1)] for i=0,1,...,n-1. So the assignment Z <- toeplitz(y) *does* create a symmetric circulant matrix. It is diagonalizable but does not have distinct eigenvalues, hence the eigenspaces may be more than one-dimensional, so you can't just pick a unit vector and call it "the" eigenvector for that eigenvalue. You choose a basis for each eigenspace. R detects the symmetry: ... symmetric: if `TRUE', the matrix is assumed to be symmetric (or Hermitian if complex) and only its lower triangle is used. If `symmetric' is not specified, the matrix is inspected for symmetry. (from help(eigen)) and knows that computations can be done with real arithmetic. As for why you get NaN, you should submit --along with your example-- details of your platform (machine, R version, how R was built and installed, etc). Reid Huntsinger -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of [EMAIL PROTECTED] Sent: Monday, May 02, 2005 5:28 PM To: r-help@stat.math.ethz.ch Subject: Re: [R] eigenvalues of a circulant matrix On 02-May-05 Ted Harding wrote: > On 02-May-05 Rolf Turner wrote: >> I just Googled around a bit and found definitions of Toeplitz and >> circulant matrices as follows: >> [...] >> A circulant matrix is an n x n matrix whose rows are composed of >> cyclically shifted versions of a length-n vector. For example, the >> circulant matrix on the vector (1, 2, 3, 4) is >> >> 4 1 2 3 >> 3 4 1 2 >> 2 3 4 1 >> 1 2 3 4 >> >> So circulant matrices are a special case of Toeplitz matrices. >> However a circulant matrix cannot be symmetric. > > I suspect the confusion may lie in what's meant by "cyclically > shifted". In Rolf's example above, each row is shifted right by 1 > and the one that falls off the end is put at the beginning. This > cannot be symmetric for general values in the fist row. > > However, if you shift left instead, then you get > > 4 1 2 3 > 1 2 3 4 > 2 3 4 1 > 3 4 1 2 > > and this *is* symmetric (and indeed will always be so, for > general values in the first row). I just had a look at ?toeplitz (We should have done that earlier!) toeplitzpackage:statsR Documentation Form Symmetric Toeplitz Matrix * Description: Forms a symmetric Toeplitz matrix given its first row. * [...] Examples: x <- 1:5 toeplitz (x) > x <- 1:5 > toeplitz (x) [,1] [,2] [,3] [,4] [,5] [1,]12345 [2,]21234 [3,]32123 [4,]43212 [5,]54321 Since "Globe Trotter's" construction was Y<-toeplitz(x) it's not surprising what he got (and it *certainly* wasn't a circulant!!!). Everybody barking up the wring tree here! Best wishes to all, Ted. E-Mail: (Ted Harding) <[EMAIL PROTECTED]> Fax-to-email: +44 (0)870 094 0861 Date: 02-May-05 Time: 22:27:32 -- XFMail -- __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch 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] Trying to understand kpss.test() in tseries package
I'm trying to understand how to use kpss.test() properly. If I have a level stationary series like rnorm() in the help page, shouldn't I get a small p-value with the null hypothesis set to "Trend"? The (condensed) output from kpss.test() for the two possible null hypotheses is given below. I don't see any significant difference between these results. > x <- rnorm(1000) # is level stationary > kpss.test(x, null="Level") KPSS Test for Level Stationarity KPSS Level = 0.0638, Truncation lag parameter = 7, p-value = 0.1 Warning: p-value greater than printed p-value > kpss.test(x, null="Trend") KPSS Test for Trend Stationarity KPSS Trend = 0.0275, Truncation lag parameter = 7, p-value = 0.1 Warning: p-value greater than printed p-value I can't get the original reference easily. Scott Waichler Pacific Northwest National Laboratory [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] making a log file with error messages?
One way is to define your own function to use as the "error" option value. E.g. my.error.fun <- function() { cat(geterrmessage(), file="rerr.txt", append=T) } then, options("error"=my.error.fun) Norm -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Chaehyung Ahn Sent: Monday, May 02, 2005 1:32 PM To: r-help@stat.math.ethz.ch Subject: [R] making a log file with error messages? Dear all: Is there anyone who would help me to generate a log file reporting error messages? I can see error messages on the screen, but I would like to report them as a file. thanks Sincerely, cahn --- ChaeHyung Ahn (cahn) Ph.D. Department of Biostatistics, School of Public Health CB #7420, University of North Carolina Chapel Hill, NC 27599-7420 phone: 919-966-7284 home: http://cyworld.nate.com/cahn88 __ R-help@stat.math.ethz.ch 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] summary(as.factor(x) - force to not sort the result according factor levels
Christoph Lehmann wrote: Hi The result of a summary(as.factor(x)) (see example below) call is sorted according to the factor level. How can I get the result not sorted but in the original order of the levels in x? by creating the factor with the levels in the order you want: > test <- c(120402, 120402, 120402, 1323, 1323,200393, 200393, 200393, 200393, 200393) > summary(factor(test, levels=unique(test))) 120402 1323 200393 3 2 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] eigenvalues of a circulant matrix
On 02-May-05 Ted Harding wrote: > On 02-May-05 Rolf Turner wrote: >> I just Googled around a bit and found definitions of Toeplitz and >> circulant matrices as follows: >> [...] >> A circulant matrix is an n x n matrix whose rows are composed of >> cyclically shifted versions of a length-n vector. For example, the >> circulant matrix on the vector (1, 2, 3, 4) is >> >> 4 1 2 3 >> 3 4 1 2 >> 2 3 4 1 >> 1 2 3 4 >> >> So circulant matrices are a special case of Toeplitz matrices. >> However a circulant matrix cannot be symmetric. > > I suspect the confusion may lie in what's meant by "cyclically > shifted". In Rolf's example above, each row is shifted right by 1 > and the one that falls off the end is put at the beginning. This > cannot be symmetric for general values in the fist row. > > However, if you shift left instead, then you get > > 4 1 2 3 > 1 2 3 4 > 2 3 4 1 > 3 4 1 2 > > and this *is* symmetric (and indeed will always be so, for > general values in the first row). I just had a look at ?toeplitz (We should have done that earlier!) toeplitzpackage:statsR Documentation Form Symmetric Toeplitz Matrix * Description: Forms a symmetric Toeplitz matrix given its first row. * [...] Examples: x <- 1:5 toeplitz (x) > x <- 1:5 > toeplitz (x) [,1] [,2] [,3] [,4] [,5] [1,]12345 [2,]21234 [3,]32123 [4,]43212 [5,]54321 Since "Globe Trotter's" construction was Y<-toeplitz(x) it's not surprising what he got (and it *certainly* wasn't a circulant!!!). Everybody barking up the wring tree here! Best wishes to all, Ted. E-Mail: (Ted Harding) <[EMAIL PROTECTED]> Fax-to-email: +44 (0)870 094 0861 Date: 02-May-05 Time: 22:27:32 -- XFMail -- __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] summary(as.factor(x) - force to not sort the result according factor levels
Hi The result of a summary(as.factor(x)) (see example below) call is sorted according to the factor level. How can I get the result not sorted but in the original order of the levels in x? > test <- c(120402, 120402, 120402, 1323, 1323,200393, 200393, 200393, 200393, 200393) > summary(as.factor(test)) 1323 120402 200393 2 3 5 I need: 120402 1323 200393 32 5 thanks for a hint christoph __ R-help@stat.math.ethz.ch 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] eigenvalues of a circulant matrix
The example that I submitted earlier in the day. Would you like me to send again? Thanks! --- "Huntsinger, Reid" <[EMAIL PROTECTED]> wrote: > For which X? > > Reid Huntsinger > > -Original Message- > From: Globe Trotter [mailto:[EMAIL PROTECTED] > Sent: Monday, May 02, 2005 2:34 PM > To: Huntsinger, Reid; Rolf Turner > Cc: r-help@stat.math.ethz.ch > Subject: RE: [R] eigenvalues of a circulant matrix > > > By the way, I just noticed that eigen(X) returns eigenvectors, at least two > of > which are NaN's. > > Best wishes! > > --- "Huntsinger, Reid" <[EMAIL PROTECTED]> wrote: > > > When the matrix is symmetric and omega is not real, omega and its > conjugate > > (= inverse) give the same eigenvalue, so you have a 2-dimensional > > eigenspace. R chooses a real basis of this, which is perfectly fine since > > it's not looking for circulant structure. > > > > For example, > > > > > m > > [,1] [,2] [,3] [,4] [,5] > > [1,]12332 > > [2,]21233 > > [3,]32123 > > [4,]33212 > > [5,]23321 > > > > > eigen(m) > > $values > > [1] 11.00 -0.381966 -0.381966 -2.618034 -2.618034 > > > > $vectors > > [,1] [,2] [,3] [,4] [,5] > > [1,] 0.4472136 0.00 -0.6324555 0.6324555 0.00 > > [2,] 0.4472136 0.371748 0.5116673 0.1954395 0.601501 > > [3,] 0.4472136 -0.601501 -0.1954395 -0.5116673 0.371748 > > [4,] 0.4472136 0.601501 -0.1954395 -0.5116673 -0.371748 > > [5,] 0.4472136 -0.371748 0.5116673 0.1954395 -0.601501 > > > > and you can match these columns up with the "canonical" eigenvectors > > exp(2*pi*1i*(0:4)*j/5) for j = 0,1,2,3,4. E.g., > > > > > Im(exp(2*pi*1i*(0:4)*3/5)) > > [1] 0.000 -0.5877853 0.9510565 -0.9510565 0.5877853 > > > > which can be seen to be a scalar multiple of column 2. > > > > Reid Huntsinger > > > > Reid Huntsinger > > > > -Original Message- > > From: [EMAIL PROTECTED] > > [mailto:[EMAIL PROTECTED] On Behalf Of Huntsinger, Reid > > Sent: Monday, May 02, 2005 10:43 AM > > To: 'Globe Trotter'; Rolf Turner > > Cc: r-help@stat.math.ethz.ch > > Subject: RE: [R] eigenvalues of a circulant matrix > > > > > > It's hard to argue against the fact that a real symmetric matrix has real > > eigenvalues. The eigenvalues of the circulant matrix with first row v are > > *polynomials* (not the roots of 1 themselves, unless as Rolf suggested you > > start with a vector with all zeros except one 1) in the roots of 1, with > > coefficients equal to the entries in v. This is the finite Fourier > transform > > of v, by the way, and takes real values when the coefficients are real and > > symmetric, ie when the matrix is symmetric. > > > > Reid Huntsinger > > > > -Original Message- > > From: [EMAIL PROTECTED] > > [mailto:[EMAIL PROTECTED] On Behalf Of Globe Trotter > > Sent: Monday, May 02, 2005 10:23 AM > > To: Rolf Turner > > Cc: r-help@stat.math.ethz.ch > > Subject: Re: [R] eigenvalues of a circulant matrix > > > > > > > > --- Rolf Turner <[EMAIL PROTECTED]> wrote: > > > I just Googled around a bit and found definitions of Toeplitz and > > > circulant matrices as follows: > > > > > > A Toeplitz matrix is any n x n matrix with values constant along each > > > (top-left to lower-right) diagonal. matrix has the form > > > > > > a_0 a_1 . . . . ... a_{n-1} > > > a_{-1} a_0 a_1... a_{n-2} > > > a_{-2} a_{-1} a_0 a_1 .... > > > . .. . . . > > > . .. . . . > > > . .. . . . > > > a_{-(n-1)} a_{-(n-2)} ... a_1 a_0 > > > > > > (A Toeplitz matrix ***may*** be symmetric.) > > > > Agreed. As may a circulant matrix if a_i = a_{p-i+2} > > > > > > > > A circulant matrix is an n x n matrix whose rows are composed of > > > cyclically shifted versions of a length-n vector. For example, the > > > circulant matrix on the vector (1, 2, 3, 4) is > > > > > > 4 1 2 3 > > > 3 4 1 2 > > > 2 3 4 1 > > > 1 2 3 4 > > > > > > So circulant matrices are a special case of Toeplitz matrices. > > > However a circulant matrix cannot be symmetric. > > > > > > The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i, > > > 2 - 2i, and 2 --- certainly not roots of unity. > > > > The eigenvalues are 4+1*omega+2*omega^2+3*omega^3. > > omega=cos(2*pi*k/4)+isin(2*pi*k/4) as k ranges over 1, 2, 3, 4, so the > above > > holds. > > > > Bellman may have > > > been talking about the particular (important) case of a circulant > > > matrix where the vector from which it is constructed is a canonical > > > basis vector e_i with a 1 in the i-th slot and zeroes elsewhere. > > > > No, that is not true: his result can be verified for any circulant matrix, > > directly. > > > > > Such a matrix is in fact a unitary matrix (operator), whence its > > > spectrum is contained in the unit circle; its eigenvalues are indeed > > > n-th roots of unity. > >
[R] making a log file with error messages?
Dear all: Is there anyone who would help me to generate a log file reporting error messages? I can see error messages on the screen, but I would like to report them as a file. thanks Sincerely, cahn --- ChaeHyung Ahn (cahn) Ph.D. Department of Biostatistics, School of Public Health CB #7420, University of North Carolina Chapel Hill, NC 27599-7420 phone: 919-966-7284 home: http://cyworld.nate.com/cahn88 __ R-help@stat.math.ethz.ch 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] Reading in a dataset with uneven variable lengths
One possibility is to use read.fwf(). Andy > From: Owen, Jason > > Suppose I have a text file that I want to read into R like the > following: > > X Y > 649 699 > 657 891 > 714 632 > 849 727 > 721 597 > 791 868 > 874 652 > 405 978 > 733 > 549 > 790 > > This is a simple example -- I could have a huge file with many > columns of unequal lengths. > > What is the best way to do it? I can't see how a data frame > can be used. I checked the FAQ and did a web search on the topic > but I came up empty. > > Jason > > __ > R-help@stat.math.ethz.ch 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] Reading in a dataset with uneven variable lengths
Suppose I have a text file that I want to read into R like the following: X Y 649 699 657 891 714 632 849 727 721 597 791 868 874 652 405 978 733 549 790 This is a simple example -- I could have a huge file with many columns of unequal lengths. What is the best way to do it? I can't see how a data frame can be used. I checked the FAQ and did a web search on the topic but I came up empty. Jason __ R-help@stat.math.ethz.ch 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] opimization problem
The form of the problem looks like you are trying to do a mean-variance portfolio optimization. If that is the case, you should not be dealing with variance as a restriction, but as part of the objective function: max (r'*w - rho*w'*V*w) s.t. sum(w) == 1 where rho is a risk aversion parameter. You can solve this as a quadratic programming problem using either 1) solve.QP from the quadprog package; 2) portfolio.optim in package tseries see http://tolstoy.newcastle.edu.au/R/help/05/01/10505.html for details on how to use the two. On Sunday 01 May 2005 19:21, Gottfried Gruber wrote: hi, i want to execute the following opimization problem: max r*w s.t.: w*z=1 # sum of w is 1 r, w are [nx1] vectors, z is a [nx1] vector consisting of 1 so far so good, works fine with lp the problem arises with the additional restriction w' * V * w where V is a [nxn] matrix how can i include this restriction since w arises twice? __ R-help@stat.math.ethz.ch 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] formula in fixed-effects part of GLMM
> From: Douglas Bates > > Henric Nilsson wrote: > > On Mon, 2005-05-02, 17:24, Douglas Bates wrote: > > > > > >>weihong wrote: > >> > >>>Can GLMM take formula derived from another object? > >>> > >>>foo <- glm (OVEN ~ h + h2, poisson, dataset) > >>> > >>># ok > >>>bar <- GLMM (OVEN ~ h + h2, poisson, dataset, random = > list (yr = ~1)) > >>> > >>>#error > >>>bar <- GLMM (foo$formula, poisson, dataset, random = list > (yr = ~1)) > >>>#Error in foo$("formula" + yr + 1) : invalid subscript type > >> > >>That won't work without some tweaking of the GLMM function. In the > >>0.95-1 and later versions of the lme4 package the > capabilities of GLMM > >>have been folded in to the lmer function and you would need > to fit that > >>model as > >> > >> lmer(OVEN ~ h + h2 + (1|yr), dataset, poisson) > >> > >>anyway. > > > > > > I don't have access to the "White Book" right now (or R, > for that matter), > > but doesn't it say that something like the following works? > > > > fit.lm <- lm(y ~ x) > > fit.glm <- update(fit.lm, class = "glm") > > > > But this isn't implemented in R, right? If one can make a > wish, it would > > be really nice being able to (using weihong's example): > > > > bar <- update(foo, . ~ . + (1|yr), class = "lmer") > > > > > > //Henric > > It's not very convenient from the point of view of dispatch. The > dispatch rules are such that the update method for the lm or glm class > would need to be aware of the lmer class to be able to do > this. The lm > and glm classes predate lmer by a long time. In that case, wouldn't lmer also need to be in "stats" (or somewhere "stats" can import from)? Best, Andy __ R-help@stat.math.ethz.ch 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] formula in fixed-effects part of GLMM
Henric Nilsson wrote: > On Mon, 2005-05-02, 17:24, Douglas Bates wrote: > > >>weihong wrote: >> >>>Can GLMM take formula derived from another object? >>> >>>foo <- glm (OVEN ~ h + h2, poisson, dataset) >>> >>># ok >>>bar <- GLMM (OVEN ~ h + h2, poisson, dataset, random = list (yr = ~1)) >>> >>>#error >>>bar <- GLMM (foo$formula, poisson, dataset, random = list (yr = ~1)) >>>#Error in foo$("formula" + yr + 1) : invalid subscript type >> >>That won't work without some tweaking of the GLMM function. In the >>0.95-1 and later versions of the lme4 package the capabilities of GLMM >>have been folded in to the lmer function and you would need to fit that >>model as >> >> lmer(OVEN ~ h + h2 + (1|yr), dataset, poisson) >> >>anyway. > > > I don't have access to the "White Book" right now (or R, for that matter), > but doesn't it say that something like the following works? > > fit.lm <- lm(y ~ x) > fit.glm <- update(fit.lm, class = "glm") > > But this isn't implemented in R, right? If one can make a wish, it would > be really nice being able to (using weihong's example): > > bar <- update(foo, . ~ . + (1|yr), class = "lmer") > > > //Henric It's not very convenient from the point of view of dispatch. The dispatch rules are such that the update method for the lm or glm class would need to be aware of the lmer class to be able to do this. The lm and glm classes predate lmer by a long time. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] eigenvalues of a circulant matrix
By the way, I just noticed that eigen(X) returns eigenvectors, at least two of which are NaN's. Best wishes! --- "Huntsinger, Reid" <[EMAIL PROTECTED]> wrote: > When the matrix is symmetric and omega is not real, omega and its conjugate > (= inverse) give the same eigenvalue, so you have a 2-dimensional > eigenspace. R chooses a real basis of this, which is perfectly fine since > it's not looking for circulant structure. > > For example, > > > m > [,1] [,2] [,3] [,4] [,5] > [1,]12332 > [2,]21233 > [3,]32123 > [4,]33212 > [5,]23321 > > > eigen(m) > $values > [1] 11.00 -0.381966 -0.381966 -2.618034 -2.618034 > > $vectors > [,1] [,2] [,3] [,4] [,5] > [1,] 0.4472136 0.00 -0.6324555 0.6324555 0.00 > [2,] 0.4472136 0.371748 0.5116673 0.1954395 0.601501 > [3,] 0.4472136 -0.601501 -0.1954395 -0.5116673 0.371748 > [4,] 0.4472136 0.601501 -0.1954395 -0.5116673 -0.371748 > [5,] 0.4472136 -0.371748 0.5116673 0.1954395 -0.601501 > > and you can match these columns up with the "canonical" eigenvectors > exp(2*pi*1i*(0:4)*j/5) for j = 0,1,2,3,4. E.g., > > > Im(exp(2*pi*1i*(0:4)*3/5)) > [1] 0.000 -0.5877853 0.9510565 -0.9510565 0.5877853 > > which can be seen to be a scalar multiple of column 2. > > Reid Huntsinger > > Reid Huntsinger > > -Original Message- > From: [EMAIL PROTECTED] > [mailto:[EMAIL PROTECTED] On Behalf Of Huntsinger, Reid > Sent: Monday, May 02, 2005 10:43 AM > To: 'Globe Trotter'; Rolf Turner > Cc: r-help@stat.math.ethz.ch > Subject: RE: [R] eigenvalues of a circulant matrix > > > It's hard to argue against the fact that a real symmetric matrix has real > eigenvalues. The eigenvalues of the circulant matrix with first row v are > *polynomials* (not the roots of 1 themselves, unless as Rolf suggested you > start with a vector with all zeros except one 1) in the roots of 1, with > coefficients equal to the entries in v. This is the finite Fourier transform > of v, by the way, and takes real values when the coefficients are real and > symmetric, ie when the matrix is symmetric. > > Reid Huntsinger > > -Original Message- > From: [EMAIL PROTECTED] > [mailto:[EMAIL PROTECTED] On Behalf Of Globe Trotter > Sent: Monday, May 02, 2005 10:23 AM > To: Rolf Turner > Cc: r-help@stat.math.ethz.ch > Subject: Re: [R] eigenvalues of a circulant matrix > > > > --- Rolf Turner <[EMAIL PROTECTED]> wrote: > > I just Googled around a bit and found definitions of Toeplitz and > > circulant matrices as follows: > > > > A Toeplitz matrix is any n x n matrix with values constant along each > > (top-left to lower-right) diagonal. matrix has the form > > > > a_0 a_1 . . . . ... a_{n-1} > > a_{-1} a_0 a_1... a_{n-2} > > a_{-2} a_{-1} a_0 a_1 .... > >. .. . . . > >. .. . . . > >. .. . . . > > a_{-(n-1)} a_{-(n-2)} ... a_1 a_0 > > > > (A Toeplitz matrix ***may*** be symmetric.) > > Agreed. As may a circulant matrix if a_i = a_{p-i+2} > > > > > A circulant matrix is an n x n matrix whose rows are composed of > > cyclically shifted versions of a length-n vector. For example, the > > circulant matrix on the vector (1, 2, 3, 4) is > > > > 4 1 2 3 > > 3 4 1 2 > > 2 3 4 1 > > 1 2 3 4 > > > > So circulant matrices are a special case of Toeplitz matrices. > > However a circulant matrix cannot be symmetric. > > > > The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i, > > 2 - 2i, and 2 --- certainly not roots of unity. > > The eigenvalues are 4+1*omega+2*omega^2+3*omega^3. > omega=cos(2*pi*k/4)+isin(2*pi*k/4) as k ranges over 1, 2, 3, 4, so the above > holds. > > Bellman may have > > been talking about the particular (important) case of a circulant > > matrix where the vector from which it is constructed is a canonical > > basis vector e_i with a 1 in the i-th slot and zeroes elsewhere. > > No, that is not true: his result can be verified for any circulant matrix, > directly. > > > Such a matrix is in fact a unitary matrix (operator), whence its > > spectrum is contained in the unit circle; its eigenvalues are indeed > > n-th roots of unity. > > > > Such matrices are related to the unilateral shift operator on > > Hilbert space (which is the ``primordial'' Toeplitz operator). > > It arises as multiplication by z on H^2 --- the ``analytic'' > > elements of L^2 of the unit circle. > > > > On (infinite dimensional) Hilbert space the unilateral shift > > looks like > > > > 0 0 0 0 0 ... > > 1 0 0 0 0 ... > > 0 1 0 0 0 ... > > 0 0 1 0 0 ... > > . . . . . ... > > . . . . . ... > > > > which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ... on and on > > forever. On (say) 4 dimensional space we can have a unilateral > > shift operator/matrix > > > > 0 0 0 0 > >
RE: [R] congratulations to the JGR developers
Well put, Gabor! P.S. Sorry for wasting the bandwidth. Just couldn't resist, so much true it is. > -Original Message- > From: [EMAIL PROTECTED] > [mailto:[EMAIL PROTECTED] On Behalf Of Gabor > Grothendieck > Sent: Monday, May 02, 2005 8:06 AM > To: Andy Bunn > Cc: R-Help > Subject: Re: [R] congratulations to the JGR developers > > On 5/2/05, Andy Bunn <[EMAIL PROTECTED]> wrote: > > > Just want to offer my congratulations to the JGR > developers as the > > > recepient of the 2005 Chambers Award. Great job, guys!! > > > http://stats.math.uni-augsburg.de/JGR/ > > > > This feels like the future of R to me. It's simple, powerful, and > > elegant just like R. As soon as the binary that works with 2.1 is > > released I'll use it exclusively on Linux and Windows. I'm > deeply impressed. > > I have not tried JGR but regarding your three adjective > describing R, R is very powerful but I am not sure I would > characterize it as simple and elegant -- complex and > practical seem nearer to the mark to me. > Some parts of R may be simple and elegant but when I think of > simple and elegant languages I think of ones that are > organized around a > single concept like APL (arrays), Smalltalk (objects), etc. > The underlying > Lisp roots of R may have a certain simplicity to them and S3 > (though not > S4) is relatively simple but not R as a whole. On the other > hand, the fact that it is practical, powerful and free with a > broad set of builtin and addon functionality and has become a > de facto standard for statistical research have been > sufficient reason for me to do all my computing using R. > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! > http://www.R-project.org/posting-guide.html > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] User-defined random variable
On 5/1/05, Matthias Kohl <[EMAIL PROTECTED]> wrote: > >>>I would like to know whether it is possible with R to define a > >>>discrete random variable different from the ones already defined > >>>inside R and generate random numbers from that user-defined > >>>distribution. > >>> > >>Yes. One generic way is to specify the quantile function (as in > >>qpois() etc.) and do qfun(runif(N)). > > > >If the support discrete but also finite, you can also use sample(), e.g. > > > > sample(myset, N, replace = TRUE, prob = myprob) > > one can also use our R package "distr" to generate discrete random > variables. The subsequent code provides a function which generates an > object of class "DiscreteDistribution" based on a finite support "supp". > If "prob" is missing all elements in "supp" are equally weighted. Thanks you for all your helpful replies to my question. Paul __ R-help@stat.math.ethz.ch 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] formula in fixed-effects part of GLMM
On Mon, 2005-05-02, 17:24, Douglas Bates wrote: > weihong wrote: >> Can GLMM take formula derived from another object? >> >> foo <- glm (OVEN ~ h + h2, poisson, dataset) >> >> # ok >> bar <- GLMM (OVEN ~ h + h2, poisson, dataset, random = list (yr = ~1)) >> >> #error >> bar <- GLMM (foo$formula, poisson, dataset, random = list (yr = ~1)) >> #Error in foo$("formula" + yr + 1) : invalid subscript type > > That won't work without some tweaking of the GLMM function. In the > 0.95-1 and later versions of the lme4 package the capabilities of GLMM > have been folded in to the lmer function and you would need to fit that > model as > > lmer(OVEN ~ h + h2 + (1|yr), dataset, poisson) > > anyway. I don't have access to the "White Book" right now (or R, for that matter), but doesn't it say that something like the following works? fit.lm <- lm(y ~ x) fit.glm <- update(fit.lm, class = "glm") But this isn't implemented in R, right? If one can make a wish, it would be really nice being able to (using weihong's example): bar <- update(foo, . ~ . + (1|yr), class = "lmer") //Henric __ R-help@stat.math.ethz.ch 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] opimization problem
Hi Arne, i'm sorry - the additional restriction is w' * V * w <= j where j is a [1x1] matrix. i don't know how to incorporate this restriction in the lp-object. would be grateful on any response, thanks, gg On Monday 02 May 2005 09:55, Arne Henningsen wrote: > Hi Gottfried, > > w' * V * w is not a restriction, because there is no equal sign. > Do you mean w' * V * w = 1? > > Arne > > On Sunday 01 May 2005 19:21, Gottfried Gruber wrote: > > hi, > > > > i want to execute the following opimization problem: > > max r*w > > s.t.: w*z=1 # sum of w is 1 > > r, w are [nx1] vectors, z is a [nx1] vector consisting of 1 > > so far so good, works fine with lp > > > > the problem arises with the additional restriction > > w' * V * w > > where V is a [nxn] matrix > > how can i include this restriction since w arises twice? > > > > thanks, > > gg -- --- Gottfried Gruber mailto:[EMAIL PROTECTED] www: http://gogo.sehrsupa.net __ R-help@stat.math.ethz.ch 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] xinch/yinch equivalent for log axis
I would like to draw a vertical line from a given point, in user coordinates, to x inches before from another point, also in user coordinates. This is easy enough to do for linear scales, using code based on xinch/yinch, but I do not know how to do this for logarithmic scales. This code shows an example of what I mean[1]: split.screen(c(2,1)) screen(1) # Linear scale, works fine plot(1:100, cex=0.5, pch=19) plotheight <- diff(par("usr")[3:4]) igap = 0.5 # intended gap, in inches gap <- igap/par("pin")[2]*plotheight # gap in user units lines(c(20,20),c(100,60+gap)) lines(c(20,20),c(60-gap,1)) points(20,60, pch=18) screen(2) # Logarithmic scale, point no longer centered plot(1:100, cex=0.5, pch=19, log="y") plotheight <- diff(par("usr")[3:4]) igap = 10 gap <- igap/par("pin")[2]*plotheight lines(c(20,20),c(100,10+gap)) lines(c(20,20),c(10-gap,1)) points(20,10, pch=18) close.screen(all=TRUE) The top graph is as I want. The diamond is centered in the gap, and the gap is 1 inch high (2*igap). In bottom graph, using a log scale, the diamond is no longer centered in the gap and there is a non-linear relationship between the gap height in inches and igap. I understand why this is happening, and this is why the xinch and yinch functions raise a warning, but is there a way to handle it? For example, is there a function which will return the user coordinate of a point, x inches above a given point p (in user coordinates), for a logarithmic y axis? Thanks in advance, Steven Murdoch. [1] I actually want to use this in a larger script, where I leave a gap in the axis where the median in. This needs to be small but legible to the eye, so that is why I am defining it in inches. The source code is at: http://www.cl.cam.ac.uk/users/sjm217/projects/graphics/fancyaxis.R -- w: http://www.cl.cam.ac.uk/users/sjm217/ __ R-help@stat.math.ethz.ch 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] Restricted cubic spline function ERROR?: glm(Y~rcs(x,5))
Jan Verbesselt wrote: Dear all, Is the restricted cubic spline function working properly in the glm model? We used glm(y~rcs(x,5), family=binomial) but it seems that for some theoretical reasons the rcs, restricted cubic spline function can not be fitted by a glm function. Is this correct? Regards, Jan ((Originally, we used lrm(y~ rcs(x,5)) but we couldn't find how to derive the AIC value of the fitted model. Is there a solution?)) The fit object from lrm has all the information you need to easily compute AIC. Type names(fit object) to see what's there. If you really needed glm (which you don't), Design wants you to use glmD. -- Frank E Harrell Jr Professor and Chair School of Medicine Department of Biostatistics Vanderbilt University __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Re fixed -- R2.1.0: X11 font at size 14 could not be loaded
Thanks again to Prof. Ripley, Peter Dalgaard, and Seth Falcon for their prompt responses to my X11 font problem. After restarting the font server, the problem has been fixed. Xiang-Jun __ R-help@stat.math.ethz.ch 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] RMySQL query: why result takes so much memory in R ?
Hi, IMHO you need only when your columns are numeric X rows /100/100/8 MB. >>(1200*3)/100/100/8 [1] 450 But one of your columns is group char. I'm suffering in the past in lot of things with massive data and R and recognize doing how many as possible in the database, or you have to upgrade your computer to 2-4GB like a database machine!? regards, christian Christoph Lehmann schrieb: Hi I just started with RMySQL. I have a database with roughly 12 millions rows/records and 8 columns/fields. From all 12 millions of records I want to import 3 fields only. The fields are specified as:id int(11), group char(15), measurement float(4,2). Why does this take > 1G RAM? I run R on suse linux, with 1G RAM and with the code below it even fills the whole 1G of swap. I just don't understand how 12e6 * 3 can fill such a huge range of RAM? Thanks for clarification and potential solutions. ## my code library(RMySQL) drv <- dbDriver("MySQL") ch <- dbConnect(drv,dbname="testdb", user="root",password="mysql") testdb <- dbGetQuery(ch, "select id, group, measurement from mydata") dbDisconnect(ch) dbUnloadDriver(drv) ## end of my code Cheers Christoph __ R-help@stat.math.ethz.ch 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] ANOSIM
Hello, I am interested in using ANOSIM (analysis of similarities) which is part of the vegan package. When I look at the programmers notes they mention that "They don't quite trust this method, and someone should study it closely." I was just wondering if this had been done, and whether it had been found reliable or not? Thanks, Kerry Reeves __ R-help@stat.math.ethz.ch 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 ifelse matrix question
try something like this: dat <- matrix(rnorm(20*205), 20, 205) sds <- sd(t(dat)) ## dat[, 6:205][abs(dat[, 6:205]) > sds] <- NA Best, Dimitris Dimitris Rizopoulos Ph.D. Student Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven, Belgium Tel: +32/16/336899 Fax: +32/16/337015 Web: http://www.med.kuleuven.ac.be/biostat/ http://www.student.kuleuven.ac.be/~m0390867/dimitris.htm - Original Message - From: "Christine Krisky" <[EMAIL PROTECTED]> To: Sent: Monday, May 02, 2005 5:13 PM Subject: [R] newbie ifelse matrix question Hi all, I have time series data in a matrix format 20 rows x 205 columns and have been trying to replace outliers with NA. My first column contains the outliers threshold (3 Standard deviations) for each ow - here's a bit of the first row sd3 V1 V2 V3 V4 V5V6V7V8 V9 1 13.03267 1797157 75 84 58 -1.958649 0.048775 2.056198 8.063622 3.071045 What I want is a statment that says if the absolute value of [,5:205] <= column 1 keep the value, otherwise replace with NA I've tried this without success - I do know this line works with s-plus data1[,6:205]<-ifelse(abs(data1[,6:205])<=data1[,1],data1[,6:205], NA) But what I get with R is only partially correct. I get NA in places where there shouldn't be. About every 7 columns or so I end up with 3 columns of NA and no replacement with NA where there should be (V11) in other places. The majority of the matrix though is correct. sd3 V1 V2 V3 V4 V5V6V7V8 V9 V10 V11 V12 V13 V14 V15 1 13.03267 1797157 75 84 58 -1.958649 0.048775 2.056198 8.063622 3.071045 0.078468 -21.9141NA NA NA 3.115585 I've tried searching the manual and user list and countless changes to the syntax over the last week, but no luck so far. Here is my syntax and error data1<-read.table("testR_data",na.string="0.00") dim(data1) [1] 20 204 # Find sd of rows, then multiply * 3 sd_data3<-apply(data1[,5:204],1,sd,na.rm=TRUE)*3 # Attach sd_data3 to data1 sd_and_data1<-cbind(sd_data3,data1) # Replace values >= 3 SDs with NA sd_and_data1[,6:205]<-ifelse(abs(sd_and_data1[,6:205])<=sd_and_data1[,1],sd_and_data1[,6:205], NA) Warning message: provided 4000 variables to replace 200 variables in: "[<-.data.frame"(`*tmp*`, , 6:205, value = list(c(-1.958649, I can't figure this out and would be very grateful for your help. -chris __ R-help@stat.math.ethz.ch 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] eigenvalues of a circulant matrix
On 02-May-05 Rolf Turner wrote: > I just Googled around a bit and found definitions of Toeplitz and > circulant matrices as follows: > [...] > A circulant matrix is an n x n matrix whose rows are composed of > cyclically shifted versions of a length-n vector. For example, the > circulant matrix on the vector (1, 2, 3, 4) is > > 4 1 2 3 > 3 4 1 2 > 2 3 4 1 > 1 2 3 4 > > So circulant matrices are a special case of Toeplitz matrices. > However a circulant matrix cannot be symmetric. I suspect the confusion may lie in what's meant by "cyclically shifted". In Rolf's example above, each row is shifted right by 1 and the one that falls off the end is put at the beginning. This cannot be symmetric for general values in the fist row. However, if you shift left instead, then you get 4 1 2 3 1 2 3 4 2 3 4 1 3 4 1 2 and this *is* symmetric (and indeed will always be so, for general values in the first row). All the formal definitions of "circulant" which I have seen use Rolf's definition ("shift right"). Suppose we call the left-shifted one "anti-circulant" ("AC"). The vector (which *is* and eigenvector of a circulant matrix): c(1, w, w^2, ... , w^(p-1)) where w is *any* complex p-th root of unity (including 1), is not in general an eigenvector of a pxp AC matrix. Example: M<-matrix(c(1,2,3,2,3,1,3,1,2),nrow=3) M [,1] [,2] [,3] [1,]123 [2,]231 [3,]312 p1 <- 1 ; p2 <- (-1-sqrt(3))/2 ; p3 <- (-1+sqrt(3))/2 e1 <- c(1,p1,p1^2) e2 <- c(1,p2,p2^2) e3 <- c(1,p3,p3^2) v1<-(M%*%e1)[1] ; v1 [1] 6 cbind(round(M%*%e1/v1,15), e1) e1 [1,] 1 1 [2,] 1 1 [3,] 1 1 v2<-(M%*%e2)[1] ; v2 [1] -0.2320508 cbind(round(M%*%e2/v2,15), e2) e2 [1,] 1.0+0.000i 1.0+0.000i [2,] -0.5+0.8660254i -0.5-0.8660254i [3,] -0.5-0.8660254i -0.5+0.8660254i v3<-(M%*%e3)[1] ; v3 [1] 2.133975 cbind(round(M%*%e3/v3,15), e3) e3 [1,] 1.0+0.000i 1.0+0.000i [2,] -0.5-0.8660254i -0.5+0.8660254i [3,] -0.5+0.8660254i -0.5-0.8660254i (I've out in rounding because of nasty little specks of "i" that keep dropping out, as in p2^3 [1] 1 - 6.432571e-16i ) So (except for e1) e2 and e3 are not eigenvectors of M (note the switching of signs in the imaginary parts of row 2 and in row 3). The AC matrix, being symmetric, heas real eigenvalues and real eigevectors, as can be easily verified using 'eigen'. Therefore I suspect that "Globe Trotter"'s "circulant matrix" was in fact an "anti-circulant" ("AC") matrix. However, from the information he gave I'm not clear how to verify that this is the case. Hoping this helps, Ted. E-Mail: (Ted Harding) <[EMAIL PROTECTED]> Fax-to-email: +44 (0)870 094 0861 Date: 02-May-05 Time: 16:32:27 -- XFMail -- __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] eigenvalues of a circulant matrix
When the matrix is symmetric and omega is not real, omega and its conjugate (= inverse) give the same eigenvalue, so you have a 2-dimensional eigenspace. R chooses a real basis of this, which is perfectly fine since it's not looking for circulant structure. For example, > m [,1] [,2] [,3] [,4] [,5] [1,]12332 [2,]21233 [3,]32123 [4,]33212 [5,]23321 > eigen(m) $values [1] 11.00 -0.381966 -0.381966 -2.618034 -2.618034 $vectors [,1] [,2] [,3] [,4] [,5] [1,] 0.4472136 0.00 -0.6324555 0.6324555 0.00 [2,] 0.4472136 0.371748 0.5116673 0.1954395 0.601501 [3,] 0.4472136 -0.601501 -0.1954395 -0.5116673 0.371748 [4,] 0.4472136 0.601501 -0.1954395 -0.5116673 -0.371748 [5,] 0.4472136 -0.371748 0.5116673 0.1954395 -0.601501 and you can match these columns up with the "canonical" eigenvectors exp(2*pi*1i*(0:4)*j/5) for j = 0,1,2,3,4. E.g., > Im(exp(2*pi*1i*(0:4)*3/5)) [1] 0.000 -0.5877853 0.9510565 -0.9510565 0.5877853 which can be seen to be a scalar multiple of column 2. Reid Huntsinger Reid Huntsinger -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Huntsinger, Reid Sent: Monday, May 02, 2005 10:43 AM To: 'Globe Trotter'; Rolf Turner Cc: r-help@stat.math.ethz.ch Subject: RE: [R] eigenvalues of a circulant matrix It's hard to argue against the fact that a real symmetric matrix has real eigenvalues. The eigenvalues of the circulant matrix with first row v are *polynomials* (not the roots of 1 themselves, unless as Rolf suggested you start with a vector with all zeros except one 1) in the roots of 1, with coefficients equal to the entries in v. This is the finite Fourier transform of v, by the way, and takes real values when the coefficients are real and symmetric, ie when the matrix is symmetric. Reid Huntsinger -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Globe Trotter Sent: Monday, May 02, 2005 10:23 AM To: Rolf Turner Cc: r-help@stat.math.ethz.ch Subject: Re: [R] eigenvalues of a circulant matrix --- Rolf Turner <[EMAIL PROTECTED]> wrote: > I just Googled around a bit and found definitions of Toeplitz and > circulant matrices as follows: > > A Toeplitz matrix is any n x n matrix with values constant along each > (top-left to lower-right) diagonal. matrix has the form > > a_0 a_1 . . . . ... a_{n-1} > a_{-1} a_0 a_1... a_{n-2} > a_{-2} a_{-1} a_0 a_1 .... > . .. . . . > . .. . . . > . .. . . . > a_{-(n-1)} a_{-(n-2)} ... a_1 a_0 > > (A Toeplitz matrix ***may*** be symmetric.) Agreed. As may a circulant matrix if a_i = a_{p-i+2} > > A circulant matrix is an n x n matrix whose rows are composed of > cyclically shifted versions of a length-n vector. For example, the > circulant matrix on the vector (1, 2, 3, 4) is > > 4 1 2 3 > 3 4 1 2 > 2 3 4 1 > 1 2 3 4 > > So circulant matrices are a special case of Toeplitz matrices. > However a circulant matrix cannot be symmetric. > > The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i, > 2 - 2i, and 2 --- certainly not roots of unity. The eigenvalues are 4+1*omega+2*omega^2+3*omega^3. omega=cos(2*pi*k/4)+isin(2*pi*k/4) as k ranges over 1, 2, 3, 4, so the above holds. Bellman may have > been talking about the particular (important) case of a circulant > matrix where the vector from which it is constructed is a canonical > basis vector e_i with a 1 in the i-th slot and zeroes elsewhere. No, that is not true: his result can be verified for any circulant matrix, directly. > Such a matrix is in fact a unitary matrix (operator), whence its > spectrum is contained in the unit circle; its eigenvalues are indeed > n-th roots of unity. > > Such matrices are related to the unilateral shift operator on > Hilbert space (which is the ``primordial'' Toeplitz operator). > It arises as multiplication by z on H^2 --- the ``analytic'' > elements of L^2 of the unit circle. > > On (infinite dimensional) Hilbert space the unilateral shift > looks like > > 0 0 0 0 0 ... > 1 0 0 0 0 ... > 0 1 0 0 0 ... > 0 0 1 0 0 ... > . . . . . ... > . . . . . ... > > which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ... on and on > forever. On (say) 4 dimensional space we can have a unilateral > shift operator/matrix > > 0 0 0 0 > 1 0 0 0 > 0 1 0 0 > 0 0 1 0 > > but its range is a 3 dimensional subspace (e_4 gets ``killed''). > > The ``corresponding'' circulant matrix is > > 0 0 0 1 > 1 0 0 0 > 0 1 0 0 > 0 0 1 0 > > which is an onto mapping --- e_4 gets sent back to e_1. > > I hope this clears up some of the confusion. > > cheers, > >
Re: [R] congratulations to the JGR developers
On 5/2/05, Andy Bunn <[EMAIL PROTECTED]> wrote: > > I have not tried JGR but regarding your three adjective describing R, > > R is very powerful but I am not sure I would characterize it as simple > > and elegant -- complex and practical seem nearer to the mark to me. > > I take umbrage (and not in the sense of affording shade). > > Take 'simple' to mean plain, or having few ornamentations. Then it sounds > like R. As for elegance, R is refined, tasteful, and beautiful. When I grow > up, I want to marry R. This seems to me to be like the Oscar Wilde sketch where they are criticizing the King behind his back and every time they get caught redefine their words to turn the insult into praise in order to please him. They keep coming up with more and more outrageous insults making it harder and harder to figure out how to redefine them as praise. __ R-help@stat.math.ethz.ch 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] formula in fixed-effects part of GLMM
weihong wrote: > Can GLMM take formula derived from another object? > > foo <- glm (OVEN ~ h + h2, poisson, dataset) > > # ok > bar <- GLMM (OVEN ~ h + h2, poisson, dataset, random = list (yr = ~1)) > > #error > bar <- GLMM (foo$formula, poisson, dataset, random = list (yr = ~1)) > #Error in foo$("formula" + yr + 1) : invalid subscript type That won't work without some tweaking of the GLMM function. In the 0.95-1 and later versions of the lme4 package the capabilities of GLMM have been folded in to the lmer function and you would need to fit that model as lmer(OVEN ~ h + h2 + (1|yr), dataset, poisson) anyway. __ R-help@stat.math.ethz.ch 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] congratulations to the JGR developers
> I have not tried JGR but regarding your three adjective describing R, > R is very powerful but I am not sure I would characterize it as simple > and elegant -- complex and practical seem nearer to the mark to me. I take umbrage (and not in the sense of affording shade). Take 'simple' to mean plain, or having few ornamentations. Then it sounds like R. As for elegance, R is refined, tasteful, and beautiful. When I grow up, I want to marry 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] "apply" question
On 5/2/05, Christoph Scherber <[EMAIL PROTECTED]> wrote: > Dear R users, > > I´ve got a simple question but somehow I can´t find the solution: > > I have a data frame with columns 1-5 containing one set of integer > values, and columns 6-10 containing another set of integer values. > Columns 6-10 contain NA´s at some places. > > I now want to calculate > (1) the number of values in each row of columns 6-10 that were NA´s Supposing our data is called DF, rowSums(!is.na(DF[,6:10])) > (2) the sum of all values on columns 1-5 for which there were no missing > values in the corresponding cells of columns 6-10. In the expression below 1 + 0 *DF[,6:10] is like DF[,6:10] except all non-NAs are replaced by 1. Multiplying DF[,1:5] by that effectively replaces each element in DF[,1:5] with an NA if the corresponding DF[,6:10] contained an NA. rowSums( DF[,1:5] * (1 + 0 * DF[,6:10]), na.rm = TRUE ) > > Example: (let´s call the data frame "data") > > Col1 Col2 Col3 Col4 Col5 Col6 Col7 Col8 Col9 Col10 > 1 2 5 2 3 NA 5 NA1 4 > 3 1 4 5 2 6 NA 4 NA 1 > > The result would then be (for the first row) > (1) "There were 2 NA´s in columns 6-10." > (2) The mean of Columns 1-5 was 2+2+3=7" (because there were NA´s in the > 1st and 3rd position in rows 6-10) I guess you meant sum when you referred to mean in (2). If you really do want the mean replace rowSums with rowMeans in the expression given above in the answer to (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] RMySQL query: why result takes so much memory in R ?
Christoph Lehmann wrote: Hi I just started with RMySQL. I have a database with roughly 12 millions rows/records and 8 columns/fields. From all 12 millions of records I want to import 3 fields only. The fields are specified as:id int(11), group char(15), measurement float(4,2). Why does this take > 1G RAM? I run R on suse linux, with 1G RAM and with the code below it even fills the whole 1G of swap. I just don't understand how 12e6 * 3 can fill such a huge range of RAM? Thanks for clarification and potential solutions. Those fields are each 8 or 20 bytes in size, so you're talking 12e6 times 36 or about nearly half a Gig for each copy. Presumably the code is storing more than one or two copies of the data. Why don't you use fetch() to get your records in more manageable chunks? Duncan Murdoch ## my code library(RMySQL) drv <- dbDriver("MySQL") ch <- dbConnect(drv,dbname="testdb", user="root",password="mysql") testdb <- dbGetQuery(ch, "select id, group, measurement from mydata") dbDisconnect(ch) dbUnloadDriver(drv) ## end of my code Cheers Christoph __ R-help@stat.math.ethz.ch 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 ifelse matrix question
Hi all, I have time series data in a matrix format 20 rows x 205 columns and have been trying to replace outliers with NA. My first column contains the outliers threshold (3 Standard deviations) for each row - here's a bit of the first row sd3 V1 V2 V3 V4 V5V6V7V8V9 1 13.03267 1797157 75 84 58 -1.958649 0.048775 2.056198 8.063622 3.071045 What I want is a statment that says if the absolute value of [,5:205] <= column 1 keep the value, otherwise replace with NA I've tried this without success - I do know this line works with s-plus data1[,6:205]<-ifelse(abs(data1[,6:205])<=data1[,1],data1[,6:205], NA) But what I get with R is only partially correct. I get NA in places where there shouldn't be. About every 7 columns or so I end up with 3 columns of NA and no replacement with NA where there should be (V11) in other places. The majority of the matrix though is correct. sd3 V1 V2 V3 V4 V5V6V7V8V9 V10 V11 V12 V13V14 V15 1 13.03267 1797157 75 84 58 -1.958649 0.048775 2.056198 8.063622 3.071045 0.078468 -21.9141NA NA NA 3.115585 I've tried searching the manual and user list and countless changes to the syntax over the last week, but no luck so far. Here is my syntax and error > data1<-read.table("testR_data",na.string="0.00") > dim(data1) [1] 20 204 # Find sd of rows, then multiply * 3 > sd_data3<-apply(data1[,5:204],1,sd,na.rm=TRUE)*3 # Attach sd_data3 to data1 > sd_and_data1<-cbind(sd_data3,data1) # Replace values >= 3 SDs with NA > sd_and_data1[,6:205]<-ifelse(abs(sd_and_data1[,6:205])<=sd_and_data1[,1],sd_and_data1[,6:205], > NA) Warning message: provided 4000 variables to replace 200 variables in: "[<-.data.frame"(`*tmp*`, , 6:205, value = list(c(-1.958649, I can't figure this out and would be very grateful for your help. -chris __ R-help@stat.math.ethz.ch 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] Roots of quadratic system.
Hi, I used polyroot() and it works fine. Sebastian On Monday 02 May 2005 02:51, John Janmaat wrote: > Hello Bill, > > I have used the optimization approach you suggest in past. I was hoping > that someone had written something specifically for solving a system of > nonlinear equations, as the fsolve function does in MatLab. The Octave > version is somewhat limited compared to the MatLab version, and I like > working in R. > > Thanks, > > John. > > ps: I would like the system to have a unique solution, but there is > nothing about the system that precludes multiple equilibria. Of course, > the L(x) = ... approach can search for multiple equilibria if I try > enough different starting points. > > [EMAIL PROTECTED] wrote: > > Are you looking for a unique solution or families of solutions? > > > > Can't you turn a root-finding problem for a system of equations > > with a unique solution into an optimisation problem, anyway? > > > > E.g. You want to solve > > > > f1(x) = g1 > > f2(x) = g2 > > ... > > > > Why not optimise L(x) = (f1(x) - g1)^2 + (f2(x) - g2)^2 + ... > > with respect to x? If the minimum value is zero, then you are > > done; if it is greater than zero your original system does not > > have a solution. > > > > If you are in the complex domain the changes needed are obvious. > > > > V. > > > > : -Original Message- > > : From: [EMAIL PROTECTED] > > : [mailto:[EMAIL PROTECTED] On Behalf Of John Janmaat > > : Sent: Monday, 2 May 2005 12:48 AM > > : To: r-help@stat.math.ethz.ch > > : Subject: [R] Roots of quadratic system. > > : > > : > > : Hello, > > : > > : I have a system of quadratic equations (results of a > > : Hamiltonian optimization) > > : which I need to find the roots for. Is there a package > > : and/or function which > > : will find the roots for a quadratic system? Note that I am > > : not opimizing, but > > : rather solving the first order conditions which come from a > > : Hamiltonian. I am > > : basically looking for something in R that will do the same > > : thing as fsolve in > > : Matlab. > > : > > : Thanks, > > : > > : John. > > : > > : == > > : Dr. John Janmaat > > : Department of Economics > > : Acadia University > > : Tel: 902-585-1461 > > : > > : __ > > : R-help@stat.math.ethz.ch mailing list > > : https://stat.ethz.ch/mailman/listinfo/r-help > > : PLEASE do read the posting guide! > > : http://www.R-project.org/posting-guide.html -- Sebastian Leuzinger Institute of Botany, University of Basel Schönbeinstr. 6 CH-4056 Basel ph0041 (0) 61 2673511 fax 0041 (0) 61 2673504 email [EMAIL PROTECTED] web http://pages.unibas.ch/botschoen/leuzinger __ R-help@stat.math.ethz.ch 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] "apply" question
you could try something like this: dat <- rbind(c(1, 2, 5, 2, 3, NA, 5, NA, 1, 4), c(3, 1, 4, 5, 2, 6, NA, 4, NA, 1)) ## # (1) rowSums(is.na(dat[, 6:10])) ## (2) dat. <- dat[, 1:5] dat.[is.na(dat[, 6:10])] <- NA rowSums(dat., na.rm=TRUE) rowMeans(dat., na.rm=TRUE) 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/16/336899 Fax: +32/16/337015 Web: http://www.med.kuleuven.ac.be/biostat/ http://www.student.kuleuven.ac.be/~m0390867/dimitris.htm - Original Message - From: "Christoph Scherber" <[EMAIL PROTECTED]> To: Sent: Monday, May 02, 2005 4:52 PM Subject: [R] "apply" question Dear R users, I´ve got a simple question but somehow I can´t find the solution: I have a data frame with columns 1-5 containing one set of integer values, and columns 6-10 containing another set of integer values. Columns 6-10 contain NA´s at some places. I now want to calculate (1) the number of values in each row of columns 6-10 that were NA´s (2) the sum of all values on columns 1-5 for which there were no missing values in the corresponding cells of columns 6-10. Example: (let´s call the data frame "data") Col1 Col2 Col3 Col4 Col5 Col6 Col7 Col8 Col9 Col10 1 2 5 2 3 NA 5 NA1 4 3 1 4 5 2 6 NA 4 NA 1 The result would then be (for the first row) (1) "There were 2 NA´s in columns 6-10." (2) The mean of Columns 1-5 was 2+2+3=7" (because there were NA´s in the 1st and 3rd position in rows 6-10) So far, I know how to calculate the rowSums for the data.frame, but I don´t know how to condition these on the values of columns 6-10 rowSums(data[,1:5]) #that´s straightforward apply(data[,6:19],1,function(x)sum(is.na(x))) #this also works fine But I don´t know how to select just the desired values of columns 1-5 (as described above) Can anyone help me? Thanks a lot in advance! Best regards Christoph __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Re: simulate zero-truncated Poisson distribution
Dear Ted and Peter, Thank you so much for your answers. These solutions you suggested are really helpful. Best regards, Galina [[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] "apply" question
Try: > ## Number of NAs in columns 6-10. > colSums(is.na(data[6:10])) Col6 Col7 Col8 Col9 Col10 1 1 1 1 0 > > ## Number of NAs in each row of columns 6-10. > rowSums(is.na(data[6:10])) 1 2 2 2 > > ## Sums of rows 1-5 omitting corresponding NAs in cols 6-10. > rowSums(data[,1:5] * !is.na(data[,6:10])) 1 2 7 9 If all entries are numeric, it'd be easier to use matrices instead of data frames. HTH, Andy > From: Christoph Scherber > > Dear R users, > > I´ve got a simple question but somehow I can´t find the solution: > > I have a data frame with columns 1-5 containing one set of integer > values, and columns 6-10 containing another set of integer values. > Columns 6-10 contain NA´s at some places. > > I now want to calculate > (1) the number of values in each row of columns 6-10 that were NA´s > (2) the sum of all values on columns 1-5 for which there were > no missing > values in the corresponding cells of columns 6-10. > > > Example: (let´s call the data frame "data") > > Col1 Col2 Col3 Col4 Col5 Col6 Col7 Col8 Col9 Col10 > 1 2 5 2 3 NA 5 NA1 4 > 3 1 4 5 2 6 NA 4 NA 1 > > The result would then be (for the first row) > (1) "There were 2 NA´s in columns 6-10." > (2) The mean of Columns 1-5 was 2+2+3=7" (because there were > NA´s in the > 1st and 3rd position in rows 6-10) > > So far, I know how to calculate the rowSums for the data.frame, but I > don´t know how to condition these on the values of columns 6-10 > > rowSums(data[,1:5]) #that´s straightforward > apply(data[,6:19],1,function(x)sum(is.na(x))) #this also works fine > > But I don´t know how to select just the desired values of columns 1-5 > (as described above) > > > Can anyone help me? Thanks a lot in advance! > > Best regards > Christoph > > __ > R-help@stat.math.ethz.ch 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] eigenvalues of a circulant matrix
Please, please, we are talking about eigenvectors in the original post, not eigenvalues. They are indeed all real when they are symmetric (even for circulant matrices because the imaginary parts sum up to zero). Many thanks! --- "Huntsinger, Reid" <[EMAIL PROTECTED]> wrote: > It's hard to argue against the fact that a real symmetric matrix has real > eigenvalues. The eigenvalues of the circulant matrix with first row v are > *polynomials* (not the roots of 1 themselves, unless as Rolf suggested you > start with a vector with all zeros except one 1) in the roots of 1, with > coefficients equal to the entries in v. This is the finite Fourier transform > of v, by the way, and takes real values when the coefficients are real and > symmetric, ie when the matrix is symmetric. > > Reid Huntsinger > > -Original Message- > From: [EMAIL PROTECTED] > [mailto:[EMAIL PROTECTED] On Behalf Of Globe Trotter > Sent: Monday, May 02, 2005 10:23 AM > To: Rolf Turner > Cc: r-help@stat.math.ethz.ch > Subject: Re: [R] eigenvalues of a circulant matrix > > > > --- Rolf Turner <[EMAIL PROTECTED]> wrote: > > I just Googled around a bit and found definitions of Toeplitz and > > circulant matrices as follows: > > > > A Toeplitz matrix is any n x n matrix with values constant along each > > (top-left to lower-right) diagonal. matrix has the form > > > > a_0 a_1 . . . . ... a_{n-1} > > a_{-1} a_0 a_1... a_{n-2} > > a_{-2} a_{-1} a_0 a_1 .... > >. .. . . . > >. .. . . . > >. .. . . . > > a_{-(n-1)} a_{-(n-2)} ... a_1 a_0 > > > > (A Toeplitz matrix ***may*** be symmetric.) > > Agreed. As may a circulant matrix if a_i = a_{p-i+2} > > > > > A circulant matrix is an n x n matrix whose rows are composed of > > cyclically shifted versions of a length-n vector. For example, the > > circulant matrix on the vector (1, 2, 3, 4) is > > > > 4 1 2 3 > > 3 4 1 2 > > 2 3 4 1 > > 1 2 3 4 > > > > So circulant matrices are a special case of Toeplitz matrices. > > However a circulant matrix cannot be symmetric. > > > > The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i, > > 2 - 2i, and 2 --- certainly not roots of unity. > > The eigenvalues are 4+1*omega+2*omega^2+3*omega^3. > omega=cos(2*pi*k/4)+isin(2*pi*k/4) as k ranges over 1, 2, 3, 4, so the above > holds. > > Bellman may have > > been talking about the particular (important) case of a circulant > > matrix where the vector from which it is constructed is a canonical > > basis vector e_i with a 1 in the i-th slot and zeroes elsewhere. > > No, that is not true: his result can be verified for any circulant matrix, > directly. > > > Such a matrix is in fact a unitary matrix (operator), whence its > > spectrum is contained in the unit circle; its eigenvalues are indeed > > n-th roots of unity. > > > > Such matrices are related to the unilateral shift operator on > > Hilbert space (which is the ``primordial'' Toeplitz operator). > > It arises as multiplication by z on H^2 --- the ``analytic'' > > elements of L^2 of the unit circle. > > > > On (infinite dimensional) Hilbert space the unilateral shift > > looks like > > > > 0 0 0 0 0 ... > > 1 0 0 0 0 ... > > 0 1 0 0 0 ... > > 0 0 1 0 0 ... > > . . . . . ... > > . . . . . ... > > > > which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ... on and on > > forever. On (say) 4 dimensional space we can have a unilateral > > shift operator/matrix > > > > 0 0 0 0 > > 1 0 0 0 > > 0 1 0 0 > > 0 0 1 0 > > > > but its range is a 3 dimensional subspace (e_4 gets ``killed''). > > > > The ``corresponding'' circulant matrix is > > > > 0 0 0 1 > > 1 0 0 0 > > 0 1 0 0 > > 0 0 1 0 > > > > which is an onto mapping --- e_4 gets sent back to e_1. > > > > I hope this clears up some of the confusion. > > > > cheers, > > > > Rolf Turner > > [EMAIL PROTECTED] > > Many thanks and best wishes! > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! > http://www.R-project.org/posting-guide.html > > > > > > -- > Notice: This e-mail message, together with any attachment...{{dropped}} __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Variance Homogenity
Working with living models sometimes the response is reflected in the variability more than the middle parameters, i.e., the variability can reflected the system response, p.ej., in toxicology models the variability in the response means something about the selective pressure, buy how to evaluated this variability?. The ANOVA Test evaluated difference between the medians using the variance between groups and behind groups, but i can´t use ANOVA if the two suposes don´t exists, normalitiy and variance homogeniety, then how to evaluated the differences in variance?, Barttlet Test evaluated the homoscedascity using the square sum and not the variance as a response variable. How could I evaluated the difference in groups with differents variances?. I don´t want use Kruskall-Wallis Test because my variable is continue and numerical. Do someone know what do must do it? José Herrera Bazán Science Faculty UNAM http://search.t1msn.com.mx/ __ R-help@stat.math.ethz.ch 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] congratulations to the JGR developers
On 5/2/05, Andy Bunn <[EMAIL PROTECTED]> wrote: > > Just want to offer my congratulations to the JGR developers as > > the recepient > > of the 2005 Chambers Award. Great job, guys!! > > http://stats.math.uni-augsburg.de/JGR/ > > This feels like the future of R to me. It's simple, powerful, and elegant > just like R. As soon as the binary that works with 2.1 is released I'll use > it exclusively on Linux and Windows. I'm deeply impressed. I have not tried JGR but regarding your three adjective describing R, R is very powerful but I am not sure I would characterize it as simple and elegant -- complex and practical seem nearer to the mark to me. Some parts of R may be simple and elegant but when I think of simple and elegant languages I think of ones that are organized around a single concept like APL (arrays), Smalltalk (objects), etc. The underlying Lisp roots of R may have a certain simplicity to them and S3 (though not S4) is relatively simple but not R as a whole. On the other hand, the fact that it is practical, powerful and free with a broad set of builtin and addon functionality and has become a de facto standard for statistical research have been sufficient reason for me to do all my computing using R. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] RMySQL query: why result takes so much memory in R ?
Hi I just started with RMySQL. I have a database with roughly 12 millions rows/records and 8 columns/fields. From all 12 millions of records I want to import 3 fields only. The fields are specified as:id int(11), group char(15), measurement float(4,2). Why does this take > 1G RAM? I run R on suse linux, with 1G RAM and with the code below it even fills the whole 1G of swap. I just don't understand how 12e6 * 3 can fill such a huge range of RAM? Thanks for clarification and potential solutions. ## my code library(RMySQL) drv <- dbDriver("MySQL") ch <- dbConnect(drv,dbname="testdb", user="root",password="mysql") testdb <- dbGetQuery(ch, "select id, group, measurement from mydata") dbDisconnect(ch) dbUnloadDriver(drv) ## end of my code Cheers Christoph __ R-help@stat.math.ethz.ch 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] "apply" question
- Original Message - From: "Christoph Scherber" <[EMAIL PROTECTED]> To: Sent: Monday, May 02, 2005 10:52 AM Subject: [R] "apply" question Dear R users, I´ve got a simple question but somehow I can´t find the solution: I have a data frame with columns 1-5 containing one set of integer values, and columns 6-10 containing another set of integer values. Columns 6-10 contain NA´s at some places. I now want to calculate (1) the number of values in each row of columns 6-10 that were NA´s (2) the sum of all values on columns 1-5 for which there were no missing values in the corresponding cells of columns 6-10. Example: (let´s call the data frame "data") Col1 Col2 Col3 Col4 Col5 Col6 Col7 Col8 Col9 Col10 1 2 5 2 3 NA 5 NA1 4 3 1 4 5 2 6 NA 4 NA 1 The result would then be (for the first row) (1) "There were 2 NA´s in columns 6-10." (2) The mean of Columns 1-5 was 2+2+3=7" (because there were NA´s in the 1st and 3rd position in rows 6-10) So far, I know how to calculate the rowSums for the data.frame, but I don´t know how to condition these on the values of columns 6-10 rowSums(data[,1:5]) #that´s straightforward apply(data[,6:19],1,function(x)sum(is.na(x))) #this also works fine But I don´t know how to select just the desired values of columns 1-5 (as described above) tmp <- rowSums(data[apply(data[,6:19],1,function(x) sum(is.na(x)))==0,1:5]) Now, tmp contains only the rowsums for the rows with no NAs in the other columns. Sean __ R-help@stat.math.ethz.ch 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] "apply" question
Dear R users, I´ve got a simple question but somehow I can´t find the solution: I have a data frame with columns 1-5 containing one set of integer values, and columns 6-10 containing another set of integer values. Columns 6-10 contain NA´s at some places. I now want to calculate (1) the number of values in each row of columns 6-10 that were NA´s (2) the sum of all values on columns 1-5 for which there were no missing values in the corresponding cells of columns 6-10. Example: (let´s call the data frame "data") Col1 Col2 Col3 Col4 Col5 Col6 Col7 Col8 Col9 Col10 1 2 5 2 3 NA 5 NA1 4 3 1 4 5 2 6 NA 4 NA 1 The result would then be (for the first row) (1) "There were 2 NA´s in columns 6-10." (2) The mean of Columns 1-5 was 2+2+3=7" (because there were NA´s in the 1st and 3rd position in rows 6-10) So far, I know how to calculate the rowSums for the data.frame, but I don´t know how to condition these on the values of columns 6-10 rowSums(data[,1:5]) #that´s straightforward apply(data[,6:19],1,function(x)sum(is.na(x))) #this also works fine But I don´t know how to select just the desired values of columns 1-5 (as described above) Can anyone help me? Thanks a lot in advance! Best regards Christoph __ R-help@stat.math.ethz.ch 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] eigenvalues of a circulant matrix
It's hard to argue against the fact that a real symmetric matrix has real eigenvalues. The eigenvalues of the circulant matrix with first row v are *polynomials* (not the roots of 1 themselves, unless as Rolf suggested you start with a vector with all zeros except one 1) in the roots of 1, with coefficients equal to the entries in v. This is the finite Fourier transform of v, by the way, and takes real values when the coefficients are real and symmetric, ie when the matrix is symmetric. Reid Huntsinger -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Globe Trotter Sent: Monday, May 02, 2005 10:23 AM To: Rolf Turner Cc: r-help@stat.math.ethz.ch Subject: Re: [R] eigenvalues of a circulant matrix --- Rolf Turner <[EMAIL PROTECTED]> wrote: > I just Googled around a bit and found definitions of Toeplitz and > circulant matrices as follows: > > A Toeplitz matrix is any n x n matrix with values constant along each > (top-left to lower-right) diagonal. matrix has the form > > a_0 a_1 . . . . ... a_{n-1} > a_{-1} a_0 a_1... a_{n-2} > a_{-2} a_{-1} a_0 a_1 .... > . .. . . . > . .. . . . > . .. . . . > a_{-(n-1)} a_{-(n-2)} ... a_1 a_0 > > (A Toeplitz matrix ***may*** be symmetric.) Agreed. As may a circulant matrix if a_i = a_{p-i+2} > > A circulant matrix is an n x n matrix whose rows are composed of > cyclically shifted versions of a length-n vector. For example, the > circulant matrix on the vector (1, 2, 3, 4) is > > 4 1 2 3 > 3 4 1 2 > 2 3 4 1 > 1 2 3 4 > > So circulant matrices are a special case of Toeplitz matrices. > However a circulant matrix cannot be symmetric. > > The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i, > 2 - 2i, and 2 --- certainly not roots of unity. The eigenvalues are 4+1*omega+2*omega^2+3*omega^3. omega=cos(2*pi*k/4)+isin(2*pi*k/4) as k ranges over 1, 2, 3, 4, so the above holds. Bellman may have > been talking about the particular (important) case of a circulant > matrix where the vector from which it is constructed is a canonical > basis vector e_i with a 1 in the i-th slot and zeroes elsewhere. No, that is not true: his result can be verified for any circulant matrix, directly. > Such a matrix is in fact a unitary matrix (operator), whence its > spectrum is contained in the unit circle; its eigenvalues are indeed > n-th roots of unity. > > Such matrices are related to the unilateral shift operator on > Hilbert space (which is the ``primordial'' Toeplitz operator). > It arises as multiplication by z on H^2 --- the ``analytic'' > elements of L^2 of the unit circle. > > On (infinite dimensional) Hilbert space the unilateral shift > looks like > > 0 0 0 0 0 ... > 1 0 0 0 0 ... > 0 1 0 0 0 ... > 0 0 1 0 0 ... > . . . . . ... > . . . . . ... > > which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ... on and on > forever. On (say) 4 dimensional space we can have a unilateral > shift operator/matrix > > 0 0 0 0 > 1 0 0 0 > 0 1 0 0 > 0 0 1 0 > > but its range is a 3 dimensional subspace (e_4 gets ``killed''). > > The ``corresponding'' circulant matrix is > > 0 0 0 1 > 1 0 0 0 > 0 1 0 0 > 0 0 1 0 > > which is an onto mapping --- e_4 gets sent back to e_1. > > I hope this clears up some of the confusion. > > cheers, > > Rolf Turner > [EMAIL PROTECTED] Many thanks and best wishes! __ R-help@stat.math.ethz.ch 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] congratulations to the JGR developers
> Just want to offer my congratulations to the JGR developers as > the recepient > of the 2005 Chambers Award. Great job, guys!! > http://stats.math.uni-augsburg.de/JGR/ This feels like the future of R to me. It's simple, powerful, and elegant just like R. As soon as the binary that works with 2.1 is released I'll use it exclusively on Linux and Windows. I'm deeply impressed. -Andy __ R-help@stat.math.ethz.ch 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] eigenvalues of a circulant matrix
--- Rolf Turner <[EMAIL PROTECTED]> wrote: > I just Googled around a bit and found definitions of Toeplitz and > circulant matrices as follows: > > A Toeplitz matrix is any n x n matrix with values constant along each > (top-left to lower-right) diagonal. matrix has the form > > a_0 a_1 . . . . ... a_{n-1} > a_{-1} a_0 a_1... a_{n-2} > a_{-2} a_{-1} a_0 a_1 .... > . .. . . . > . .. . . . > . .. . . . > a_{-(n-1)} a_{-(n-2)} ... a_1 a_0 > > (A Toeplitz matrix ***may*** be symmetric.) Agreed. As may a circulant matrix if a_i = a_{p-i+2} > > A circulant matrix is an n x n matrix whose rows are composed of > cyclically shifted versions of a length-n vector. For example, the > circulant matrix on the vector (1, 2, 3, 4) is > > 4 1 2 3 > 3 4 1 2 > 2 3 4 1 > 1 2 3 4 > > So circulant matrices are a special case of Toeplitz matrices. > However a circulant matrix cannot be symmetric. > > The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i, > 2 - 2i, and 2 --- certainly not roots of unity. The eigenvalues are 4+1*omega+2*omega^2+3*omega^3. omega=cos(2*pi*k/4)+isin(2*pi*k/4) as k ranges over 1, 2, 3, 4, so the above holds. Bellman may have > been talking about the particular (important) case of a circulant > matrix where the vector from which it is constructed is a canonical > basis vector e_i with a 1 in the i-th slot and zeroes elsewhere. No, that is not true: his result can be verified for any circulant matrix, directly. > Such a matrix is in fact a unitary matrix (operator), whence its > spectrum is contained in the unit circle; its eigenvalues are indeed > n-th roots of unity. > > Such matrices are related to the unilateral shift operator on > Hilbert space (which is the ``primordial'' Toeplitz operator). > It arises as multiplication by z on H^2 --- the ``analytic'' > elements of L^2 of the unit circle. > > On (infinite dimensional) Hilbert space the unilateral shift > looks like > > 0 0 0 0 0 ... > 1 0 0 0 0 ... > 0 1 0 0 0 ... > 0 0 1 0 0 ... > . . . . . ... > . . . . . ... > > which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ... on and on > forever. On (say) 4 dimensional space we can have a unilateral > shift operator/matrix > > 0 0 0 0 > 1 0 0 0 > 0 1 0 0 > 0 0 1 0 > > but its range is a 3 dimensional subspace (e_4 gets ``killed''). > > The ``corresponding'' circulant matrix is > > 0 0 0 1 > 1 0 0 0 > 0 1 0 0 > 0 0 1 0 > > which is an onto mapping --- e_4 gets sent back to e_1. > > I hope this clears up some of the confusion. > > cheers, > > Rolf Turner > [EMAIL PROTECTED] Many thanks and best wishes! __ R-help@stat.math.ethz.ch 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] eigenvalues of a circulant matrix
I just Googled around a bit and found definitions of Toeplitz and circulant matrices as follows: A Toeplitz matrix is any n x n matrix with values constant along each (top-left to lower-right) diagonal. matrix has the form a_0 a_1 . . . . ... a_{n-1} a_{-1} a_0 a_1... a_{n-2} a_{-2} a_{-1} a_0 a_1 .... . .. . . . . .. . . . . .. . . . a_{-(n-1)} a_{-(n-2)} ... a_1 a_0 (A Toeplitz matrix ***may*** be symmetric.) A circulant matrix is an n x n matrix whose rows are composed of cyclically shifted versions of a length-n vector. For example, the circulant matrix on the vector (1, 2, 3, 4) is 4 1 2 3 3 4 1 2 2 3 4 1 1 2 3 4 So circulant matrices are a special case of Toeplitz matrices. However a circulant matrix cannot be symmetric. The eigenvalues of the forgoing circulant matrix are 10, 2 + 2i, 2 - 2i, and 2 --- certainly not roots of unity. Bellman may have been talking about the particular (important) case of a circulant matrix where the vector from which it is constructed is a canonical basis vector e_i with a 1 in the i-th slot and zeroes elsewhere. Such a matrix is in fact a unitary matrix (operator), whence its spectrum is contained in the unit circle; its eigenvalues are indeed n-th roots of unity. Such matrices are related to the unilateral shift operator on Hilbert space (which is the ``primordial'' Toeplitz operator). It arises as multiplication by z on H^2 --- the ``analytic'' elements of L^2 of the unit circle. On (infinite dimensional) Hilbert space the unilateral shift looks like 0 0 0 0 0 ... 1 0 0 0 0 ... 0 1 0 0 0 ... 0 0 1 0 0 ... . . . . . ... . . . . . ... which maps e_0 to e_1, e_1 to e_2, e_2 to e_3, ... on and on forever. On (say) 4 dimensional space we can have a unilateral shift operator/matrix 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 0 but its range is a 3 dimensional subspace (e_4 gets ``killed''). The ``corresponding'' circulant matrix is 0 0 0 1 1 0 0 0 0 1 0 0 0 0 1 0 which is an onto mapping --- e_4 gets sent back to e_1. I hope this clears up some of the confusion. 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] eigenvalues of a circulant matrix
Again, what's your "kinv"? On Mon, 2 May 2005, Globe Trotter wrote: > OK, lets redo this again, and ensure that we start with a row that will indeed > lead to a symmetric matrix for the circulant matrix: > > x<-scan("kinv") > y<-x[c(109:1,2:108)] > > X=toeplitz(y) > Z=y > for (i in 2:216) Z=rbind(Z,y[c((216-i+2):216,1:(216-i+1))]) > > range(X-Z) > [1] 0 0 > > eigen(X) is the same as eigen(Z), but we know that Z is a circulant matrix so > the eigenvectors are complex > > Any thoughts/screams? > > > --- Prof Brian Ripley <[EMAIL PROTECTED]> wrote: > > On Sun, 1 May 2005, someone who didn't give his name wrote: > > > > > It is my understanding that the eigenvectors of a circulant matrix are > > > given as follows: > > > > > > 1,omega,omega^2,,omega^{p-1} > > > > > > where the matrix has dimension given by p x p and omega is one of p > > > complex > > > roots of unity. (See Bellman for an excellent discussion on this). > > > > What is the relevance of this? Also, your reference is useless to us, > > which is important as this all hinges on your definitions. > > > > > The matrix created by the attached row and obtained using the following > > > commands indicates no imaginary parts for the eigenvectors. It appears > > > that the real values are close, but not exactly so, and there is no > > > imaginary part whatsoever. > > > > > > x<-scan("kinv.dat") #length(x) = 216 > > > y<-x[c(109:216,1:108)] > > > X<-toeplitz(y) > > > eigen(X)$vectors > > > > We don't have "kinv.dat", but X is not circulant as usually defined. > > > > > Note that the eigenvectors are correct, and they are indeed real, > > > because X is symmetric. > > > > > > Is this a bug in R? Any insight if not, please! > > > > Well, first R calls LAPACK or EISPACK, so it would be a bug in one of > > those. But in so far as I understand you, X is a real symmetric matrix, > > and those have real eigenvalues and eigenvectors. > > > > I think you are confused about the meaning of Toeplitz and circulant. > > Compare > > > > http://mathworld.wolfram.com/CirculantMatrix.html > > http://mathworld.wolfram.com/ToeplitzMatrix.html > > > > and note that ?toeplitz says it computes the *symmetric* Toeplitz matrix. > > > > There is a very regretable tendency here for people to assume their > > lack of understanding is `a bug in R'. > > > > -- > > 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-help@stat.math.ethz.ch 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] eigenvalues of a circulant matrix
OK, lets redo this again, and ensure that we start with a row that will indeed lead to a symmetric matrix for the circulant matrix: x<-scan("kinv") y<-x[c(109:1,2:108)] X=toeplitz(y) Z=y for (i in 2:216) Z=rbind(Z,y[c((216-i+2):216,1:(216-i+1))]) range(X-Z) [1] 0 0 eigen(X) is the same as eigen(Z), but we know that Z is a circulant matrix so the eigenvectors are complex Any thoughts/screams? --- Prof Brian Ripley <[EMAIL PROTECTED]> wrote: > On Sun, 1 May 2005, someone who didn't give his name wrote: > > > It is my understanding that the eigenvectors of a circulant matrix are > > given as follows: > > > > 1,omega,omega^2,,omega^{p-1} > > > > where the matrix has dimension given by p x p and omega is one of p complex > > roots of unity. (See Bellman for an excellent discussion on this). > > What is the relevance of this? Also, your reference is useless to us, > which is important as this all hinges on your definitions. > > > The matrix created by the attached row and obtained using the following > > commands indicates no imaginary parts for the eigenvectors. It appears > > that the real values are close, but not exactly so, and there is no > > imaginary part whatsoever. > > > > x<-scan("kinv.dat") #length(x) = 216 > > y<-x[c(109:216,1:108)] > > X<-toeplitz(y) > > eigen(X)$vectors > > We don't have "kinv.dat", but X is not circulant as usually defined. > > > Note that the eigenvectors are correct, and they are indeed real, > > because X is symmetric. > > > > Is this a bug in R? Any insight if not, please! > > Well, first R calls LAPACK or EISPACK, so it would be a bug in one of > those. But in so far as I understand you, X is a real symmetric matrix, > and those have real eigenvalues and eigenvectors. > > I think you are confused about the meaning of Toeplitz and circulant. > Compare > > http://mathworld.wolfram.com/CirculantMatrix.html > http://mathworld.wolfram.com/ToeplitzMatrix.html > > and note that ?toeplitz says it computes the *symmetric* Toeplitz matrix. > > There is a very regretable tendency here for people to assume their > lack of understanding is `a bug in R'. > > -- > 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] Multivariate kernel density estimation
locfit() in the `locfit' package should be able to handle at least 5 dimensions, but you need tons of data and either lots of time or a very fast computer. ssden() in the `gss' package will handle up to 4 dimensions. Those are the only two I know about. Others may well know better. HTH, Andy > From: Stephan Tolksdorf > > Hi, > > I need to estimate the density at the mean of a sample of a few > thousands data points with a dimesion up to 5. The data is > uni-modal and > regularly shaped. > > I couldn't find any kernel density package for R which supports more > than 3 dimensions. Have I overlooked a package or does somebody have > code for this purpose? Any other advice? > > Regards, >Stephan > > __ > R-help@stat.math.ethz.ch 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] eigenvalues of a circulant matrix
Dear Professor Ripley: Lets do this professionally, shall we? --- Prof Brian Ripley <[EMAIL PROTECTED]> wrote: > On Sun, 1 May 2005, someone who didn't give his name wrote: > > > It is my understanding that the eigenvectors of a circulant matrix are > > given as follows: > > > > 1,omega,omega^2,,omega^{p-1} > > > > where the matrix has dimension given by p x p and omega is one of p complex > > roots of unity. (See Bellman for an excellent discussion on this). > > What is the relevance of this? Also, your reference is useless to us, > which is important as this all hinges on your definitions. Bellman is an excellent book on the topic and that is what I was alluding to, that you can calculate the eigendecomposition by hand without any costly computations, actually. > > > The matrix created by the attached row and obtained using the following > > commands indicates no imaginary parts for the eigenvectors. It appears > > that the real values are close, but not exactly so, and there is no > > imaginary part whatsoever. > > > > x<-scan("kinv.dat") #length(x) = 216 > > y<-x[c(109:216,1:108)] > > X<-toeplitz(y) > > eigen(X)$vectors > > We don't have "kinv.dat", but X is not circulant as usually defined. Sorry about kinv.dat -- in the e-mail that came back to me, it read "kinv", btw. I am unclear why you say that "X is not circulant as usually defined" -- do you think you could clarify? It is true I use a Toeplitz matrix to set this up, but how does that matter? The end result in this case is still a circulant matrix that is symmetric, is it not? I would like to know why my result is not circulant here. > > Note that the eigenvectors are correct, and they are indeed real, > > because X is symmetric. > > > > Is this a bug in R? Any insight if not, please! > > Well, first R calls LAPACK or EISPACK, so it would be a bug in one of > those. But in so far as I understand you, X is a real symmetric matrix, > and those have real eigenvalues and eigenvectors. Yes, I know that R calls LAPACK (which now contains EISPACK, btw). But I also know that LAPACK contains complex eigendecomposition routines in addition to double precision ones and it would need to be used if there is reason to believe that the result is complex valued. (In particular ZGESDD would do it.) The eigendecomposition of a matrix is unique. Whatever you think of Bellman, the book does show how the eigenvectors of a circulant matrix are given by the complex roots of unity as given above. We have therefore exhibited an eigendecomposition without actually going through major computations (which is good, because statistical computing is best when you use it sparingly). Why then does the result differ from that in R, and why by so much? (After all, the eigendecomposition is unique, or is that only fpr real matrices?) > I think you are confused about the meaning of Toeplitz and circulant. Unclear, but would like to hear about your views on the actual differences in this specific example. > Compare > > http://mathworld.wolfram.com/CirculantMatrix.html > http://mathworld.wolfram.com/ToeplitzMatrix.html > > and note that ?toeplitz says it computes the *symmetric* Toeplitz matrix. In my case, my matrix is symmetric and the result is a circulant matrix. > There is a very regretable tendency here for people to assume their > lack of understanding is `a bug in R'. True, but bugs in software are not exactly rare. Though R does have very few bugs and which is why I recommend the software to every Tom. Besides I asked a question here because I was confused Many thanks and best wishes! > -- > 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] Take each one cell
see ?sapply and ?for Sean - Original Message - From: "Jan Sabee" <[EMAIL PROTECTED]> To: Sent: Monday, May 02, 2005 8:23 AM Subject: [R] Take each one cell Are there any way to take x <- c("0", "large", "medium", "small") x [1] "0" "large" "medium" "small" like x="0" x="large" x="medium" x="small" Best regards Jan Sabee __ R-help@stat.math.ethz.ch 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] Multivariate kernel density estimation
Hi, I need to estimate the density at the mean of a sample of a few thousands data points with a dimesion up to 5. The data is uni-modal and regularly shaped. I couldn't find any kernel density package for R which supports more than 3 dimensions. Have I overlooked a package or does somebody have code for this purpose? Any other advice? Regards, Stephan __ R-help@stat.math.ethz.ch 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] Take each one cell
Are there any way to take > x <- c("0", "large", "medium", "small") > x [1] "0" "large" "medium" "small" > like x="0" x="large" x="medium" x="small" Best regards Jan Sabee __ R-help@stat.math.ethz.ch 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] Roots of quadratic system.
Hello Bill, I have used the optimization approach you suggest in past. I was hoping that someone had written something specifically for solving a system of nonlinear equations, as the fsolve function does in MatLab. The Octave version is somewhat limited compared to the MatLab version, and I like working in R. Thanks, John. ps: I would like the system to have a unique solution, but there is nothing about the system that precludes multiple equilibria. Of course, the L(x) = ... approach can search for multiple equilibria if I try enough different starting points. [EMAIL PROTECTED] wrote: Are you looking for a unique solution or families of solutions? Can't you turn a root-finding problem for a system of equations with a unique solution into an optimisation problem, anyway? E.g. You want to solve f1(x) = g1 f2(x) = g2 ... Why not optimise L(x) = (f1(x) - g1)^2 + (f2(x) - g2)^2 + ... with respect to x? If the minimum value is zero, then you are done; if it is greater than zero your original system does not have a solution. If you are in the complex domain the changes needed are obvious. V. : -Original Message- : From: [EMAIL PROTECTED] : [mailto:[EMAIL PROTECTED] On Behalf Of John Janmaat : Sent: Monday, 2 May 2005 12:48 AM : To: r-help@stat.math.ethz.ch : Subject: [R] Roots of quadratic system. : : : Hello, : : I have a system of quadratic equations (results of a : Hamiltonian optimization) : which I need to find the roots for. Is there a package : and/or function which : will find the roots for a quadratic system? Note that I am : not opimizing, but : rather solving the first order conditions which come from a : Hamiltonian. I am : basically looking for something in R that will do the same : thing as fsolve in : Matlab. : : Thanks, : : John. : : == : Dr. John Janmaat : Department of Economics : Acadia University : Tel: 902-585-1461 : : __ : R-help@stat.math.ethz.ch mailing list : https://stat.ethz.ch/mailman/listinfo/r-help : PLEASE do read the posting guide! : http://www.R-project.org/posting-guide.html : -- == Dr. John Janmaat Department of Economics Acadia University Wolfville, Nova Scotia, Canada. __ R-help@stat.math.ethz.ch 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] Roots of quadratic system.
Kjetil Brinchmann Halvorsen wrote: John Janmaat wrote: Hello, I have a system of quadratic equations (results of a Hamiltonian optimization) which I need to find the roots for. Is there a package and/or function which will find the roots for a quadratic system? Certainly you cxould use solve, see ?solve Alternatively you could go for a computer algebra system with an implemantation of groebner basis, and use an symbolic method. I have looked into using solve. However, solve works on a system of linear equations (at least that is how I read it). I have a system of quadratic equations. They can be written to equate to zero, so that a non-linear system solver should do the trick. John. Kjetil Note that I am not opimizing, but rather solving the first order conditions which come from a Hamiltonian. I am basically looking for something in R that will do the same thing as fsolve in Matlab. Thanks, John. == Dr. John Janmaat Department of Economics Acadia University Tel: 902-585-1461 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- == Dr. John Janmaat Department of Economics Acadia University Wolfville, Nova Scotia, Canada. __ R-help@stat.math.ethz.ch 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] Keep R-mirror updated
Hi, I have a R mirror on www.termix.ufv.br/CRAN I try to keep it updated using rsync and cron. rsync -rzuav --progress --delete cran.r-project.org::CRAN /var/www/termix/mirrors/linux/CRAN/ & But all times the rsync abort with error: 95.99kB/s0:00:02 548960 77% 110.50kB/s0:00:01rsync: read error: Connection reset by peer (104) rsync error: error in rsync protocol data stream (code 12) at io.c(515) rsync: connection unexpectedly closed (391811 bytes received so far) [generator] rsync error: error in rsync protocol data stream (code 12) at io.c(359) And it dont made my mirror updated. Exist any method more efficient to made it updated? Thanks Ronaldo -- O cachorro abana o rabo porque o rabo, enquanto segmento discriminado do corpo da sociedade, ainda não conquistou o direito de abanar o cachorro -- |> // | \\ [***] | ( õ õ ) [Ronaldo Reis Júnior] |> V [UFV/DBA-Entomologia] |/ \ [36571-000 Viçosa - MG ] |> /(.''`.)\ [Fone: 31-3899-4007 ] | /(: :' :)\ [EMAIL PROTECTED]] |>/ (`. `'` ) \[ICQ#: 5692561 | LinuxUser#: 205366 ] |( `- ) [***] |>> _/ \_Powered by GNU/Debian Woody/Sarge __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] [R-pkgs] RWinEdt_1.7-1 for R>=2.1.0
Dear RWinEdt users, unfortunately, it seems to be impossible to support localized (i.e. translated) versions of RGui in MDI mode. The support of this feature has been removed now. All the workarounds I tried caused bugs (as many of you told me re. RWinEdt_1.7-0). Thanks for all the bug reports! It is recommended to use the new version of RWinEdt (1.7-1) together with RGui in SDI mode (which should work with arbitrary translations). Further on, only english versions of RGui in MDI mode are supported by RWinEdt. Uwe Ligges ___ R-packages mailing list [EMAIL PROTECTED] https://stat.ethz.ch/mailman/listinfo/r-packages __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Multiple regression
None: basic R comes with a package `stats' that does both (very well!). It is loaded by default. You do need to read `An Introduction to R', as the posting guide asks. On Mon, 2 May 2005, alessandro carletti wrote: what package could I install to best perform a Multiple Linear Regression (and what for PCA)? 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] Multiple regression
Hi, what package could I install to best perform a Multiple Linear Regression (and what for PCA)? Thanks Alessandro Carletti __ R-help@stat.math.ethz.ch 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] Restricted cubic spline function ERROR?: glm(Y~rcs(x,5))
R itself does not have a function rcs(), and help.search("rcs") does not find one. My guess is that you are using package Design and failing to tell us. [Aside here: this is not the first time I have seen help.search() results that give no indication why the result was found. To wit Design.trans(Design)Design Special Transformation Functions rcspline.eval(Hmisc)Restricted Cubic Spline Design Matrix rcspline.plot(Hmisc)Plot Restricted Cubic Spline Function rcspline.restate(Hmisc) Re-state Restricted Cubic Spline Function GtkVisibility(RGtk) Automatically generated S functions for bindings to the RGtk library mysqlInitDriver(RMySQL) Support Functions sqliteInitDriver(RSQLite) Support Functions prob.frcs.dat(verification) Probablisitic Forecast Dataset. is not much help and suggested to me to look at the second to fourth entries.] If you read the help page ?Design.trans you will see that rcs is part of the Design system built on top of S/R. What makes you think it is intended to work with glm()? On Mon, 2 May 2005, Jan Verbesselt wrote: Is the restricted cubic spline function working properly in the glm model? We used glm(y~rcs(x,5), family=binomial) but it seems that for some theoretical reasons the rcs, restricted cubic spline function can not be fitted by a glm function. Is this correct? ((Originally, we used lrm(y~ rcs(x,5)) but we couldn't find how to derive the AIC value of the fitted model. Is there a solution?)) Yes. Make use of your theory to write an AIC() method. Note that lrm> fit <- lrm(y ~ blood.pressure + sex * (age + rcs(cholesterol, 4)), x = TRUE, y = TRUE) class(fit) [1] "lrm""Design" "glm" is apparently incorrectly asserting inheritance and so the current logLik and hence AIC methods selected do not work. (An "lrm" object is missing the "family" component that is documented for a "glm" object, and in any case I believe lrm also fits non-GLM models.) -- 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] Restricted cubic spline function ERROR?: glm(Y~rcs(x,5))
Dear all, Is the restricted cubic spline function working properly in the glm model? We used glm(y~rcs(x,5), family=binomial) but it seems that for some theoretical reasons the rcs, restricted cubic spline function can not be fitted by a glm function. Is this correct? Regards, Jan ((Originally, we used lrm(y~ rcs(x,5)) but we couldn't find how to derive the AIC value of the fitted model. Is there a solution?)) __ R-help@stat.math.ethz.ch 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] opimization problem
Hi Gottfried, w' * V * w is not a restriction, because there is no equal sign. Do you mean w' * V * w = 1? Arne On Sunday 01 May 2005 19:21, Gottfried Gruber wrote: > hi, > > i want to execute the following opimization problem: > max r*w > s.t.: w*z=1 # sum of w is 1 > r, w are [nx1] vectors, z is a [nx1] vector consisting of 1 > so far so good, works fine with lp > > the problem arises with the additional restriction > w' * V * w > where V is a [nxn] matrix > how can i include this restriction since w arises twice? > > thanks, > gg -- Arne Henningsen Department of Agricultural Economics University of Kiel Olshausenstr. 40 D-24098 Kiel (Germany) Tel: +49-431-880 4445 Fax: +49-431-880 1397 [EMAIL PROTECTED] http://www.uni-kiel.de/agrarpol/ahenningsen/ __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html