[R] nls() and numerical integration (e.g. integrate()) working together?
Dear List-Members, since 3 weeks I have been heavily working on reproducing the results of an economic paper. The method there uses the numerical solution of an integral within nonlinear least squares. Within the integrand there is also some parameter to estimate. Is that in the end possible to implement in R [Originally it was done in GAUSS]? I'm nearly into giving up. I constucted an example to showing the problems I face. I have three questions - related to three errors shown below: 1) How to make clear that in the integrand z is the integration variable and b1 is a parameter while x1 is a data variable 2) and 3) How to set up a correct estimation of the integral? library(stats) y <- c(2,15,24,21,5,6,) x1 <- c(2.21,5,3,5,2,1) x2 <- c(4.51,6,2,11,0.4,3) f <- function(z) {z + b1*x1} vf <- Vectorize(f) g <- function(z) {z + x1} vg <- Vectorize(f) Error 1: > nls(y ~ integrate(vf,0,1)+b2*x2,start=list(b1=0.5,b2=2)) Error in function (z) : object "b1" not found Error 2: > nls(y ~ integrate(vg,0,1)+b2*x2,start=list(b1=0.5,b2=2)) Error in integrate(vg, 0, 1) : REAL() can only be applied to a 'numeric', not a 'list' Error 3: > nls(y ~ integrate(g,0,1)+b2*x2,start=list(b1=0.5,b2=2)) Error in integrate(g, 0, 1) + b2 * x2 : non-numeric argument to binary operator In addition: Warning messages: 1: longer object length is not a multiple of shorter object length in: z + x1 With a lot of thanks in advance, Michael __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] two labels on x-axis (year and month)
hej i'm plotting time-series and label the x-axis as follows: r <- as.POSIXct(round(range(p1$time), "month")) to define the time range for labeling the xaxis plot(p1$time,p1$ p1, type="l", xaxt="n") plots p1 against time axis.POSIXct(1, at=seq(r[1], r[2], by="month"), format="%m") labels the axis in mothly steps. what I want do do now, is a second label for the x-axis, that stands lower and indicates the years. like 05 06 07 08 09 10 11 12 01 02 03 04 05 etc... 2004 2005 I don't know how to proceed with all the possibilities in axis, par and plot thanks for your help __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] gWidgets (tcltk): problem extracting values f rom widgets in glayout grid
Dear Evan, epamail.epa.gov> writes: > > > Hello, > > I haven't been able to find an example for the second case below -- or > perhaps I didn't recognize it when I saw it. > Is there a value for x such that svalue(x) will return "bbb", either by > itself or as part of an array? Or do I need to do something else > entirely? > (R2.5.1; Windows XP) > You will need to store the widget. See below. > > ### Case 2: this makes an identical window & widget, with "b" > replacing "a" > > ### but does't return the value > > wtesta = gwindow("wtestb",visible=TRUE) > > wtestb = glayout(visible=TRUE,container=wtesta) > > wtestb[1,1]=gedit("bbb",container=wtestb) > > print(svalue(wtestb[1,1])) the glayout container doesn't work that way. You don't get a [ method, just a [<- one. Try something like this instead to store the widget: wtestb[1,1]=(bbb <- gedit("bbb",container=wtestb) ) print(svalue(bbb)) --John > Evan Englund > U.S. EPA > 702-798-2248 > > __ > R-help stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] degrees of freedom question
R2.3, WinXP Dear all, I am using the following functions: f1 = Phi1+(Phi2-Phi1)/(1+exp((log(Phi3)-log(x))/exp(log(Phi4))) f2 = Phi1+(Phi2-Phi1)/(1+exp((log(Phi3)-log(r)-log(x))/exp(log(Phi4))) subject to the residual weighting Var(e[i]) = sigma^2 * abs( E(y) )^(2*Delta) Here is my question, in steps: 1. Function f1 is separately fitted to two different datasets corresponding to two different dose response curves. These fits are unweighted. 2. Function f2 is fitted to the pooled data such that the two dose response curves are assumed to differ _only_ in log(r). This fit is also unweighted. 3. The residuals from #2 are used to estimate an appropriate sigma^2 and Delta to use in weighting. 4. The functions described in #1 and #2 are refitted, but this time weighted using the information gathered in #3. 5. How many degrees of freedom should be allocated to the weighted residual sums of squares? (There are three such SSE's, one for each individual model, and one for the overall joint model) Much thanks in advance, Greg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Need a variant of rbind for datasets with different numbers of columns
Where is the data coming from since it has a variable number of columns in each row? Is it coming from a text file? If so, you can use the "fill=TRUE" option when reading to fill out empty columns. You need to provide at least a subset of the data so we can see what you are working with. On 8/22/07, Kirsten Beyer <[EMAIL PROTECTED]> wrote: > Hello. I am looking for a function that will allow me to paste rows > together without regard for the numbers of columns in the datasets to > be joined. The only columns where it matters if they are aligned > correctly are at the beginning - the rest of the columns represent > differing numbers of ICD9 (disease) codes reported by each > person(record) at a health visit. They are in no particular order. > > For example, a result would look like this: > > patient ICD91 ICD92 ICD93 > patient A 12345 67891543 > patient B3469 9090 > patient C 1234 > > I am trying to accomplish this inside a loop which first identifies > the codes associated with the person and then joins them to the > person. I have the code working so that it can create a row for each > person, but I can't figure out how to join these rows together! FYI, > my dataset has 200,000+ people. > > 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 > and provide commented, minimal, self-contained, reproducible code. > -- Jim Holtman Cincinnati, OH +1 513 646 9390 What is the problem you are trying to solve? __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] gWidgets (tcltk): problem extracting values from widgets in glayout grid
Hello, I haven't been able to find an example for the second case below -- or perhaps I didn't recognize it when I saw it. Is there a value for x such that svalue(x) will return "bbb", either by itself or as part of an array? Or do I need to do something else entirely? (R2.5.1; Windows XP) > gWidgets test > options("guiToolkit"="tcltk") > require(gWidgets) [1] TRUE > > ## Case 1: this works > wtesta = gwindow("wtesta",visible=TRUE) > wtesta1 = gedit("aaa",container=wtesta) > print(svalue(wtesta1)) [1] "aaa" > > ### Case 2: this makes an identical window & widget, with "b" replacing "a" > ### but does't return the value > wtesta = gwindow("wtestb",visible=TRUE) > wtestb = glayout(visible=TRUE,container=wtesta) > wtestb[1,1]=gedit("bbb",container=wtestb) > print(svalue(wtestb[1,1])) Error in function (classes, fdef, mtable) : unable to find an inherited method for function ".leftBracket", for signature "gLayouttcltk", "guiWidgetsToolkittcltk" Error in svalue(wtestb[1, 1]) : error in evaluating the argument 'obj' in selecting a method for function 'svalue' > print(svalue(wtestb)) Error in function (classes, fdef, mtable) : unable to find an inherited method for function ".svalue", for signature "gLayouttcltk", "NULL", "NULL", "guiWidgetsToolkittcltk" Evan Englund U.S. EPA 702-798-2248 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] integrate
You can divide your domain of integration into smaller intervals and then add up the individual contributions. This could improve the speed of adaptive Gauss-Kronrod quadrature used in integrate(). Ravi. --- Ravi Varadhan, Ph.D. Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University Ph: (410) 502-2619 Fax: (410) 614-9625 Email: [EMAIL PROTECTED] Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Santanu Pramanik Sent: Wednesday, August 22, 2007 6:41 PM To: [EMAIL PROTECTED] Cc: r-help@stat.math.ethz.ch; [EMAIL PROTECTED] Subject: Re: [R] integrate Hi Andy, Thank your very much for your input. I also tried something like that which gives a value close to 20, basically using the same trapezoidal rule. > sum(apply(as.matrix(seq(-10,10,by=0.1)),1,my.fcn))*0.1 [1] 20.17385 Actually my function is much more complicated than the posted example and it is taking a long time... Anyway, thanks again, Santanu JPSM, 1218J Lefrak Hall University of Maryland, College Park Phone: 301-314-9916 >>> <[EMAIL PROTECTED]> 8/22/2007 5:20:05 PM >>> As Duncan Murdoch mentioned in his reply, the problem is with the fact that your function is not really a properly defined function in the sense of "assigning a unique y to each x". The integrate function uses an adaptive quadrature routine which probably makes multiple calls to the function being integrated and expects to get the same y's for the same x's every time. If you want to get a number close to 20 (for your example) you need an integration routine which will use single evaluation of your "function" at each value of x. A simple method like rectangular approximation on a grid or the so-called trapezoidal rule will do just that. Here is a very crude prototype of such an integrator: integrate1 <- function(f, lower, upper){ f <- match.fun(f) xx <- seq(lower, upper, length=100) del <- xx[2] - xx[1] yy <- f(xx[-100]) return(del*sum(yy)) Now when you run integrate1(my.fun, -10, 10) you will get a number close to 20 but, of course, every time you do it you will get a different value. Hope this helps, Andy __ Andy Jaworski 518-1-01 Process Laboratory 3M Corporate Research Laboratory - E-mail: [EMAIL PROTECTED] Tel: (651) 733-6092 Fax: (651) 736-3122 "Santanu Pramanik" <[EMAIL PROTECTED] To .umd.edu> Sent by: cc [EMAIL PROTECTED] at.math.ethz.ch Subject [R] integrate 08/22/2007 02:56 PM Hi, I am trying to integrate a function which is approximately constant over the range of the integration. The function is as follows: > my.fcn = function(mu){ + m = 1000 + z = 0 + z.mse = 0 + for(i in 1:m){ + z[i] = rnorm(1, mu, 1) + z.mse = z.mse + (z[i] - mu)^2 + } + return(z.mse/m) + } > my.fcn(-10) [1] 1.021711 > my.fcn(10) [1] 0.9995235 > my.fcn(-5) [1] 1.012727 > my.fcn(5) [1] 1.033595 > my.fcn(0) [1] 1.106282 > The function takes the value (approx) 1 over the range of mu. So the integration result should be close to 20 if we integrate over (-10, 10). But R gives: > integrate(my.fcn, -10, 10) 685.4941 with absolute error < 7.6e-12 > integrate(Vectorize(my.fcn), -10, 10) # this code never terminate I have seen in the "?integrate" page it is clearly written: If the function is approximately constant (in particular, zero) over nearly all its range it is possible that the result and error estimate may be seriously wrong. But this doesn't help in solving the problem. Thanks, Santanu JPSM, 1218J Lefrak Hall University of Maryland, College Park Phone: 301-314-9916 [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list
Re: [R] integrate
Hi Andy, Thank your very much for your input. I also tried something like that which gives a value close to 20, basically using the same trapezoidal rule. > sum(apply(as.matrix(seq(-10,10,by=0.1)),1,my.fcn))*0.1 [1] 20.17385 Actually my function is much more complicated than the posted example and it is taking a long time... Anyway, thanks again, Santanu JPSM, 1218J Lefrak Hall University of Maryland, College Park Phone: 301-314-9916 >>> <[EMAIL PROTECTED]> 8/22/2007 5:20:05 PM >>> As Duncan Murdoch mentioned in his reply, the problem is with the fact that your function is not really a properly defined function in the sense of "assigning a unique y to each x". The integrate function uses an adaptive quadrature routine which probably makes multiple calls to the function being integrated and expects to get the same y's for the same x's every time. If you want to get a number close to 20 (for your example) you need an integration routine which will use single evaluation of your "function" at each value of x. A simple method like rectangular approximation on a grid or the so-called trapezoidal rule will do just that. Here is a very crude prototype of such an integrator: integrate1 <- function(f, lower, upper){ f <- match.fun(f) xx <- seq(lower, upper, length=100) del <- xx[2] - xx[1] yy <- f(xx[-100]) return(del*sum(yy)) Now when you run integrate1(my.fun, -10, 10) you will get a number close to 20 but, of course, every time you do it you will get a different value. Hope this helps, Andy __ Andy Jaworski 518-1-01 Process Laboratory 3M Corporate Research Laboratory - E-mail: [EMAIL PROTECTED] Tel: (651) 733-6092 Fax: (651) 736-3122 "Santanu Pramanik" <[EMAIL PROTECTED] To .umd.edu> Sent by: cc [EMAIL PROTECTED] at.math.ethz.ch Subject [R] integrate 08/22/2007 02:56 PM Hi, I am trying to integrate a function which is approximately constant over the range of the integration. The function is as follows: > my.fcn = function(mu){ + m = 1000 + z = 0 + z.mse = 0 + for(i in 1:m){ + z[i] = rnorm(1, mu, 1) + z.mse = z.mse + (z[i] - mu)^2 + } + return(z.mse/m) + } > my.fcn(-10) [1] 1.021711 > my.fcn(10) [1] 0.9995235 > my.fcn(-5) [1] 1.012727 > my.fcn(5) [1] 1.033595 > my.fcn(0) [1] 1.106282 > The function takes the value (approx) 1 over the range of mu. So the integration result should be close to 20 if we integrate over (-10, 10). But R gives: > integrate(my.fcn, -10, 10) 685.4941 with absolute error < 7.6e-12 > integrate(Vectorize(my.fcn), -10, 10) # this code never terminate I have seen in the "?integrate" page it is clearly written: If the function is approximately constant (in particular, zero) over nearly all its range it is possible that the result and error estimate may be seriously wrong. But this doesn't help in solving the problem. Thanks, Santanu JPSM, 1218J Lefrak Hall University of Maryland, College Park Phone: 301-314-9916 [[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 and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] triviality solved: Unwanted but added additional first column when using write.csv2()
Hello! When I export a dataframe to a csv-file there is always an unwanted column added as the following demonstrates: > system("cat example.csv") ID;Name;Value1;Value2 100;Tom;11,90;32,90 101;Martha;15,90;49,00 102;Ute;20,00;300,00 > write.csv2(read.csv2 ("example.csv"), "examplecopy.csv") > system("cat examplecopy.csv") "";"ID";"Name";"Value1";"Value2" "1";100;"Tom";11,9;32,9 "2";101;"Martha";15,9;49 "3";102;"Ute";20;300 > I never want the column "" "1" "2" "3" to be written into a csv-file. I could find any hint neither in the manual nor in the archive, but then I found the solution: row.names = FALSE > write.csv2(read.csv2 ("example.csv"), "examplecopy.csv", row.names = FALSE) > system("cat examplecopy.csv") "ID";"Name";"Value1";"Value2" 100;"Tom";11,9;32,9 101;"Martha";15,9;49 102;"Ute";20;300 But (perhaps because I am not a native in English) it was more a solution by trial and chance than by reading and understanding the manual. I did not expect this unwanted column to be called "row.names" (the values are numbers which at last becomes after a subset() completely unusuable in most cases). So I did never search for row names but for "unwanted column write." I did not dare to ask before because of the triviality of the question (even I stumble over this problem for weeks) but I want the next beginner who searches this list archive for "unwanted column" to get a result. Greetings Delcour __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] integrate
As Duncan Murdoch mentioned in his reply, the problem is with the fact that your function is not really a properly defined function in the sense of "assigning a unique y to each x". The integrate function uses an adaptive quadrature routine which probably makes multiple calls to the function being integrated and expects to get the same y's for the same x's every time. If you want to get a number close to 20 (for your example) you need an integration routine which will use single evaluation of your "function" at each value of x. A simple method like rectangular approximation on a grid or the so-called trapezoidal rule will do just that. Here is a very crude prototype of such an integrator: integrate1 <- function(f, lower, upper){ f <- match.fun(f) xx <- seq(lower, upper, length=100) del <- xx[2] - xx[1] yy <- f(xx[-100]) return(del*sum(yy)) Now when you run integrate1(my.fun, -10, 10) you will get a number close to 20 but, of course, every time you do it you will get a different value. Hope this helps, Andy __ Andy Jaworski 518-1-01 Process Laboratory 3M Corporate Research Laboratory - E-mail: [EMAIL PROTECTED] Tel: (651) 733-6092 Fax: (651) 736-3122 "Santanu Pramanik" <[EMAIL PROTECTED] To .umd.edu> Sent by: cc [EMAIL PROTECTED] at.math.ethz.ch Subject [R] integrate 08/22/2007 02:56 PM Hi, I am trying to integrate a function which is approximately constant over the range of the integration. The function is as follows: > my.fcn = function(mu){ + m = 1000 + z = 0 + z.mse = 0 + for(i in 1:m){ + z[i] = rnorm(1, mu, 1) + z.mse = z.mse + (z[i] - mu)^2 + } + return(z.mse/m) + } > my.fcn(-10) [1] 1.021711 > my.fcn(10) [1] 0.9995235 > my.fcn(-5) [1] 1.012727 > my.fcn(5) [1] 1.033595 > my.fcn(0) [1] 1.106282 > The function takes the value (approx) 1 over the range of mu. So the integration result should be close to 20 if we integrate over (-10, 10). But R gives: > integrate(my.fcn, -10, 10) 685.4941 with absolute error < 7.6e-12 > integrate(Vectorize(my.fcn), -10, 10) # this code never terminate I have seen in the "?integrate" page it is clearly written: If the function is approximately constant (in particular, zero) over nearly all its range it is possible that the result and error estimate may be seriously wrong. But this doesn't help in solving the problem. Thanks, Santanu JPSM, 1218J Lefrak Hall University of Maryland, College Park Phone: 301-314-9916 [[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 and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Need a variant of rbind for datasets with different numbers of columns
Hello. I am looking for a function that will allow me to paste rows together without regard for the numbers of columns in the datasets to be joined. The only columns where it matters if they are aligned correctly are at the beginning - the rest of the columns represent differing numbers of ICD9 (disease) codes reported by each person(record) at a health visit. They are in no particular order. For example, a result would look like this: patient ICD91 ICD92 ICD93 patient A 12345 67891543 patient B3469 9090 patient C 1234 I am trying to accomplish this inside a loop which first identifies the codes associated with the person and then joins them to the person. I have the code working so that it can create a row for each person, but I can't figure out how to join these rows together! FYI, my dataset has 200,000+ people. 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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting lda results
Thanks for trying to answer even if my question isn't very clear. Please allow me to try again with an example: I want to predict group membership by the three variables below, and then plot the results in one figure that compares the frequency of scores of both groups on the LD obtained: Group<- c(1,1,1,1,1,1,1,2,2,2,2,2,2,2,2) var1<- c(2,3,1,2,3,2,3,3,5,6,7,6,8,5,5) var2<- c(9,9,9,8,7,8,9,3,2,2,1,1,2,3,3) var3<- c(6,7,6,6,5,6,7,1,2,1,2,3,1,1,2) data.df <- as.data.frame (cbind(Group, var1, var2, var3)) data.lda <- lda(Group~., data=data.df) predict(data.lda) plot(data.lda) As you said, I get two separate figures, one for each group, but I need them both in one figure. I tried ldahist by reading the help because it has an argument sep, but I'm not managing to get it going. I hope this code is reproducible but I'll certainly read the book suggested to make my questions easier to answer. Thanks! Silvia. Prof Brian Ripley wrote: > > Read ?plot.lda, which tells you the ... arguments are (for dimen=1, the > only option for two groups) passed to ldahist, so then read its help page. > > I tried what is suggested plot.lda but > > I don't know what you want (and your example is not reproducible): I would > expect you to get a single plot with two panels (figures), but there are > options to have a single panel. (Reading 'An Introduction to R' may help > you to use standard terminology that others will be able to follow.) > > On Wed, 22 Aug 2007, Silvia Lomascolo wrote: > >> >> Hi all, >> I am trying to plot the results of a discriminant analysis done with >> lda(MASS) but my groups appear in two different plots (in the same >> graphics >> device) and I want to combine them in one plot. My code looks like: >> >> BirdTrain.lda <- lda(Bdisperser~., data=BirdTrain.mx) >> predict(BirdTrain.lda) >> plot(BirdTrain.lda) >> >> I have two types of Bdisperser, so I only get one linear discriminant >> function. Can anyone please tell me how to combine the data in one plot? >> >> I work with R 2.4.1 using Windows. > > But the version of MASS is what is relevant, and it would have been in > the sessionInfo() output the R posting guide asked you for. > > -- > 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 > and provide commented, minimal, self-contained, reproducible code. > > -- View this message in context: http://www.nabble.com/plotting-lda-results-tf4312870.html#a12282028 Sent from the R help mailing list archive at Nabble.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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] distance by vegan
Hi, In vegan the function to calculate distance is vegdist in which the default distance is Sorensen (Bray-Curtis). So, try: vegdist(x, method="bray") Also see ?vegdist victor On 8/22/07, [EMAIL PROTECTED] <[EMAIL PROTECTED]> wrote: > > How to calculate sorensen (bray-curtis) distance by dist function > within the vegan package? > > Cheers > Duccio > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > -- Victor Lemes Landeiro Instituto Nacional de Pesquisas da Amazônia - INPA Manaus, Amazonas, Brasil Skypename: landeiro [[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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] integrate
On 8/22/2007 3:54 PM, Santanu Pramanik wrote: > Hi, > I am trying to integrate a function which is approximately constant > over the range of the integration. The function is as follows: That's not a function of the input mu. It includes a random component: > my.fcn(10) [1] 0.9786558 > my.fcn(10) [1] 1.022467 You can't expect integrate() to return a sensible answer if you don't give it a function that returns consistent results. Duncan Murdoch > >> my.fcn = function(mu){ > + m = 1000 > + z = 0 > + z.mse = 0 > + for(i in 1:m){ > + z[i] = rnorm(1, mu, 1) > + z.mse = z.mse + (z[i] - mu)^2 > + } > + return(z.mse/m) > + } > >> my.fcn(-10) > [1] 1.021711 >> my.fcn(10) > [1] 0.9995235 >> my.fcn(-5) > [1] 1.012727 >> my.fcn(5) > [1] 1.033595 >> my.fcn(0) > [1] 1.106282 >> > The function takes the value (approx) 1 over the range of mu. So the > integration result should be close to 20 if we integrate over (-10, 10). > But R gives: > >> integrate(my.fcn, -10, 10) > 685.4941 with absolute error < 7.6e-12 > >> integrate(Vectorize(my.fcn), -10, 10) # this code never terminate > > I have seen in the "?integrate" page it is clearly written: > > If the function is approximately constant (in particular, zero) over > nearly all its range it is possible that the result and error estimate > may be seriously wrong. > > But this doesn't help in solving the problem. > Thanks, > Santanu > > > > JPSM, 1218J Lefrak Hall > University of Maryland, College Park > Phone: 301-314-9916 > > > > [[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 > and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] integrate
Hi, I am trying to integrate a function which is approximately constant over the range of the integration. The function is as follows: > my.fcn = function(mu){ + m = 1000 + z = 0 + z.mse = 0 + for(i in 1:m){ + z[i] = rnorm(1, mu, 1) + z.mse = z.mse + (z[i] - mu)^2 + } + return(z.mse/m) + } > my.fcn(-10) [1] 1.021711 > my.fcn(10) [1] 0.9995235 > my.fcn(-5) [1] 1.012727 > my.fcn(5) [1] 1.033595 > my.fcn(0) [1] 1.106282 > The function takes the value (approx) 1 over the range of mu. So the integration result should be close to 20 if we integrate over (-10, 10). But R gives: > integrate(my.fcn, -10, 10) 685.4941 with absolute error < 7.6e-12 > integrate(Vectorize(my.fcn), -10, 10) # this code never terminate I have seen in the "?integrate" page it is clearly written: If the function is approximately constant (in particular, zero) over nearly all its range it is possible that the result and error estimate may be seriously wrong. But this doesn't help in solving the problem. Thanks, Santanu JPSM, 1218J Lefrak Hall University of Maryland, College Park Phone: 301-314-9916 [[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 and provide commented, minimal, self-contained, reproducible code.
[R] integrate
Hi, I am trying to integrate a function which is approximately constant over the range of the integration. The function is as follows: > my.fcn = function(mu){ + m = 1000 + z = 0 + z.mse = 0 + for(i in 1:m){ + z[i] = rnorm(1, mu, 1) + z.mse = z.mse + (z[i] - mu)^2 + } + return(z.mse/m) + } > my.fcn(-10) [1] 1.021711 > my.fcn(10) [1] 0.9995235 > my.fcn(-5) [1] 1.012727 > my.fcn(5) [1] 1.033595 > my.fcn(0) [1] 1.106282 > The function takes the value (approx) 1 over the range of mu. So the integration result should be close to 20 if we integrate over (-10, 10). But R gives: > integrate(my.fcn, -10, 10) 685.4941 with absolute error < 7.6e-12 > integrate(Vectorize(my.fcn), -10, 10) # this code never terminate I have seen in the "?integrate" page it is clearly written: If the function is approximately constant (in particular, zero) over nearly all its range it is possible that the result and error estimate may be seriously wrong. But this doesn't help in solving the problem. Thanks, Santanu JPSM, 1218J Lefrak Hall University of Maryland, College Park Phone: 301-314-9916 [[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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] All subsets regression
On Wed, 22 Aug 2007, Alan Harrison wrote: > Hey folks, > > I'm trying to do all subsets on a zero-inflated poisson regression. > I'm aware of the leaps and regsubsets functions but I don't know if they > work for ZIP regressions or how the syntax fits in for them. They don't. -thomas Thomas Lumley Assoc. Professor, Biostatistics [EMAIL PROTECTED] University of Washington, Seattle __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] tackle memory insufficiency for large dataset using save() & load()??
On Tue, 21 Aug 2007, Jessica Z wrote: [snip] I did not notice a comment on this bit in the other replies: >> >> newdata <- load ("compact_d.Rdata") >> >> summary(newdata) > Length Class Mode >1 character character newdata is a string whose value is 'd' try print( newdata ) ls() should tell you there are two objects - 'd' and 'newdata' So just continue using 'd', e.g. summary( d ) HTH, Chuck [snip] Charles C. Berry(858) 534-2098 Dept of Family/Preventive Medicine E mailto:[EMAIL PROTECTED] UC San Diego http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Output from while and for loop
On Tue, 21 Aug 2007, Ryan Briscoe Runquist wrote: > > Hello, > > I am new and am having a hard time getting the proper syntax for output > from loops. I am working on a simulation to generate a null expectation of > bee behavior. Pieces of it work. The part that I am having specific > difficulty is in output of a vector from within the while loop that I am > using. Hmmm. I think I can guess at your problem. You need to assign the result to a vector or list element along these lines: my.list <- list() k <- 0 while( keep.going ){ k <- k+1 ## do lots of stuff ## details omitted my.list[[ k ]] <- result.of.doing.lots.of.stuff } unlist( my.list ) returns a vector if result.of.doing.lots.of.stuff is atomic. There are lots of variations on this. If the result ('my.list') will be very long and each element is atomic, it helps speed to set it up a vector ahead of time and then fill one element per pass thru the loop. Basically the simulation works as such: I have a starting point > and a neighbor matrix and a certain threshold distance for travel. In the > while loop the "bee" moves to a randomly chosen neighbor location. I want > to be able to record the elevations of these points (including the starting > point) so that I can look at variance in elevation and mean elevation. The > loop itself works as does the calculation of the final elevation list, > change in elevation list, and true total distance traveled. I have looked > in all of the email archives but have not come across a correct way of > doing it. Code below: > > start.elev.list<-list() > final.mean.elev.list<-list() > final.elev.list<-list() > final.distance.list<-list() > final.delta.elev.list<-list() > final.var.elev<-list() > > > b<-length(Bees.Day.1$bee) > for (bee in 1:b){ > #this is for number of bees that are trackable in the day with starting > points and threshold distances > elev.current.vector<-vector(mode="numeric", length=0) > count<-1 > ElevSS<-0 > d.traveled<-0 > thresh<-Bees.Day.1$cum.dist[bee] > n<-Bees.Day.1$grid.pt[bee] > #I'm making this up for the threshold, want to be bee specific > #current.point<-round(runif(1,1,n)) #random starting point > current.point<-Day.1.neighbor.matrix[1,n] > #I want to specify the first point in the matrix > Elev.Sum<-Day.1.elev.vector[current.point] > > > while(d.traveled #which of the four options will be selection > transition<-round(runif(1,1,4)) > > #so, what's the new point? > new.point<- Day.1.neighbor.matrix[transition,n] > > #what is the variance in elevation changed > Elev.current<-Day.1.elev.vector[current.point] > elev.current.vector[i]<-Elev.current > Elev.new<-Day.1.elev.vector[new.point] > Elev.Sum<-(Elev.Sum+Elev.new) > > #how far will bee travelled > current.travel<- Day.1.distance.matrix[current.point, new.point] > d.traveled<- current.travel + d.traveled > current.point<- new.point > > #Number of iterations until we reach the threshold > count<-count+1 > > } > > print(count) > print(elev.current.vector) > mean.elev<-Elev.Sum/count > print(paste("Final mean elev for bee", bee, "is", mean.elev, sep=" ")) > final.mean.elev.list[bee]<-list(mean.elev) > > #What was the start elevation? > start.elev<-Day.1.elev.vector[n] > print(paste("Start elev for bee",bee,"is",start.elev, sep=" ")) > start.elev.list[bee]<-list(start.elev) > > #what is the final elevation? > final.elev<-Day.1.elev.vector[current.point] > print(paste("Final elev for bee",bee,"is", final.elev,sep=" ")) > final.elev.list[bee]<-list(final.elev) > > print(paste("Final travel distance for bee", bee,"is", d.traveled, sep=" ")) > final.distance.list[bee]<-list(d.traveled) > > net.delta.elev<-(final.elev-Day.1.elev.vector[n]) > print(paste("Final net change in elevation for bee",bee,"is", > net.delta.elev,sep=" ")) > final.delta.elev.list[bee]<-list(net.delta.elev) > > } > ~~ > Ryan D. Briscoe Runquist > Population Biology Graduate Group > University of California, Davis > [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 > and provide commented, minimal, self-contained, reproducible code. > Charles C. Berry(858) 534-2098 Dept of Family/Preventive Medicine E mailto:[EMAIL PROTECTED] UC San Diego http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] plotting lda results
Read ?plot.lda, which tells you the ... arguments are (for dimen=1, the only option for two groups) passed to ldahist, so then read its help page. I don't know what you want (and your example is not reproducible): I would expect you to get a single plot with two panels (figures), but there are options to have a single panel. (Reading 'An Introduction to R' may help you to use standard terminology that others will be able to follow.) On Wed, 22 Aug 2007, Silvia Lomascolo wrote: > > Hi all, > I am trying to plot the results of a discriminant analysis done with > lda(MASS) but my groups appear in two different plots (in the same graphics > device) and I want to combine them in one plot. My code looks like: > > BirdTrain.lda <- lda(Bdisperser~., data=BirdTrain.mx) > predict(BirdTrain.lda) > plot(BirdTrain.lda) > > I have two types of Bdisperser, so I only get one linear discriminant > function. Can anyone please tell me how to combine the data in one plot? > > I work with R 2.4.1 using Windows. But the version of MASS is what is relevant, and it would have been in the sessionInfo() output the R posting guide asked you for. -- 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 and provide commented, minimal, self-contained, reproducible code.
[R] All subsets regression
Hey folks, I'm trying to do all subsets on a zero-inflated poisson regression. I'm aware of the leaps and regsubsets functions but I don't know if they work for ZIP regressions or how the syntax fits in for them. Can anyone help? The model syntax is: > zip.zc <- zicounts(resp=Scars~.,x=~Location + Lar + Mass + Lac + Lar:Mass + > Location:Mass, data = alan2, distr="ZIP") Many thanks Alan Harrison Quercus Queen's University Belfast MBC, 97 Lisburn Road Belfast BT9 7BL T: 02890 972219 M: 07798615682 [[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 and provide commented, minimal, self-contained, reproducible code.
[R] plotting lda results
Hi all, I am trying to plot the results of a discriminant analysis done with lda(MASS) but my groups appear in two different plots (in the same graphics device) and I want to combine them in one plot. My code looks like: BirdTrain.lda <- lda(Bdisperser~., data=BirdTrain.mx) predict(BirdTrain.lda) plot(BirdTrain.lda) I have two types of Bdisperser, so I only get one linear discriminant function. Can anyone please tell me how to combine the data in one plot? I work with R 2.4.1 using Windows. -- View this message in context: http://www.nabble.com/plotting-lda-results-tf4312870.html#a12279005 Sent from the R help mailing list archive at Nabble.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 and provide commented, minimal, self-contained, reproducible code.
[R] distance by vegan
How to calculate sorensen (bray-curtis) distance by dist function within the vegan package? Cheers Duccio __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Rdonlp2 control parameter question
Hi. I have been trying to work with the Rdonlp2 package. There is a control parameter called epsx that has the following description: epsx(1e-5) - successful temination is indicated if the Kuhn-Tucker criteria are satisfied within the value I have tried to use the donlp2.control() function to adjust the value of epsx, e.g.: >cntlx1_3<-donlp2.control(epsx=1e-03) >cntlx1_3$epsx [1] 0.001 > I have adjusted it larger (as above) and smaller, but it has not had any effect on the solution when I run donlp2(...,control=cntlx1_3). And the console output from running the donlp2() function reports a value of epsx, which is always "epsx= 1.000e-05". Am I doing something wrong? Is there a way to adjust the value of epsx? I've been comparing donlp2() with quadprog() on a quadratic programming problem; the objective function values are pretty close, but Id like to see if I can get them any closer by adjusting epsx. Any help is appreciated. Thanks! -- TMK -- 212-460-5430home 917-656-5351cell __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] 2nd label category on x axis
i'm plotting time-series and label the x-axis as follows: r <- as.POSIXct(round(range(p1$time), "month")) plot(p1$time,p1$ p1, type="l", xaxt="n") axis.POSIXct(1, at=seq(r[1], r[2], by="month"), format="%Y.%m") what I want do do now, is a second label for the x-axis, that stands lower and indicates the years. I don't know how to proceed with all the possibilities in axis, par and plot thanks for your help __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Help with vector gymnastics
Philip - I don't know if this is the "best" way, but it gives you the output you want. Using your tf, vals <- rle(ifelse(tf, 5*which(tf), 0)) vals$values[vals$values == 0] <- vals$values[which(vals$values==0) - 1] inverse.rle(vals) [1] 5 5 5 5 25 30 30 Gladwin, Philip wrote: > Hello, > > What is the best way of solving this problem? > > answer <- ifelse(tf=TRUE, i * 5, previous answer) > where as an initial condition > tf[1] <- TRUE > > > For example if, > tf <- c(T,F,F,F,T,T,F) > over i = 1 to 7 > then the output of the function will be > answer = 5 5 5 5 25 30 30 > > Thank you. > > Phil, > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Help with vector gymnastics
library(zoo) tf <- c(T,F,F,F,T,T,F) i <- seq(7) answer <- ifelse(tf, i*5, NA) answer <- na.locf(answer) On 8/23/07, Gladwin, Philip <[EMAIL PROTECTED]> wrote: > Hello, > > What is the best way of solving this problem? > > answer <- ifelse(tf=TRUE, i * 5, previous answer) > where as an initial condition > tf[1] <- TRUE > > > For example if, > tf <- c(T,F,F,F,T,T,F) > over i = 1 to 7 > then the output of the function will be > answer = 5 5 5 5 25 30 30 > > Thank you. > > Phil, > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > -- Felix Andrews / 安福立 PhD candidate Integrated Catchment Assessment and Management Centre The Fenner School of Environment and Society The Australian National University (Building 48A), ACT 0200 Beijing Bag, Locked Bag 40, Kingston ACT 2604 http://www.neurofractal.org/felix/ voice:+86_1051404394 (in China) mobile:+86_13522529265 (in China) mobile:+61_410400963 (in Australia) xmpp:[EMAIL PROTECTED] 3358 543D AAC6 22C2 D336 80D9 360B 72DD 3E4C F5D8 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] within-subject factors in lme
> > If I understand correctly, you want to include the interactions > > between the random and fixed terms? > > Yes that is exactly I wanted to model. > > > This is done like: > > > > model.lme <- lme(Beta ~ Trust*Sex*Freq, > > random = ~Trust*Sex*Freq|Subj, Model) > > > > But this needs a lot of observations as quite a few > > parameters need to be estimated! > > Well, I tried this as well, but it seems R kept hanging there and > never finished the modeling. It is very likely due to some > singularity as you suspected about the large number of parameters > needed to estimate. But this is not a problem with aov. So does > it mean that I can't run a similar model to that in aov with lme? It depends what you mean by 'similar'. You could still include some of the interactions, e.g. by random = ~(Trust+Sex+Freq)^2|Subj, or even further reduced such as ~Trust+Sex+Freq|Subj. I am not very familiar with aov, but I would suspect that the model you calcualted in aov is not really the same than the one with all possible interactions in lme. In any case, I would personally trust lme much more than aov. > but I feel this is not good enough to account for cross-subject > variations for those interactions. Why wouldn't those patterned > variance-covariance matrix specifications work as I mentioned in > my previous mail? Any more thoughts and suggestions? Sorry, I have never really worked with those. Lorenz - Lorenz Gygax Centre for proper housing of ruminants and pigs Agroscope Reckenholz-Tänikon Research Station ART __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Help with vector gymnastics
Hello, What is the best way of solving this problem? answer <- ifelse(tf=TRUE, i * 5, previous answer) where as an initial condition tf[1] <- TRUE For example if, tf <- c(T,F,F,F,T,T,F) over i = 1 to 7 then the output of the function will be answer = 5 5 5 5 25 30 30 Thank you. Phil, __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] how do i use the get function to obtain an element from alist...
To Mark Leeds and the others, thank you for solving my problem, and so quickly, JM __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] within-subject factors in lme
Hi Lorenz, I really appreciate your comments. > If I understand correctly, you want to include the interactions > between the random and fixed terms? Yes that is exactly I wanted to model. > This is done like: > > model.lme <- lme(Beta ~ Trust*Sex*Freq, > random = ~Trust*Sex*Freq|Subj, Model) > > But this needs a lot of observations as quite a few parameters need > to be estimated! Well, I tried this as well, but it seems R kept hanging there and never finished the modeling. It is very likely due to some singularity as you suspected about the large number of parameters needed to estimate. But this is not a problem with aov. So does it mean that I can't run a similar model to that in aov with lme? Sure I can simply run the following model fit.lme <- lme(Beta ~ Trust*Sex*Freq, random = ~1|Subj, Model) but I feel this is not good enough to account for cross-subject variations for those interactions. Why wouldn't those patterned variance-covariance matrix specifications work as I mentioned in my previous mail? Any more thoughts and suggestions? > Possibly, you can not include the variable Sex in this, because I > assume that Subj is nested within Sex. No, Sex is NOT a subject classifying factor. Instead it is a task- related within-subject factor. Again thanks a lot for your help, Gang On Aug 22, 2007, at 6:52 AM, <[EMAIL PROTECTED]> <[EMAIL PROTECTED]> wrote: > I don't think, this has been answered: > >> I'm trying to run a 3-way within-subject anova in lme with 3 >> fixed factors (Trust, Sex, and Freq), but get stuck with handling >> the random effects. As I want to include all the possible random >> effects in the model, it would be something more or less >> equivalent to using aov >> >>> fit.aov <- aov(Beta ~ >> Trust*Sex*Freq+Error(Subj/(Trust*Sex*Freq)), >> Model) >> >> However I'm not so sure what I should do in lme. Sure >> >>> lme(Beta ~ Trust*Sex*Freq, random = ~1|Subj, Model) >> >> works fine, but it only models the random effect of the >> intercept. I tried the following 4 possibilities: > > If I understand correctly, you want to include the interactions > between the random and fixed terms? This is done like: > > model.lme <- lme(Beta ~ Trust*Sex*Freq, > random = ~Trust*Sex*Freq|Subj, Model) > > But this needs a lot of observations as quite a few parameters need > to be estimated! Possibly, you can not include the variable Sex in > this, because I assume that Subj is nested within Sex. If you just > refer to within and between subject effects and their corresponding > degrees of freedom: you should see this being handled automatically > and correctly by lme e.g. in the output of anova (model.lme) > > Lorenz > - > Lorenz Gygax > Centre for proper housing of ruminants and pigs > Agroscope Reckenholz-Tänikon Research Station ART > Tänikon, CH-8356 Ettenhausen / Switzerland __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Evaluating f(x(2,3)) on a function f<- function(a,b){a+b}
On 8/22/07, Søren Højsgaard <[EMAIL PROTECTED]> wrote: > I have a function and a vector, say > f <- function(a,b){a+b} > x <- c(2,3) > I want to "evaluate f on x" in the sense of computing f(x[1],x[2]). I would > like it to be so that I can write f(x). (I know I can write a wrapper > function g <- function(x){f(x[1],x[2])}, but this is not really what I am > looking for). Is there a general way doing this (programmatically)? (E.g. by > "unpacking" the elements of x and putting them in the "right places" when > calling f...) > I've looked under formals, alist etc. but so far without luck. I hope that the following helps: > f <- function(x) {sum(x)} > f(c(2,3)) [1] 5 > f(c(2,3,5)) [1] 10 > 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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] Evaluating f(x(2,3)) on a function f<- function(a,b){a+b}
Try this: do.call(f, as.list(x)) On 8/22/07, Søren Højsgaard <[EMAIL PROTECTED]> wrote: > Dear list > I have a function and a vector, say >f <- function(a,b){a+b} >x <- c(2,3) > I want to "evaluate f on x" in the sense of computing f(x[1],x[2]). I would > like it to be so that I can write f(x). (I know I can write a wrapper > function g <- function(x){f(x[1],x[2])}, but this is not really what I am > looking for). Is there a general way doing this (programmatically)? (E.g. by > "unpacking" the elements of x and putting them in the "right places" when > calling f...) > I've looked under formals, alist etc. but so far without luck. > > Regards > Søren > > >[[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 > and provide commented, minimal, self-contained, reproducible code. > > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Evaluating f(x(2,3)) on a function f<- function(a,b){a+b}
Dear list I have a function and a vector, say f <- function(a,b){a+b} x <- c(2,3) I want to "evaluate f on x" in the sense of computing f(x[1],x[2]). I would like it to be so that I can write f(x). (I know I can write a wrapper function g <- function(x){f(x[1],x[2])}, but this is not really what I am looking for). Is there a general way doing this (programmatically)? (E.g. by "unpacking" the elements of x and putting them in the "right places" when calling f...) I've looked under formals, alist etc. but so far without luck. Regards Søren [[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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] Formatting Sweave in R-News
Hi Duncan, Duncan Murdoch wrote: >> Slightly off Arjuns topic is a problem with Schunk, Sinput and Soutput >> environments. They just use to much inline space. I have tweade this >> problem with sed (see bellow for Makefile content), but wonder if there >> is a better solution. > > You can change the spacing by redefining those environments. I use this: > > % This removes the extra spacing after code and output chunks in Sweave, > % but keeps the spacing around the whole block. > > \fvset{listparameters={\setlength{\topsep}{0pt}}} > \renewenvironment{Schunk}{\vspace{\topsep}}{\vspace{\topsep}} Nice trick. Thank you. I am CCing to a friend that had the same issue. -- Lep pozdrav / With regards, Gregor Gorjanc -- University of Ljubljana PhD student Biotechnical Facultywww: http://www.bfro.uni-lj.si/MR/ggorjan Zootechnical Department blog: http://ggorjan.blogspot.com Groblje 3 mail: gregor.gorjanc bfro.uni-lj.si SI-1230 Domzale fax: +386 (0)1 72 17 888 Slovenia, Europetel: +386 (0)1 72 17 861 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Synchronzing workspaces
Eric Turkheimer wrote: > How do people go about synchronizing multiple workspaces on different > workstations? I tend to wind up with projects spread around the various > machines I work on. I find that placing the directories on a server and > reading them remotely tends to slow things down. If R were to store all its workspace data objects in individual files instead of one big .RData file, then you could use a revision control system like SVN. Check out the data, work on it, check it in, then on another machine just update to get the changes. However SVN doesn't work too well for binary files - conflicts being hard to resolve without someone backing down - so maybe its not such a good idea anyway... On unix boxes and derivatives, you can keep things in sync efficiently with the 'rsync' command. I think there are GUI addons for it, and Windows ports. Barry __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Synchronzing workspaces
On 8/22/2007 9:20 AM, Eric Turkheimer wrote: > How do people go about synchronizing multiple workspaces on different > workstations? I tend to wind up with projects spread around the various > machines I work on. I find that placing the directories on a server and > reading them remotely tends to slow things down. I use Subversion. It works best with text, but can handle binary files too. So the idea would be to write scripts that either generate or load the data I want to work on, and source those at the beginning of a session. I always try to start with a clean workspace. Separate the code that does complex calculations into functions (and consider organizing it into a package, if it's big enough); put those functions in separate files from the scripts that run them on particular examples. This way you can run source("fns.R") to update the function definitions without re-doing all the simulations in your project. Duncan Murdoch __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Synchronzing workspaces
You could use a versioning system like Subversion or Git on the server but work with the local repository. I believe Git is the easier of the two to set up. Best, Jim Eric Turkheimer wrote: > How do people go about synchronizing multiple workspaces on different > workstations? I tend to wind up with projects spread around the various > machines I work on. I find that placing the directories on a server and > reading them remotely tends to slow things down. > > thanks, > Eric > -- James W. MacDonald, M.S. Biostatistician Affymetrix and cDNA Microarray Core University of Michigan Cancer Center 1500 E. Medical Center Drive 7410 CCGC Ann Arbor MI 48109 734-647-5623 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Synchronzing workspaces
How do people go about synchronizing multiple workspaces on different workstations? I tend to wind up with projects spread around the various machines I work on. I find that placing the directories on a server and reading them remotely tends to slow things down. thanks, Eric -- Eric Turkheimer, PhD Department of Psychology University of Virginia PO Box 400400 Charlottesville, VA 22904-4400 434-982-4732 434-982-4766 (FAX) [[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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] Processing Sweave documents without Sweave
On 8/22/2007 6:40 AM, Renaud Lancelot wrote: > Sweave produces a TeX file which can be exchanged like any other such > file, providing that you give the possible graphic files (and bib, > etc.) called from the TeX code. Yes, and since Sweave doesn't touch the formatting of the TeX code, it's not hard to recognize and merge changes to the .tex file back into the Rnw file. Duncan Murdoch > > Best, > > Renaud > > 2007/8/22, Werner Wernersen <[EMAIL PROTECTED]>: >> Hi, >> >> I am very intrigued by the idea of integrating >> statistical analysis directly with a paper as Sweave >> does it. But as I am collaborating with several people >> and they won't have set up R and also they're very >> unlikely interested in learning it, I am concerned >> about how such an Latex/Sweave document will work on >> their side. >> >> Can they Latex-compile the document without having R / >> Sweave installed and still see the numbers and figures >> I produced by an earlier run of Latex / R / Sweave? Is >> there any option or so planned for this case? Or how >> are you coping with this problem otherwise? >> >> Many thanks, >> Werner >> >> __ >> R-help@stat.math.ethz.ch mailing list >> https://stat.ethz.ch/mailman/listinfo/r-help >> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html >> and provide commented, minimal, self-contained, reproducible code. >> > > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Formatting Sweave in R-News
On 8/22/2007 3:20 AM, Gregor Gorjanc wrote: > Arjun Ravi Narayan arjunnarayan.com> writes: >> I am editing a document for submission to the R-news newsletter, and >> in my article my Sweave code inserts a dynamically generated PDF >> report that my R program generates. >> > > Slightly off Arjuns topic is a problem with Schunk, Sinput and Soutput > environments. They just use to much inline space. I have tweade this > problem with sed (see bellow for Makefile content), but wonder if there > is a better solution. You can change the spacing by redefining those environments. I use this: % This removes the extra spacing after code and output chunks in Sweave, % but keeps the spacing around the whole block. \fvset{listparameters={\setlength{\topsep}{0pt}}} \renewenvironment{Schunk}{\vspace{\topsep}}{\vspace{\topsep}} > Thanks, Gregor > > default: Sweave fixTex > texi2pdf --clean wrapper.tex && evince wrapper.pdf > > Sweave: # Sweave myPaper.Rnw > R CMD Sweave myPaper.Rnw > > fixTex: Sweave # Change all S* environments to smallverbatim > @cat myPaper.tex | sed -e '/\\begin{Sinput}/d' \ > -e '/\\end{Sinput}/d' \ > -e '/\\begin{Soutput}/d' \ > -e '/\\end{Soutput}/d' \ > -e 's/\\begin{Schunk}/\\begin{smallverbatim}/g' \ > -e 's/\\end{Schunk}/\\end{smallverbatim}/g' \ > > myPaper2.tex > @mv -f myPaper2.tex myPaper.tex > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] C code generators
On 8/21/2007 10:20 PM, [EMAIL PROTECTED] wrote: > Dear R-helpers > > Are there any established R packages that include a C code generator -- > that generates new C language files and compiles them? Oleg Sklyar's inline package does that. It takes input in the form of a fragment of C (or C++, Fortran, etc.), completes the source code, compiles it, loads it, and creates an R function to call it. See the ?cfunction man page for examples. I wouldn't recommend this for extensive use: you're better off putting your own code in a package instead. But it's good for quick tests, or writing the first version of the code that goes in the package. > To be precise what I'm looking for is a process that takes text input in > some format (it might be pseudocode, fragments of C code, etc) and creates > a valid C language source file that can be compiled by R CMD COMPILE. > Ideally the procedure should also cause the C code to be compiled and > dynamically loaded. > > To give a trivial example, suppose I want to be able to perform image > filtering operations on matrices. All filters have the same structure: > each position [i,j] in the matrix is visited in a double loop; a > calculation (depending on the filter) is performed using the value at > [i,j] and also the values at neighbouring positions [i+1,j] , [i+1,j+1] > etc; the result is written to the output matrix at the position [i,j]. The > C code for the loop is always the same; only a few lines of code that > perform the filter calculation will change. I would like a procedure that > accepts a text file containing just a few lines of C code or pseudocode, > and inserts these lines into the appropriate place in the loop, producing > a valid C language routine which can then be compiled by R CMD COMPILE and > dynamically loaded. As a matter of fact, one of Oleg's examples is an image filter. Duncan Murdoch > > (Of course, it can get more complicated than just inserting a single text > fragment!) > > I once implemented such a feature in an image processing package, so I > know it's not hard. Before dusting off this ancient code I would like to > learn whether there's an R package or other open source program that > already does it. > > Any pointers would be welcome > > thanks > Adrian Baddeley > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] IRT model nonparametric from Ramsay in R
Hi: I need to apply the IRT model from Ramsay (1991), He apply the Smoothing Kernel to multiple choice test. Is possible in R?. Thank you, Xavier G. Ordonez Doctoral Student Universidad Complutense de Madrid __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Optimization problem
Try this. 1. following Ben remove the Randalstown point and reset the levels of the Location factor. 2. then replace solve with ginv so it uses the generalized inverse to calculate the hessian: alan2 <- subset(alan, subset = Location != "Randalstown") alan2$Location <- factor(as.character(alan2$Location)) library(MASS) solve <- ginv zinb.zc <- zicounts(resp=Scars~.,x =~Location + Lar + Mass + Lar:Mass + Location:Mass,z =~Location + Lar + Mass + Lar:Mass + Location:Mass, data = alan2) rm(solve) On 8/21/07, Ben Bolker <[EMAIL PROTECTED]> wrote: > > (Hope this gets threaded properly. Sorry if it doesn't.) > > Gabor: Lac and Lacfac being the same is irrelevant, wouldn't > produce NAs (but would produce something like a singular Hessian > and maybe other problems) -- but they're not even specified in this > model. > > The bottom line is that you have a location with a single > observation, so the GLM that zicounts runs to get the initial > parameter values has an unestimable location:mass interaction > for one location, so it gives an NA, so optim complains. > > In gruesome detail: > > ## set up data > scardat = read.table("scars.dat",header=TRUE) > library(zicounts) > ## try to run model > zinb.zc <- zicounts(resp=Scars~., >x =~Location + Lar + Mass + Lar:Mass + Location:Mass, >z =~Location + Lar + Mass + Lar:Mass + Location:Mass, >data=scardat) > ## tried to debug this by dumping zicounts.R to a file, modifying > ## it to put a "trace" argument in that would print out the parameters > ## and log-likelihood for every call to the log-likelihood function. > dump("zicounts",file="zicounts.R") > source("zicounts.R") > zinb.zc <- zicounts(resp=Scars~., >x =~Location + Lar + Mass + Lar:Mass + Location:Mass, >z =~Location + Lar + Mass + Lar:Mass + Location:Mass, >data=scardat,trace=TRUE) > ## this actually didn't do any good because the negative log-likelihood > ## function never gets called -- as it turns out optim() barfs when it > ## gets its initial values, before it ever gets to evaluating the > log-likelihood > > ## check the glm -- this is the equivalent of what zicounts does to > ## get the initial values of the x parameters > p1 <- glm(Scars~Location + Lar + Mass + Lar:Mass + Location:Mass, > data=scardat,family="poisson") > which(is.na(coef(p1))) > > ## find out what the deal is > table(scardat$Location) > > scar2 = subset(scardat,Location!="Randalstown") > ## first step to removing the bad point from the data set -- but ... > table(scar2$Location) > ## it leaves the Location factor with the same levels, so > ## now we have ZERO counts for one location: > ## redefine the factor to drop unused levels > scar2$Location <- factor(scar2$Location) > ## OK, looks fine now > table(scar2$Location) > > zinb.zc <- zicounts(resp=Scars~., >x =~Location + Lar + Mass + Lar:Mass + Location:Mass, >z =~Location + Lar + Mass + Lar:Mass + Location:Mass, >data=scar2) > ## now we get another error ("system is computationally singular" when > ## trying to compute Hessian -- overparameterized?) Not in any > ## trivial way that I can see. It would be nice to get into the guts > ## of zicounts and stop it from trying to invert the Hessian, which is > ## I think where this happens. > > In the meanwhile, I have some other ideas about this analysis (sorry, > but you started it ...) > > Looking at the data in a few different ways: > > library(lattice) > xyplot(Scars~Mass,groups=Location,data=scar2,jitter=TRUE, > auto.key=list(columns=3)) > xyplot(Scars~Mass|Location,data=scar2,jitter=TRUE) > > xyplot(Scars~Lar,groups=Location,data=scar2, > auto.key=list(columns=3)) > xyplot(Scars~Mass|Lar,data=scar2) > xyplot(Scars~Lar|Location,data=scar2) > > Some thoughts: (1) I'm not at all sure that > zero-inflation is necessary (see Warton 2005, Environmentrics). > This is a fairly small, noisy data set without huge numbers > of zeros -- a plain old negative binomial might be fine. > > I don't actually see a lot of signal here, period (although there may > be some) ... > there's not a huge range in Lar (whatever it is -- the rest of the > covariates I > think I can interpret). It would be tempting to try to fit location as > a random > effect, because fitting all those extra degrees of freedom is going to > kill you. > On the other hand, GLMMs are a bit hairy. > > cheers > Ben > > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Ans on: How do i print a "main title" on a win.graph with several plots?
Thanks Mr.Pr. Ripley, For pointing me in the good direction, i use win.graph() so i get an overview of the different graphs. Until now, i v had no problems with it, hope it stays that way. for those whom have the same porblem, this is what i do: win.graph(); op <- par(mfrow = c(1,2) , oma= c(1,1,1,1)) boxplot(lg_value~labo, main="Test 1 at day 1",ylab="log(x) ", xlab="different lab's", data=dataset,ylim=c(-0.05,5)) mtext("At day 1",line=-1,outer=TRUE) boxplot(lg_value~labo, main="Test 2 at day 1",ylab="log(x)", xlab="different lab's", data=dataset,ylim=c(-0.05,5)) par(op); Kind regards, Tom W. Disclaimer: click here [[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 and provide commented, minimal, self-contained, reproducible code.
[R] ordering of variables in plots
I'm making plots where the x-axis variables are called l hl s (for "long" "half-long" and "short"). I'd like them to appear in that order, but they turn up in alphabetical order: hl l s. Except for fiddling with the names of the variables, how can I fix this? That is, how can I specify what order I want? Thanks, Elin __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] within-subject factors in lme
I don't think, this has been answered: > I'm trying to run a 3-way within-subject anova in lme with 3 > fixed factors (Trust, Sex, and Freq), but get stuck with handling > the random effects. As I want to include all the possible random > effects in the model, it would be something more or less > equivalent to using aov > > > fit.aov <- aov(Beta ~ > Trust*Sex*Freq+Error(Subj/(Trust*Sex*Freq)), > Model) > > However I'm not so sure what I should do in lme. Sure > > > lme(Beta ~ Trust*Sex*Freq, random = ~1|Subj, Model) > > works fine, but it only models the random effect of the > intercept. I tried the following 4 possibilities: If I understand correctly, you want to include the interactions between the random and fixed terms? This is done like: model.lme <- lme(Beta ~ Trust*Sex*Freq, random = ~Trust*Sex*Freq|Subj, Model) But this needs a lot of observations as quite a few parameters need to be estimated! Possibly, you can not include the variable Sex in this, because I assume that Subj is nested within Sex. If you just refer to within and between subject effects and their corresponding degrees of freedom: you should see this being handled automatically and correctly by lme e.g. in the output of anova (model.lme) Lorenz - Lorenz Gygax Centre for proper housing of ruminants and pigs Agroscope Reckenholz-Tänikon Research Station ART Tänikon, CH-8356 Ettenhausen / Switzerland __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Processing Sweave documents without Sweave
That sound's perfect, that's all I needed to know. Thanks for your kind help, Werner --- Renaud Lancelot <[EMAIL PROTECTED]> schrieb: > Sweave produces a TeX file which can be exchanged > like any other such > file, providing that you give the possible graphic > files (and bib, > etc.) called from the TeX code. > > Best, > > Renaud > > 2007/8/22, Werner Wernersen > <[EMAIL PROTECTED]>: > > Hi, > > > > I am very intrigued by the idea of integrating > > statistical analysis directly with a paper as > Sweave > > does it. But as I am collaborating with several > people > > and they won't have set up R and also they're very > > unlikely interested in learning it, I am concerned > > about how such an Latex/Sweave document will work > on > > their side. > > > > Can they Latex-compile the document without having > R / > > Sweave installed and still see the numbers and > figures > > I produced by an earlier run of Latex / R / > Sweave? Is > > there any option or so planned for this case? Or > how > > are you coping with this problem otherwise? > > > > Many thanks, > > Werner > > > > __ > > R-help@stat.math.ethz.ch mailing list > > https://stat.ethz.ch/mailman/listinfo/r-help > > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > > and provide commented, minimal, self-contained, > reproducible code. > > > > > -- > Renaud LANCELOT > Département Systèmes Biologiques du CIRAD > CIRAD, Biological Systems Department > > Campus International de Baillarguet > TA 30 / B > F34398 Montpellier > Tel +33 (0)4 67 59 37 17 > Secr. +33 (0)4 67 59 37 37 > Fax +33 (0)4 67 59 37 95 > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Processing Sweave documents without Sweave
Sweave produces a TeX file which can be exchanged like any other such file, providing that you give the possible graphic files (and bib, etc.) called from the TeX code. Best, Renaud 2007/8/22, Werner Wernersen <[EMAIL PROTECTED]>: > Hi, > > I am very intrigued by the idea of integrating > statistical analysis directly with a paper as Sweave > does it. But as I am collaborating with several people > and they won't have set up R and also they're very > unlikely interested in learning it, I am concerned > about how such an Latex/Sweave document will work on > their side. > > Can they Latex-compile the document without having R / > Sweave installed and still see the numbers and figures > I produced by an earlier run of Latex / R / Sweave? Is > there any option or so planned for this case? Or how > are you coping with this problem otherwise? > > Many thanks, > Werner > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > -- Renaud LANCELOT Département Systèmes Biologiques du CIRAD CIRAD, Biological Systems Department Campus International de Baillarguet TA 30 / B F34398 Montpellier Tel +33 (0)4 67 59 37 17 Secr. +33 (0)4 67 59 37 37 Fax +33 (0)4 67 59 37 95 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How do i print a "main title" on a win.graph with several plots?
?title, look at the 'outer' argument. You can see further discussion of the outer margins in 'An Introduction to R'. I don't know why you are using win.graph(): it is a deprecated form of windows() with many of the arguments taking unchangable defaults. On Wed, 22 Aug 2007, Tom Willems wrote: > Good Mornig All, > > How R you today? ;o) > > I have lots of questions, but i l start with the simplest one, > to wich i am shy to say, i did not find the answer. > > It is the following: > > When i make a summary plot like for example plot( summary(glm)), > i get one window, one main title, and 4 graph's in that window. > > Now i do know how to get several graphs in one window, > buit i don't manage geting one main title , in the top middel of the > window. > I can give every plot a different main and subtitle, but i can't put no > title in the "win.graph()" box. > > this is what i do: > > win.graph(); op <- par(mfrow = c(1,2)) > boxplot(lg_value~labo, main="Test 1 at day 1",ylab="log(x) ", > xlab="different lab's", data=dataset,ylim=c(-0.05,5)) > boxplot(lg_value~labo, main="Test 2 at day 1",ylab="log(x)", > xlab="different lab's", data=dataset,ylim=c(-0.05,5)) > par(op); > > Then i get one graph window, two plots each with their onw main title. > What i'd like to have is, one main title saying " At day 1", and then two > plots with the different tests. > How can i do this, pls? > > Kind regards, > Tom W. > > > Disclaimer: click here > [[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 > and provide commented, minimal, self-contained, reproducible code. > -- 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 and provide commented, minimal, self-contained, reproducible code.
[R] Processing Sweave documents without Sweave
Hi, I am very intrigued by the idea of integrating statistical analysis directly with a paper as Sweave does it. But as I am collaborating with several people and they won't have set up R and also they're very unlikely interested in learning it, I am concerned about how such an Latex/Sweave document will work on their side. Can they Latex-compile the document without having R / Sweave installed and still see the numbers and figures I produced by an earlier run of Latex / R / Sweave? Is there any option or so planned for this case? Or how are you coping with this problem otherwise? Many thanks, Werner __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Recurrence plots for CGH
>CGH Hi, I guess you should start looking at the bioconductor project. http://www.bioconductor.org/packages/release/Software.html best wishes, ido __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] divided scatter plots
Caroline Nganga wrote: > I have a data set which contains two columns. The first column is a > list of countries, and the second column contains their political risk > ratings. I would like to create one large plot that contains 5 > different sections, each with a scatter plot. To clarify, I have > divided the countries into 5 groups. For each group (continent), I > would like to have the name of the continent on the x-axis, and points > representing countries and their risk rating on the y-axis. However, > I want all 5 scatter plots to be in one large plot. What function > should I use to do this? Also, is it possible to label each point? > thanks for any help! > Hi Caroline, If I understand your request, you might be able to use the axis.break function in the plotrix package. That is, you make one big plot with the points in five columns and then put gap style axis breaks between the columns. Here's a toy example: library(plotrix) prr.df<-data.frame(country=c("us","mx","ca","br","ar","pe", "ch","mn","in","nl","fr","es","na","mz","rw"), continent=c(1,1,1,2,2,2,3,3,3,4,4,4,5,5,5), prr=rnorm(15)+4) plot(prr.df$continent,prr.df$prr,main="Political risk ratings", xlim=c(0.7,5.3),xlab="Continent",ylab="Risk rating",type="n") text(prr.df$continent,prr.df$prr,prr.df$country) axis.break(1,1.5,style="gap") axis.break(1,2.5,style="gap") axis.break(1,3.5,style="gap") axis.break(1,4.5,style="gap") Jim __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] plot in for loop with "i" indexed variables does not work
you need something like this: par(mfrow = c(5, 5)) for (i in 1:10) { nam <- paste("Var.", i, sep = "") plot(x = time, y = mat[, nam], xlab = "Time", ylab = nam) } where `mat' is the matrix containing Var.1, Var.2, Var.3, etc. I hope it helps. Best, Dimitris Dimitris Rizopoulos Ph.D. Student Biostatistical Centre School of Public Health Catholic University of Leuven Address: Kapucijnenvoer 35, Leuven, Belgium Tel: +32/(0)16/336899 Fax: +32/(0)16/337015 Web: http://med.kuleuven.be/biostat/ http://www.student.kuleuven.be/~m0390867/dimitris.htm - Original Message - From: <[EMAIL PROTECTED]> To: Sent: Wednesday, August 22, 2007 10:56 AM Subject: [R] plot in for loop with "i" indexed variables does not work > Hi there > > Does anyone know, why this is not working?: > > > par(mfrow = c(5, 5)) > > for (i in 1:10){ > plot(x=time, y=Var.i) > } > > (x-variable is time > y-variable is Var out of matrix with columns Var.1 Var.2 Var.3 > etc...) > > even > > for (i in 1:10){ > plot(x=time, y=paste("Var", i, sep=".") > } > > does not work > > Thanks > > marc > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] plot in for loop with "i" indexed variables does not work
Hi there Does anyone know, why this is not working?: par(mfrow = c(5, 5)) for (i in 1:10){ plot(x=time, y=Var.i) } (x-variable is time y-variable is Var out of matrix with columns Var.1 Var.2 Var.3 etc...) even for (i in 1:10){ plot(x=time, y=paste("Var", i, sep=".") } does not work Thanks marc __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] library(fCalendar) timeDate("12.03.2005", format="%d.%m.%Y")
Martin Maechler wrote: >> "OL" == Ola Lindqvist <[EMAIL PROTECTED]> >> on Tue, 21 Aug 2007 14:32:19 +0200 writes: >> > > OL> Thanks! > OL> Seems to work fine now! > > Well, for your example. > > But sorry to say, the patch breaks other cases. > > I'm investigating further > (and will hopefully contribute to a new CRAN release of fCalendar > once Diethelm Wuertz is back from wherever; I've already made > more changes) > > Good to see experts investigating :-) Sorry, I should have mentioned, at least, that I am not an Rmetrics developer, that I only took a quick glance at the code (focused on the reported problem), and that the patched version is neither official nor tested. Now I see that I completely messed up the case for "Date" input and sometimes .midnightStandard() was not working correct any more. One further issue, which was brought to my attention by another R-user off list: I think, that all occurences of if (sum(lt$sec+lt$min+lt$hour) == 0) isoFormat = "%Y-%m-%d" in R/3A-TimeDateClass.R could safely be replaced by if (sum(lt$sec+lt$min+lt$hour,na.rm=TRUE) == 0) isoFormat = "%Y-%m-%d" in order to make fCalendar more NA-friendly. Currently, at least for some cases, a single NA entry in the input vector suffices to let timeDate() break down. Again, I am not sure, if this change breaks some special cases, or if breaking down in case of NAs is intended. Nevertheless, I have built a new, unofficial and untested patched (windows-)version of fCalendar (same weblink), which makes the "Date" input case and .midnightStandard() working again (I hope so, at least), and incorporates the (improved?) NA-handling. But of course: only use this version with care (and only if the official version doesn't work for your particular problem). Hopefully, we see an official update on CRAN soon :-) Kind regards, Martin Becker > Martin Maechler, > ETH Zurich [but different department than D.Wuertz] > > > OL> Martin Becker wrote: > >> Dear Ola, > >> > >> I think you spotted a small bug in *package* fCalendar. > >> Explicit specification should prevent "autodetection" of the date > >> format, which is not the case for fCalendar v251.70, instead > >> autodetection is done at least once (twice, if actually appropriate). > >> With the following patch, things should work ok: > >> > >> diff --recursive fCalendar.orig/R/3A-TimeDateClass.R > >> fCalendar/R/3A-TimeDateClass.R > >> 433c433 > >> < charvec = format(strptime(charvec, .whichFormat(charvec)), > >> isoFormat) > >> --- > >> > charvec = format(strptime(charvec, format), isoFormat) > >> > >> You did not provide the output of sessionInfo() (which you are asked > >> for in the posting guide). If you are using Windows and don't know how > >> to apply the patch, you can download a patched binary version here: > >> http://www.saar-gate.net/download/fCalendar_251.70.zip > >> > >> Regards, > >> > >> Martin > >> > >> PS: Maybe r-sig-finance is more appropriate for questions concerning > >> Rmetrics. > >> > >> > >> Ola Lindqvist wrote: > >>> Dear R users, > >>> I have problem with the library fCalendar. > >>> > >>> I am not using the US standard format notations. It seems like it is > >>> not possible to have different format than the US standards. > >>> Anyone how knows a way to go around this problem? > >>> > >>> Here is the code I enter: > >>> myDate = "12.03.2005" > >>> timeDate(myDate, format = "%d.%m.%Y") > >>> > >>> And I get following error message: > >>> Error in if (sum(lt$sec + lt$min + lt$hour) == 0) isoFormat = > >>> "%Y-%m-%d" : > >>> missing value where TRUE/FALSE needed > >>> > >>> Thanks, > >>> Ola > >>> > >>> __ > >>> R-help@stat.math.ethz.ch mailing list > >>> https://stat.ethz.ch/mailman/listinfo/r-help > >>> PLEASE do read the posting guide > >>> http://www.R-project.org/posting-guide.html > >>> and provide commented, minimal, self-contained, reproducible code. > >>> > >> > > OL> __ > OL> R-help@stat.math.ethz.ch mailing list > OL> https://stat.ethz.ch/mailman/listinfo/r-help > OL> PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > OL> and provide commented, minimal, self-contained, reproducible code. > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] How do i print a "main title" on a win.graph with several plots?
Good Mornig All, How R you today? ;o) I have lots of questions, but i l start with the simplest one, to wich i am shy to say, i did not find the answer. It is the following: When i make a summary plot like for example plot( summary(glm)), i get one window, one main title, and 4 graph's in that window. Now i do know how to get several graphs in one window, buit i don't manage geting one main title , in the top middel of the window. I can give every plot a different main and subtitle, but i can't put no title in the "win.graph()" box. this is what i do: win.graph(); op <- par(mfrow = c(1,2)) boxplot(lg_value~labo, main="Test 1 at day 1",ylab="log(x) ", xlab="different lab's", data=dataset,ylim=c(-0.05,5)) boxplot(lg_value~labo, main="Test 2 at day 1",ylab="log(x)", xlab="different lab's", data=dataset,ylim=c(-0.05,5)) par(op); Then i get one graph window, two plots each with their onw main title. What i'd like to have is, one main title saying " At day 1", and then two plots with the different tests. How can i do this, pls? Kind regards, Tom W. Disclaimer: click here [[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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] Random Sampling from a Matrix
Anup Nandialath wrote: > Dear Friends, > > I have a matrix of size 5000 X 20. The first two columns are indicator > variables taking the value of either 0 or 1. Let us call the first two > columns Y1 and Y2. > > I need to randomly sample 1000 rows with all the associated columns, in other > words my new matrix should be of size 1000 X 20. I realize that using this > command > > newmat <- mainmat[sample(1000,replace=F),] > > achieves this. However, I would like to make sure that both Y1 and Y2 have > more or less an equal amount of 0's and 1's. At present when I sample, I get > cases where sometimes all my Y2's are 0. Is there any way to accomodate this > problem. Yes, for example drawing samples stratified by the values of Y1 or Y2. Uwe Ligges > Thanks in advance. > > Regards > > Anup > > > - > > [[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 > and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] rectify a program of seasonal dummies matrix
Friedrich Schuster wrote: > Hello, > > the main problem seems to be the "if else", should be "else if". > > Your code is hard to read, maybe you should consider using more () {}: > > T <- 100; > br <- matrix(0,T,4); Thanks for the contribution. Please note: a) It is a bad idea to have a variable called T. Some people still use it as a logical value even if they should not. b) R does not need any ";" at the end of a line. Uwe Ligges > for (i in 1:T) { >for (j in 1:4) { > if (i==j) { > br[i,j] <- 1; >} >else if ((abs(i-j)%%4)==0) { > br[i,j] <- 1; > } > else { > br[i,j] <- 0; > } > } > } > > A simpler approach is creating a diagonal matrix and multply it : > > # create small diagonal matrix > mat = diag(x=1, nrow=4, ncol=4); > mat > # multiply diagonal matrix and re-dimension it to 4 cols > br <- rep(mat, 25); > dim(br) <- c(100, 4); > br; > > Hope this helps, > FS > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Odp: displaying the means on a scatter diagram
Hi [EMAIL PROTECTED] napsal dne 22.08.2007 09:10:26: > > Hello, > > I created a simple scatter diagram: > > x= c(53, 52, 55, 65, 71, 69, 75, 78, 88, 70) > y= c(162, 165, 165, 171, 173, 175, 179, 181, 184, 176) > plot(x, y) > > Now I would like to display the mean on that diagram. do you think points(mean(x), mean(y)) or e.g. abline(h=mean(y)) abline(v=mean(x)) Regards Petr > > Can anyone tell me what possibilities I have to do that? > Thanks in advance > > Tobias > > -- > View this message in context: http://www.nabble.com/displaying-the-means-on-a- > scatter-diagram-tf4309832.html#a12269270 > Sent from the R help mailing list archive at Nabble.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 > and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Formatting Sweave in R-News
Arjun Ravi Narayan arjunnarayan.com> writes: > I am editing a document for submission to the R-news newsletter, and > in my article my Sweave code inserts a dynamically generated PDF > report that my R program generates. > Slightly off Arjuns topic is a problem with Schunk, Sinput and Soutput environments. They just use to much inline space. I have tweade this problem with sed (see bellow for Makefile content), but wonder if there is a better solution. Thanks, Gregor default: Sweave fixTex texi2pdf --clean wrapper.tex && evince wrapper.pdf Sweave: # Sweave myPaper.Rnw R CMD Sweave myPaper.Rnw fixTex: Sweave # Change all S* environments to smallverbatim @cat myPaper.tex | sed -e '/\\begin{Sinput}/d' \ -e '/\\end{Sinput}/d' \ -e '/\\begin{Soutput}/d' \ -e '/\\end{Soutput}/d' \ -e 's/\\begin{Schunk}/\\begin{smallverbatim}/g' \ -e 's/\\end{Schunk}/\\end{smallverbatim}/g' \ > myPaper2.tex @mv -f myPaper2.tex myPaper.tex __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] displaying the means on a scatter diagram
Hello, I created a simple scatter diagram: x= c(53, 52, 55, 65, 71, 69, 75, 78, 88, 70) y= c(162, 165, 165, 171, 173, 175, 179, 181, 184, 176) plot(x, y) Now I would like to display the mean on that diagram. Can anyone tell me what possibilities I have to do that? Thanks in advance Tobias -- View this message in context: http://www.nabble.com/displaying-the-means-on-a-scatter-diagram-tf4309832.html#a12269270 Sent from the R help mailing list archive at Nabble.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 and provide commented, minimal, self-contained, reproducible code.
Re: [R] Time conversion problems
[EMAIL PROTECTED] wrote: > Hi there > > I have precipitation data from 2004 to 2006 in varying resolutions (10 to > 20min intervals) with time in seconds from beginnig of the year (summation) > and a second variable as year. > > I applied follwing code to convert the time into a date: > > times<-strptime("2004-01-01", "%Y-%m-%d", tz="GMT") + precipitation$time1 > > everytihng went well, except that every year, the seconds-counter starts by > zero, therefore I have now three 2004 series instead of going further from 04 > to 05 etc. > > I tried to sum the last seconds-values of 2004 to the first of 2005 with an > if command like: > > if (year=2005) time2=time1+632489 ;(seconds) > > but it doesn't work. > > thanks for a solution > Can't you just do strptime(paste(year, "01-01", sep="-"), ? (or use ISOdatetime(year,1,1,0,0,0,tz="GMT")) __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] [R-pkgs] brew 1.0-1
brew implements a templating framework for mixing text and R code for report generation. brew template syntax is similar to PHP, Ruby's erb module, Java Server Pages, and Python's psp module. brew is written in R with no package dependencies, and it's not just for the web. It can be used as an alternative to Sweave in a limited context. See the brew-test-1.brew file in the distribution for some salient differences between the two. brew can also complement Sweave since it can be written to do conditional inclusion of or loop over Sweave code chunks. The 1.0-1 version should show up on the CRAN mirrors shortly, but in the mean time it can be got from: http://www.rforge.net/brew/ Best, Jeff -- http://biostat.mc.vanderbilt.edu/JeffreyHorner ___ 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 and provide commented, minimal, self-contained, reproducible code.