Re: [R] converting MATLAB -> R | element-wise operation
On Tue, 27 Feb 2024 14:54:26 -0500 Evan Cooch wrote: > So, trying to convert a very long, somewhat technical bit of lin alg > MATLAB code to R. Most of it working, but raninto a stumbling block > that is probaably simple enough for someone to explain. On https://cran.r-project.org/other-docs.html the documents: “Matlab® / R Reference” by David Hiebeler (PDF, 2010-05-25, 52 pages) and “R and Octave” by Robin Hankin (Text), a reference sheet translating between the most common Octave (or Matlab) and R commands might be useful. IIRC, I once used them when I had to do something in Matlab/Octave and used them as reverse-lookup (I know how to do it in R, how is it done in Matlab??). :-) Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] converting MATLAB -> R | element-wise operation
On Tue, 27 Feb 2024 13:51:25 -0800 Jeff Newmiller via R-help wrote: > The fundamental data type in Matlab is a matrix... they don't have > vectors, they have Nx1 matrices and 1xM matrices. Also known as column vectors and row vectors. :) > Vectors don't have any concept of "row" vs. "column". They do in (numerical) linear algebra. And MATLAB was written by numerical analysts for numerical analysts. :-) So they distinguish between row and column vectors. GAUSS also distinguishes between row and column vectors. R (and S) does not distinguish between row and column vectors, and is in this aspect quite unique among the groups of vector-oriented programming languages (AFAICT). But treating every vector as a column vector, together with the recycling rule, allows for the easy coding of the matrix/vector calculations that one (mostly) comes across in statistical computing. For the rare occasion when this is not true the sweep() command is provided, typically remembered once one was bitten by the lack of distinction between row and column vectors. :) Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] [External] converting MATLAB -> R | element-wise operation
On Tue, 27 Feb 2024 21:37:52 + "Richard M. Heiberger" wrote: > > t(t(NN)/lambda) > [,1] [,2] [,3] > [1,] 0.5 0.667 0.75 > [2,] 2.0 1.667 1.50 > > > > R matrices are column-based. MATLAB matrices are row-based. It might depend on what you mean with this statement, but I would be very surprised if MATLAB is not storing matrices in column-major form, just as R does. NN = [1, 2, 3; 4, 5, 6]; and NN = [ [1, 4]' , [2, 5]', [3, 6]' ]; produce the same matrix in MATLAB. So the default filling of matrices is by row, and there is no convenient argument to change that. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Packages sometimes don't update, but no error or warning is thrown
G'day Philipp, On Fri, 16 Feb 2024 17:33:13 +0100 Philipp Schneider wrote: > Thanks for all the input. It's happening again. This time for the > packages "DBI", "parallelly", "segmented", "survival", "V8". So, > RStudio shows updates for those and updating them via RStudio leads > to this output: Perhaps, as Duncan suggested, you should use update.packages() instead of install.packages(). Also, what happens if you try to do the update in R and not from RStudio? I.e. start R and say update.packages() in R? What does the command 'getOption(pkgType)' output on your machine? It should be "both" for either binary or source package. If it is "mac.binary" you should investigate where that is set and why your system is configured to install binary packages early. Finally, you can always recall the command: > install.packages(c("DBI", "parallelly", "segmented", "survival", "V8")) and modify it to: > install.packages(c("DBI", "parallelly", "segmented", "survival", "V8"), type > = "source") This, of course, assumes that you have the tools installed to compile source packages. Best wishes, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Packages sometimes don't update, but no error or warning is thrown
G'day Philipp, On Tue, 13 Feb 2024 09:59:17 +0100 gernophil--- via R-help wrote: > this question is related to this > (https://community.rstudio.com/t/packages-are-not-updating/166214/3), > [...] > To sum it up: If I am updating packages (be it via Bioconductor or > CRAN) some packages simply don’t update, [...] > I would expect any kind of message that the package will not be > updated, since no newer binary is available or a prompt, if I want to > compile from source. RStudio is doing its own thing for some task, including 'install.packages()' (and for some reasons, at least on the platforms on which I use RStudio, RStudio calls 'install.packages()' and not 'update.packages()' when an update is requested via the GUI). See: RStudio> install.packages function (...) .rs.callAs(name, hook, original, ...) compared to: R> install.packages function (pkgs, lib, repos = getOption("repos"), contriburl = contrib.url(repos, type), method, available = NULL, destdir = NULL, dependencies = NA, type = getOption("pkgType"), configure.args = getOption("configure.args"), configure.vars = getOption("configure.vars"), clean = FALSE, Ncpus = getOption("Ncpus", 1L), verbose = getOption("verbose"), libs_only = FALSE, INSTALL_opts, quiet = FALSE, keep_outputs = FALSE, ...) { [...] So if you use Install/Update in the Packages tab of RStudio and do not experience the behaviour you are expecting, it is something that you need to discuss with Posit, not with R. :) > However, the only message I get is: > ``` > trying URL '' The package name has the version number encoded in it, so theoretical you should be able to tell at this point whether the package that is downloaded is the version that is already installed, hence no update will happen. Best wishes, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] different results between cor and ccf
G'day Patrick, On Tue, 16 Jan 2024 09:19:40 +0100 Patrick Giraudoux wrote: [...] > So far so good, but when I lag one of the series, I cannot find the > same correlation as with ccf > > > cor(x[1:(length(x)-1)],y[2:length(y)]) [1] -0.7903428 > > ... where I expect -0.668 based on ccf > > Can anyone explain why ? The difference is explained by cff() seeing the complete data on x and y and calculating the sample means only once, which are then used in the calculations for each lag. cor() sees only the data you pass down, so calculates different estimates for the means of the two sequences. To illustrate: [...first execute your code...] R> xx <- x-mean(x) R> yy <- y-mean(y) R> n <- length(x) R> vx <- sum(xx^2)/n R> vy <- sum(yy^2)/n R> (c0 <- sum(xx*yy)/n/sqrt(vx*vy)) [1] -0.5948694 R> xx <- x[1:(length(x)-1)] - mean(x) R> yy <- y[2:length(y)] - mean(y) R> (c1 <- sum(xx*yy)/n/sqrt(vx*vy)) [1] -0.6676418 The help page of cff() points to MASS, 4ed, the more specific reference is p 389ff. :) Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] non-linear regression and root finding
G'day Troels, On Tue, 7 Nov 2023 07:14:02 +0100 Troels Ring wrote: > Be as it may, I wonder if not your method might work if only we KNOW > that pK1 is either positive OR negative, in which case we have pK1 = > -exp(theta1)? If pK1 can be either negative or positive (or 0 :-) ), and it is just the ordering that you want to have enforced, then I would try the parameterisation: pK1 = pK1 pK2 = pK1 + exp(theta2) pK3 = pk2 + exp(theta3) and optimise over pK1, theta2 and theta3. As long as you want to know the estimates only. Asking for standard errors of the original estimates would open another can of worms. :-) Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] non-linear regression and root finding
G'day Troels, On Mon, 6 Nov 2023 20:43:10 +0100 Troels Ring wrote: > Thanks a lot! This was amazing. I'm not sure I see how the conditiion > pK1 < pK2 < pK3 is enforced? One way of enforcing such constraints (well, in finite computer arithemtic only "<=" can be enforced) is to rewrite the parameters as: pK1 = exp(theta1) ## only if pK1 > 0 pK2 = pK1 + exp(theta2) pK3 = pk2 + exp(theta3) And then use your optimiser to optimise over theta1, theta2 and theta3. There might be better approaches depending on the specific problem. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Issues when trying to fit a nonlinear regression model
G'day Paul, On Sun, 20 Aug 2023 12:15:08 -0500 Paul Bernal wrote: > Any idea on how to proceed in this situation? What could I do? You are fitting a simple asymptotic model for which nls() can find good starting values if you use the self starting models (SSxyz()). Well, Doug (et al.) choose to parameterise the asymptotic model differently, but you can easily change to your parameterisation if you want: ``` fm1 <- nls(y ~ SSasymp(x, Asym, R0, lrc), data=mod14data2_random) theta1 <- coef(fm1)["Asym"] theta2 <- coef(fm1)["Asym"] - coef(fm1)["R0"] theta3 <- exp(coef(fm1)["lrc"]) fm2 <- nls(y ~ theta1 - theta2 * exp(-theta3*x), start=list(theta1=theta1, theta2=theta2, theta3=theta3), data=mod14data2_random) summary(fm2) ``` Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] geom_smooth
G'day Thomas, On Sat, 12 Aug 2023 04:17:42 + (UTC) Thomas Subia via R-help wrote: > Here is my reproducible code for a graph using geom_smooth The call "library(tidyverse)" was missing. :) > I'd like to add a black boundary around the shaded area. I suspect > this can be done with geom_ribbon but I cannot figure this out. Some > advice would be welcome. This works for me: ggplot(scatter_data,aes(x=x_var,y=y_var,))+ geom_point()+ geom_smooth(se=TRUE,fill="blue",color="black",linetype="dashed") + geom_ribbon(stat="smooth", aes(ymin=after_stat(ymin), ymax=after_stat(ymax)), fill=NA, color="black")+ theme_cowplot() Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] suprising behaviour of tryCatch()
G'day Henrik, On Thu, 18 May 2023 08:35:38 -0700 Henrik Bengtsson wrote: > ... or just put the R expression inside curly brackets, e.g. I was wondering whether I was having a fortunes::fortune(106) moment. :-) Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] suprising behaviour of tryCatch()
G'day Greg, On Thu, 18 May 2023 08:57:55 -0600 Greg Snow <538...@gmail.com> wrote: [...] > `fun(x <- expr)` Assigns the value of `expr` to the variable `x` in > the frame/environment that `fun` is called from. Only if the argument 'x' is evaluated during the function call AFAIK. If it is not, due to lazy evaluation rules, then no assignment takes place. > When you run the code with `<-` it, then the ith element of the global > variable `sexsnp` is assigned the p-value. When you run the version > with `=` then R tries to find an argument to `tryCatch` that matches > (possibly partially) `sexsnp[i]` and gives the error because it does > not find a match (not sure why it is not handled by `...`, but > tryCatch may be parsed special, or non-legal argument names on the RHS > of `-` may be checked). After exact matching on names comes partial matching on names and then positional matching. If after these three rounds there are still actual arguments left over, they are assigned to "..." (if it exist as formal argument). So under the usual rules of how actual arguments are passed to formal arguments in Federico's call the argument "sexsnp[i] = fisher.test(table(data[,3], data[,i + 38]))$p" should be assigned to the formal argument "expr" of tryCatch() and "error = function(e) print(NA))" should be absorbed by "...". But it seems, as "=" has a different meaning in function call, the expression passed to the "expr" formal argument is not allowed to contain a "=". Not sure why. :-) I guess that the intended use is actually: R> for(i in 1:1750){sexsnp[i] = tryCatch(fisher.test(table(data[,3], data[,i + 38]))$p, error = function(e) print(NA))} I.e. tryCatch() returns the result of the evaluated expression (if successful, otherwise the result of an error handler) and that result can be assigned. It is not expected the that expression contains an assignment operation. But I am just guessing. :) Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] suprising behaviour of tryCatch()
G'day Federico, On Wed, 17 May 2023 10:42:17 + "Calboli Federico (LUKE)" wrote: > sexsnp = rep(NA, 1750) > for(i in 1:1750){tryCatch(sexsnp[i] = fisher.test(table(data[,3], > data[,i + 38]))$p, error = function(e) print(NA))} Error: unexpected > '=' in "for(i in 1:1750){tryCatch(sexsnp[i] =" Try: R> for(i in 1:1750){tryCatch(eval(expression("sexsnp[i] = fisher.test(table(data[,3], data[,i+38]))$p")), error=function(e)print(NA))} or R> for(i in 1:1750){tryCatch(bquote("sexsnp[i] = fisher.test(table(data[,3], data[,i+38]))$p"), error=function(e) print(NA))} or R> for(i in 1:1750){tryCatch(.("sexsnp[i] = fisher.test(table(data[,3], data[,i+38]))$p"), error=function(e) print(NA))} If you want to use the '='. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] suprising behaviour of tryCatch()
G'day Federico, On Wed, 17 May 2023 10:42:17 + "Calboli Federico (LUKE)" wrote: > I_d be obliged if someone can explain why tryCatch assigns items with _<-_ > and not _=_. It is not just tryCatch but any function. In function calls the "=" is used to assign actual arguments to formal arguments, not for assignments. That is why "plot(fm <- lm(eruptions ~ waiting, faithful)" works but "plot(fm = lm(eruptions ~ waiting, faithful)" does not. The only universal assignment operator, AFAIK, is "<-". Well, and "->". Long time ago "_" used to be one too, but that is long ago. The usual advice was that if you do not know when "=" does not work as assignment operator then always use "<-" for assignments. :) The final thing you have to take care of that in function calls you have to be sure that the argument is actually evaluated (R has lazy evaluation) to ensure that the assignment is actually executed. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] question about update()
G'day Adelchi, hope all is well with you. On Thu, 4 May 2023 10:34:00 +0200 Adelchi Azzalini via R-help wrote: > Thanks, Duncan. What you indicate is surely the ideal route. > Unfortunately, in my case this is not feasible, because the > construction of xf and the update call are within an iterative > procedure where xf is changed at each iteration, so that the steps > > obj$data <- cbind(obj$data, xf=xf) > new.obj <- update(obj, . ~ . + xf) > > must be repeated hundreds of times, each with a different xf. If memory serves correctly, update() takes the object that is passed to it, looks at what the call was that created that object, modifies that call according to the additional arguments, and finally executes the modified call. So there is a lot of manipulations going on in update(). In particular it would result each time in a call to lm(), glm() or whatever call was used to create the object. Inside any of these modelling functions a lot of symbolic manipulations/calculations are needed too (parsing the formula, creating the design matrix and response vector from the parsed formula and data frame, checking if weights are used &c). If you do the same calculation essentially over and over again, just with minor modification, all these symbolic manipulations are just time consuming. IMHO, you will be better off to bypass update() and just use lm.fit() (for which lm() is a nice front-end) and glm.fit() (for which glm() is a nice front-end), or whatever routine does the grunt work of fitting the model to the data in your application (hopefully, the package creator used a set up of XXX.fit() to fit the model, called by XXX() that does all the fancy formula handling). Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Trying to learn how to write an "advanced" function
G'day John, For these snippets to produce what I think you want them to produce it is just necessary to define doit() as follows: doit <- function(x){ lm(formula=x) } R> # run formula y~x R> JD <- doit(y~x) R> JD Call: lm(formula = x) Coefficients: (Intercept)x 0.8403 0.8558 R> # run formula y~x+z R> JD2 <- doit(y~x+z) R> JD2 Call: lm(formula = x) Coefficients: (Intercept)xz -0.7490 0.6237 1.9023 Cheers, Berwin On Thu, 16 Mar 2023 14:53:33 + "Sorkin, John" wrote: > Although I owe thanks to Ramus and Ivan, I still do not know how to > write and "advanced" function. > > My most recent try (after looking at the material Ramus and Ivan set) > still does not work. I am trying to run the lm function on two > different formulae: 1) y~x, 2) y~x+z > Any corrections would be appreciated! > > Thank you, > John > > > doit <- function(x){ > ds <- deparse(substitute(x)) > cat("1\n") > print(ds) > eval(lm(quote(ds)),parent.frame()) > } > > # define data that will be used in regression > y <- 1:10 > x <- y+rnorm(10) > z <- c(rep(1,5),rep(2,5)) > # Show what x, y and z look like > rbind(x,y,z) > > # run formula y~x > JD <- doit(y~x) > JD > > # run formula y~x+z > JD2 <- doit(y~x+z) > JD2 __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Environmental oddity --- reproducible example.
G'day Rolf, On Sun, 7 Nov 2021 19:33:40 +1300 Rolf Turner wrote: > library(Deriv) > d1 <- Deriv(dnorm,"sd") > source("d2.txt") # d2.txt is attached > > d1(1,0,3,TRUE) # [1] -0.2962963 > d2(1,0,3,TRUE) # [1] -0.889 Fascinating: R> pryr::call_tree(body(d1)) R> pryr::call_tree(body(d2)) clearly show that the two functions have a different idea to what expression the final "/sd" is applied too (as an earlier poster suggested), but I have no idea why. Deriv() seems to return the correct function, but when it is displayed, the deparser(?) somehow omits a crucial pair of braces. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] I am struggling with contrasts
G'day all, On Tue, 10 Mar 2020 11:07:13 + "Viechtbauer, Wolfgang (SP)" wrote: > The linearHypothesis() function from the 'car' package does this. The function glht() in the 'multcomp' package should also be able to do this. The 'emmeans' package might also be useful. Will be off-line for a while now, but might look at the example again and how to do it with 'multcomp' or 'emmeans' later. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] I am struggling with contrasts
G'day John, On Tue, 10 Mar 2020 01:42:46 + "Sorkin, John" wrote: > I am running a Poisson regression with a single outcome variable, > HGE, and a single independent variable, a factor, Group which can be > one of two values, Group1, or Group2. I am trying to define contrasts > that will give me the values of my outcome variable (HGE) when > group=Group1 and when group=Group2. Not sure what you mean, but I am suspecting you are after this output: R> fit0 <- glm(HGE ~ Group - 1,family=poisson,data=dataForR,offset=logFU) R> summary(fit0) Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] does pdf() work inside of a function?
G'day Anuj, On Sun, 28 Jul 2019 16:40:28 -0700 Anuj Goyal wrote: > Why does pdf() work outside a function, but not inside it? Because you use a graphics system that need printing to produce the plot. On the command line auto-printing is your friend, but in functions you have to do it explicitly yourself. > R version 3.6.0 (2019-04-26) on Mac 10.14.5 > > # R --silent --vanilla < s2.r > library(quantmod) > options("getSymbols.warning4.0"=FALSE) > options("getSymbols.yahoo.warning"=FALSE) > > gs <- function(f) { > csvText <- "s,n\nAFL,AFLAC\nAIG,AIG" > csv <- read.csv(text=csvText, stringsAsFactors = FALSE) > symVec <- getSymbols(as.vector(csv$s)) > > # create PDF > fname = paste0(f,".pdf") > pdf(file = fname) > par(mfrow = c( 4,2 ) ) > mapply (chart_Series, mget(symVec)) Change this to mapply( function(x) print(chart_Series(x)), mget(symVec)) and it should work. At least it does on my machine: -rw-r--r-- 1 berwin berwin 154281 Jul 29 12:19 t.pdf Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] R shared library (/usr/lib64/R/lib/libR.so) not found.
G'day Rolf, On Sun, 26 Aug 2018 13:00:34 +1200 Rolf Turner wrote: > On 08/26/2018 02:59 AM, Berwin A Turlach wrote: > > > I am not sure whether I agree. :) > > Huh? Black is white and up is down??? Nope, but as I said, on my machine RStudio and the R installed from the Ubuntu repositories worked out of the box. After adding the CRAN repositories to my apt configuration and upgrading R, RStudio and the new R version worked again out of the box. For the R versions that I compile from source, as I install sub architectures, I need to install symbol links from where libR.so actually is to where Rstudio expects them. And RStudio works out of the box with all these installation, i.e. it seems to follow symbolic links. > I did as advised and it made absolutely no difference. Ergo I was > right. That is your conclusion from the data, are you are perfectly entitled to it. :) I still think that this is one of the 5% of the cases where you got your inference wrong. Given that for me (and presumably many others [but I neither follow r-sig-Debian nor the Rstudio forums, so I may have a biased input]) everything works out of the box, my conclusion from the data is that you have FUBAR'd your system. Now, this can be from the way you tried to compile R from source (you never told us what you exactly specified to ./configure, but we would also have to know whether you made changes to config.site before compilation) or it could be from some environment variables (set long time ago and since forgotten [if so, where is it set ~/.bashrc? ~/.R?). > > > > I keep wondering why you have a /usr/lib64. > > > > > How did you get this? > > I have no idea. It just seems to be there. And it seems that > rstudio thinks that it should be there. As far as I can tell, RStudio expects libR.so in a certain directory relative to where the home directory (as reported by R on query) of R. But RStudio's behaviour can also be influenced by environment variable. So the fact that RStudio keeps looking at /usr/lib64 is for me evidence that your system is compromised, and that RStudio is not starting the version of R installed from the Ubuntu package manager. > Anyway, I seem to have (for the time being at least) a working > system, and I don't feel like wasting any more time on this issue. Fair enough, but I am concerned that this issue will raise it head sooner or later again. So my final advice is to stick to the old adage from sports: Never change a winning team. If you system works for you know, freeze it. No more updates/upgrades. Otherwise, good luck when you try to update R next Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] installing R 3.5.1
Dear Bogdan, On Sat, 25 Aug 2018 07:24:40 -0700 Bogdan Tanasa wrote: > installed E: Unable to correct problems, you have held broken > packages. Perhaps this is the problem, did you try "apt-get -f install"? Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] R shared library (/usr/lib64/R/lib/libR.so) not found.
G'day Rolf, On Sat, 25 Aug 2018 11:20:09 +1200 Rolf Turner wrote: > I was pretty sure that the foregoing was a complete red herring. And > I was right. I am not sure whether I agree. :) > I have been told by younger and wiser heads that installing from > source is The Right Thing to Do. Younger heads probably have more time on their hands to play around with installing from source and debugging when things go wrong. Wiser heads probably have figured out how they should do so in the first place. :) Seriously, years ago I installed quite a bit of software from source, also to always have the newest version. And while I was using Debian unstable I accepted that from time to time an update would leave me with a broken system and I would have to spend some time to fix it again. Nowadays, I rather spend my time on other things and so am happy to go with most software with whatever version the package installer of Ubuntu installs. The one exception is R, but that is for various good reasons: * I want to have the 32 bit and 64 bit architecture installed so that I can test my packages on both packages. * I would like to keep old versions around, to test my packages on them * I like to use the R version that is installed in our computer lab when preparing teaching material, and the newest version otherwise. If I ever get the feeling that all these points are easily achievable with the compiled packages (without having to play around with docker, sandboxes etc), I would stop compiling R from source. > Moreover I'd always had the impression that the version of R provided > by the package manager persistently lags one or two releases behind > the current version. >From the official repositories, yes. But CRAN provides since longer than I can remember up-to-date binary packages for the most common Linux distributions. AFAIK, Debian and Ubuntu packages are available thanks to Dirk Eddelbuettel (initially? mainly?) and others. > The process for installing R using the package manager is far from > straightforward and few people give clear instructions on this issue. 1.) Start a web browser 2.) Go to your favourite CRAN server 3.) Select 'Download R for Linux' 4.) Select the directory "ubuntu/" from the page that is served 5.) Select README.html from the page that is served 6.) Follow the instructions > (Instructions are usually incomplete and full of jargon and acronyms > that the instructors blithely assume assume that the instructees > understand. I am sure (some) of my students have the same complains about my lecturing > (They *don't*! In this instance (mirabile dictu!) I managed (using > Uncle Google) You were lucky in this case, usually Uncle Google serves me with somewhat out-dated (if not wrong) information when I run into a technical problem with linux. > to find very clear and explicit instructions at: > > https://www.digitalocean.com/community/tutorials/how-to-install-r-on-ubuntu-18-04-quickstart > > I followed these instructions, and everything went swimmingly. You were lucky that you are not sitting behind a firewall (as I seem to be, not sure whom I have to thank for that). For me the first step (installing the GPG key) fails and some further googling was necessary. Though, I was working from the CRAN instructions. > Interestingly (???) the "new" R was installed in /usr/bin and not in > /usr/local/bin. Of course, it is a system Ubuntu package, so it installs to /usr/bin. > I then tried issuing the command: > > rstudio > > Exactly the same pop-up error. No help at all, as I expected. O.k., as I got a new machine, I installed (X)ubuntu 18.04 from scratch, no update from an earlier version. The R version installed from r-base-core was R 3.4.4. In a terminal in which I manipulated my PATH variable so that /usr/bin is early on, I could start R without problem, and I could start rstudio without problem. After following the instructions on CRAN and updating my apt information, I have now R 3.5.1 installed. And it starts fine from the command line as does rstudio. And rstudio starts now with R 3.5.1. Both these packages work out of the box. Thus I wonder whether you have somehow messed up your system during the attempt to install R from source. Or have some environment variables set that provide rstudio with wrong information > Then finally, in desperation, I copied libR.so from /usr/lib/R/lib to > /usr/lib64/R/lib. Bingo!!! I can now start Rstudio!!! I keep wondering why you have a /usr/lib64. On my Ubuntu boxes, in /usr I have /lib, /lib32 and /libx32, but no lib64. As far as I know, Debian/Ubuntu 64bit implementations always used /usr/lib for its 64 bit libraries and /usr/lib32 for the 32 bit versions (unlike other distributions who used /usr/lib64 for the former and /usr/lib for the latter). In fact, a 'locate /lib64' on my Ubuntu 18.04 system (less than 2 week old fresh installation) shows: /lib64 /lib64/ld-
Re: [R] installing R 3.5.1
G'day Bogdan, On Fri, 24 Aug 2018 18:28:59 -0700 Bogdan Tanasa wrote: > I am trying to install R 3.5.1 on my Ubuntu 14.04 system; however, I > am getting the following message : > > sudo apt-get install r-base > [...] > The following packages have unmet dependencies: > r-base : Depends: r-recommended (= 3.5.1-1trusty) but it is not > going to be installed > E: Unable to correct problems, you have held broken packages. For me such problems are usually fixed by specifying the package that "is not going to be installed" but on which the package I want to install depends also to apt-get install. What does sudo apt-get install r-base r-recommended do on your system? Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] R shared library (/usr/lib64/R/lib/libR.so) not found.
G'day Peter, On Thu, 23 Aug 2018 08:45:37 -0700 Peter Langfelder wrote: > The manual, specifically > > https://cran.r-project.org/doc/manuals/r-release/R-admin.html#Installation > > documents this way of choosing the installation directory. Yes, with the caveat that one needs GNU or Solaris make (which would be the case on Ubuntu). So it is hardly the recommended way. As I read the R Administration Manual, the recommended way is to specify the location at which you want to install R via ./configure. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] R shared library (/usr/lib64/R/lib/libR.so) not found.
G'day Rolf, On Thu, 23 Aug 2018 23:34:38 +1200 Rolf Turner wrote: > I guess I should have said --- I did > > sudo make prefix=/usr install > > which puts stuff into /usr rather than into /usr/local. ??? I do not remember ever specifying "prefix=foo" at the make install stage. Not for any software that uses autoconf &c. I thought the prefix should be specified to ./configure and after that just make make check make install I am pretty sure that the location of RHOME is set by the path specified (explicitly or implicitly) to ./configure. If you then install R at another location with your construct, some problems seem to be pre-programmed. But I could be wrong. > I forget exactly why I chose (in the dim distant past) to do this ... > I have a vague recollection that my search path was more > "comfortable" that way. In my experience this is a false comfort. Set the search path so that /opt/bin or /usr/local/bin is early on and finds programs you install to those location. Installing to /usr will sooner or later lead to tears if your program "conflicts" with some Ubuntu package (which might have been installed to satisfy the requirement of another package that you needed). If that package is update during an "apt-get update", you can end up with a broken system. > > I also have a recollection that I once installed Rstudio for some > > MOOC, and ended up putting a symlink in somewhere like /usr/lib* , > > because Rstudio was only available as a binary with the location > > of the shared lib hard-baked into it. The location is only hard-coded in relation to RHOME and with the assumption that you are not using a sub-architecture on Ubuntu. AFAIK, binary Ubuntu distributions of R do not use sub-architectures and there should be no problem on an Ubuntu system as long as all software is installed via Ubuntu. But this is probably a question better discussed on r-sig-Debian > So it looks like a symlink might be the answer for me. Only if you can be sure that that libR.so is compatible with the R version that you seem to be using. The official r-base-core package from Ubuntu seems to be R 3.4.4. If you added the CRAN repository to your Ubuntu system, you might have a newer version installed. But if your installation is partly a self-compiled R-3.5.1 version that is then linked to an R 3.4.4 libR.so in /usr/lib/R/lib/ (from r-cran-base), you are inviting trouble. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] R shared library (/usr/lib64/R/lib/libR.so) not found.
G'day Rolf, On Thu, 23 Aug 2018 22:57:35 +1200 Rolf Turner wrote: > I *think* that this is an R question (and *not* an RStudio question!) Others may disagree... :) > I have, somewhat against my better judgement, decided to experiment > with using RStudio. Very good if you are still involved with teaching and need to use the same environment as your student... or want to try some new IDE... If you have a good set-up that works for you, a bit more difficult to see why you want to change... RStudio's editor allegedly has an Emacs style but I find that style more confusing than helpful... half of the short-cuts do not seem to work... But it is a nice IDE... > Then I tried to start RStudio ("rstudio" from the command line) > and got a pop-up window with the error message: > > > R shared library (/usr/lib64/R/lib/libR.so) not found. If this > > is a custom build of R, was it built with the --enable-R-shlib > > option? Yes, I regularly got that too, until I changed my R installation scripts > Oops, no, I guess it wasn't. So I carefully did a > > sudo make uninstall > make clean > make distclean To bad that you did this. There should have been a file called "config.log" in that directory, and the top lines of that file would have told you with which option ./configure was called, in particular whether you had used the --enable-R-shlib flag. > and then did > > ./R-3.5.1/configure > > making sure I added the --enable-R-shlib flag. Well, some of the other flags might also be important... > Then I did make and sudo make install. It all seemed to go ... > but then I did > > rstudio > > again and got the same popup error. > > There is indeed *no* libR.so in /usr/lib64/R/lib. I wonder why rstudio tries to look into /usr/lib64. AFAICT, rstudio queries the R that it uses for its home directory and then expects libR.so to be at a specific location relative to this home directory. And it expects that the installation does not use sub-architectures, that is what tripped me up. > There *is* a libR.so in /usr/lib/R/lib, but (weirdly) ls -l reveals > that it dates from the my previous install of R-3.5.1 for which I > *did not* configure with --enable-R-shlib. Are you sure? I am running Ubuntu 18.04 too. My system has /usr/lib/libR.so and /usr/lib/R/lib/libR.so, with the former being a link to the latter. And these were installed via `r-base-core` which seems to be a requirement for `ess`. (The long list of `ess` on Ubuntu, together with its insistence of installing r-base-core and a whole bunch of r-cran-* package is IMHO ridiculous. Nearly bad enough to make me consider installing ESS from source again.) So the /usr/lib/R/lib/libR.so could be from r-base-core (if you somehow installed that package). Obviously you have sudo rights on your machine, so I would suggest to try: sudo updatedb locate libR.so To see how many libR.so you have installed and where they are > Can anyone explain to me WTF is going on? Not with much more information, e.g. what those "" to .configure were. Also, the great strength of unix system is that you can influence the behaviour of programs via system variables... unfortunately that is also one of its greatest weaknesses when it comes to finding out why programs do not behaved the way you expect them to work. Some stray environment variable might cause all this problem. > What should I do? Just make a symbolic link > from /usr/lib/R/lib/libR.so to /usr/lib64/R/lib/libR.so? I would not recommend this. If this file is from another installation, you are just asking for trouble down the road, which would then be even harder to debug. > It bothers me that /usr/lib/R/lib/libR.so was not "refreshed" from my > most recent install of R. Me too. But I think you should never install to this location in the first place. AFAICT, /usr/lib/R/lib/libR.so is installed by r-core-base, so if you install your own version there and then a "apt-get update" updates r-core-base, you will end up with a broken system. I learned the hard way long time ago not to install any software in areas where Ubuntu packages are installed. I restrict myself to install to /usr/local or /opt (with /opt often being on a separate partition so that material installed there survive if I have to install/upgrade Ubuntu from scratch). > I plead for enlightenment. Not sure whether my comments were very helpful. But you should probably find out why your custom installed version of R tells RStudio to look at /usr/lib64. A "locate libR.pc" could be helpful. on my system this returns /usr/lib/pkgconfig/libR.pc (from r-base-core) and /opt/R/R-3.5.1/lib/pkgconfig/libR.pc (my installation from source). Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/po
Re: [R] Fwd: Quadratic programming, for loop
G'day Maija, On Wed, 27 Jun 2018 08:48:08 +0300 Maija Sirkjärvi wrote: > Thanks for your reply! Unfortunately something is still wrong. > > After the transpose, dvec and Amat are still incompatible. > > > d <- -hsmooth > > dvec <- t(d) > > c <- dvec*Amat > Error in dvec * Amat : non-conformable arrays '*' in R is element-wise multiplication and '%*%' implements matrix/matrix (matrix/vector) multiplication as defined in matrix algebra. I presume you want to use the latter operator here. > Moreover, I don't understand the following: > > > If dvec is of length *J*, then b will be of length J too. > > I believe the length of dvec comes from the number of variables and > the length of b from the number of constraints. In this case they are > not equal. As I said: > > solve.QP solves the quadratic program: > > min(-d^T b + 1/2 b^T D b) > >where A^T b >= b_0. The minimisation is with respect to b. Note that the objective function contains the inner product of d (passed to dvec) and b, so d and b must have the same dimension/length. b contains the parameters/variables over which you want to minimise. b_0 (passed to bvec) depends on the number of constraints. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Fwd: Quadratic programming, for loop
G'day all, On Tue, 26 Jun 2018 11:16:55 +0300 Maija Sirkjärvi wrote: > It seems that my Amat and dvec are incompatible. Amat is a matrix of > zeros size: *2*J-3,J* and dvec is a vector of length *J*. There > should be no problem, but apparently there is. [...] solve.QP solves the quadratic program: min(-d^T b + 1/2 b^T D b) where A^T b >= b_0. Note the transpose. :) If dvec is of length *J*, then b will be of length J too, and Amat should be Jx(2J-3) so that its transpose is (2j-3)xJ, making it compatible for matrix multiplication with b. Cheers, Berwin --- This email has been checked for viruses by Avast antivirus software. https://www.avast.com/antivirus __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Nonlinear regression
G'day Timothy, On Tue, 19 Dec 2017 18:28:00 -0600 Timothy Axberg wrote: > Should I repost the question with reply-all? Nope, we got all from Jeff's post. :) > On Tue, Dec 19, 2017 at 6:13 PM, Jeff Newmiller > wrote: > > > You also need to reply-all so the mailing list stays in the loop. > > -- > > Sent from my phone. Please excuse my brevity. > > > > On December 19, 2017 4:00:29 PM PST, Timothy Axberg < > > axbergtimo...@gmail.com> wrote: > > >Sorry about that. Here is the code typed directly on the email. > > > > > >qe = (Qmax * Kl * ce) / (1 + Kl * ce) [...] > > >##The linearized data > > >celin <- 1/ce > > >qelin <- 1/qe Plotting qelin against celin, I can see why you call this the linearized data. But fitting a linear model to these data obviously does not give you good starting values for nls(). Given your model equation, I would linearize the model to: qe = Qmax*KL * ce + KL * ce*qe and fit a non-intercept linear model to predict qe by ce and ce*qe. From this model I would then determine starting values for nls(). It seems to work with your data set: R> ce <- c(15.17, 42.15, 69.12, 237.7, 419.77) R> qe <- c(17.65, 30.07, 65.36, 81.7, 90.2) R> fit2 <- lm(qe ~ ce + I(ce*qe) - 1) R> summary(fit2) R> Kl <- - coef(fit2)[2] R> Qmax <- coef(fit2)[1]/Kl R> plot(ce, qe) R> c <- seq(min(ce), max(ce)) R> q <- (Qmax*Kl*c)/(1+(Kl*c)) R> lines(c, q) R> fit2 <- nls(qe ~ ((Qmax*Kl*ce)/(1+(Kl*ce))), start = list(Qmax = Qmax,Kl =Kl)) R> summary(fit2) Formula: qe ~ ((Qmax * Kl * ce)/(1 + (Kl * ce))) Parameters: Estimate Std. Error t value Pr(>|t|) Qmax 106.42602 12.82808 8.296 0.00367 ** Kl 0.014560.00543 2.681 0.07496 . --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 9.249 on 3 degrees of freedom Number of iterations to convergence: 2 Achieved convergence tolerance: 6.355e-06 HTH. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Generating help files for a function
G'day Erin, On Sat, 16 Dec 2017 08:00:38 -0600 Erin Hodgess wrote: > I'm in the process of writing a package, and I'm using the lovely "R > Package" book as a guideline. > > However, in the midst of my work, I discovered that I had omitted a > function and am now putting in it the package. Not a problem. But > the problem is the help file. What is the best way to generate a > help file "after the fact" like that, please? It depends on how you decided to write the documentation. If you follow the "R Package" guidelines and use roxygen2, just add the comments for the documentation at the beginning of the file and follow the procedure outline in "R Packages" book. If you are writing the documentation separate, more like the "Writing R Extensions" manual, then (1) start R, (2) source the file in which the function is so that is is in your workspace, (3) say "prompt(foo)" if the function's name is foo and (4) copy the resulting foo.Rd into the /man directory of your package. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] [FORGED] Re: tcltk problems
G'day Peter, On Sat, 18 Nov 2017 23:20:25 +0100 peter dalgaard wrote: > Hum, missed that bit. Looking at the configure script, the only way I > can see it failing to look in /usr/lib/tcl8.6 is if ${LIBnn} is not > "lib". Any chance it might be set to lib64? Well, yes, in an earlier e-mail Rolf said: Then I scanned through BldDir/config.log and found: LIBnn='lib64' (on line 19901 !!!) But that seems to be the output of the configuration as determined by ./configure. The question is why ${LIBnn} is set to this value. If not through some environment variable or config.site, the only way I can see that the relevant code snippet[1] (lines 3788-3806 of R 3.4.2's configure, quoted below) sets ${LIBnn} to 'lib64' is if /usr/lib64 exist. This should not be the case on a 64bit Ubuntu machine. My machine has /usr/lib, /usr/lib32 and /usr/libx32. As far as I know, Ubuntu (and also Debian? But I have not used Debian for a long time) has always placed its 64bit libraries into /lib and the 32bit libraries into e.g. /lib32. And nowadays, one has also /usr/lib/i386-linux-gnu for 32bit libraries and /usr/lib/x86_64-linux-gnu for 64bit libraries. As far as I can tell, RedHat used /lib for the 32bit libraries and /lib64 for the 64bit libraries. And I believe that most (all?) of R core that were using linux were using RedHat. But I could be wrong. Definitely remember that I always found the multiple-architecture part in "R Installation and Administration" confusing to read as my linux system had no /lib64 directories. :) Finally, if ${LIBnn} is set to 'lib' then the relevant code snippet[2] (lines 39527-39543 of R 3.4.2's configure, quoted below), would find tclConfig.sh; regardless on whether it is in /usr/lib or /usr/lib/tcl8.6. Of course, if ${LIBnn} is set to 'lib64', then ./configure will not find tclConfig.sh on an Ubuntu machine. So Rolf should probably check whether there is a /usr/lib64 directory on his machine. If so, why is it there? Perhaps some other software installed from source created it? The contents of the directory might give some clues. If a /usr/lib64 exist on Rolf's machine, he should probably use config.site to set ${LIBnn} to 'lib' (line 146 in R 3.4.2's config.site) so that ./configure can automagically find tclConfig.sh again.. Cheers, Berwin PS: Who is now wondering whether his 32bit compile of R is actually using the 64bit version of TCL/TK on his system One day I might find out... :) [1] ## We need to establish suitable defaults for a 64-bit OS libnn=lib case "${host_os}" in linux*) ## Not all distros use this: some choose to march out of step ## Allow for ppc64le (Debian calls ppc64el), powerpc64le ... case "${host_cpu}" in x86_64|mips64|ppc64*|powerpc64*|sparc64|s390x) if test -d /usr/lib64; then libnn=lib64 fi ;; esac ;; solaris*) ## libnn=lib/sparcv9 ## on 64-bit only, but that's compiler-specific ;; esac : ${LIBnn=$libnn} [2] if test -z "${TCL_CONFIG}"; then { $as_echo "$as_me:${as_lineno-$LINENO}: checking for tclConfig.sh in library (sub)directories" >&5 $as_echo_n "checking for tclConfig.sh in library (sub)directories... " >&6; } if ${r_cv_path_TCL_CONFIG+:} false; then : $as_echo_n "(cached) " >&6 else for ldir in /usr/local/${LIBnn} /usr/${LIBnn} /${LIBnn} /opt/lib /sw/lib /opt/csw/lib /usr/sfw/lib /opt/freeware/lib; do for dir in \ ${ldir} \ `ls -d ${ldir}/tcl[8-9].[0-9]* 2>/dev/null | sort -r`; do if test -f ${dir}/tclConfig.sh; then r_cv_path_TCL_CONFIG="${dir}/tclConfig.sh" break 2 fi done done fi __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] deviance in GLM vs. summary.glm
Dear Yi, On Wed, 31 May 2017 04:53:25 + (UTC) array chip via R-help wrote: > Hi, I am running a logistic regression on a simple dataset (attached) > using glm: [...] As you can see, the interaction term is very > insignificant (p = 0.996)! Well, all terms are not significant (actually, AFAIK, the phrase "very insignificant" does not exist; and if it does, than it ought not). But look at the estimates and the standard errors too (might be easier if you had formatted your e-mail in a way that was readable)! What does an estimate for the intercept of 19.57 on the linear predictor scale mean? What it the estimated probability if you transform back? Perhaps the following command R> xtabs(~x1+x2+y,dat) will shed some more light on what is going on in your data set, and why the interaction term is highly significant. > But if I use a anova() to compare a full > vs reduced model to evaluate the interaction term: > > anova(glm(y~x1+x2,dat,family='binomial'), > > glm(y~x1*x2,dat,family='binomial')) To you know about the (optional) argument 'test="ChiSq"' for the anova() command if you use it to compare models fitted by glm()? (see help(anova.glm)) > So I get very different p value on the interaction term, can someone > share what's going wrong here? Data separation, aka Hauck-Donner phenomenon, discussed in any good book on logistic regression. Best wishes, Berwin == Full address A/Prof Berwin A Turlach, Director Tel.: +61 (8) 6488 3338 (secr) Centre for Applied Statistics +61 (8) 6488 3383 (self) School of Maths and Stats (M019) FAX : +61 (8) 6488 1028 The University of Western Australia 35 Stirling Highway e-mail: berwin.turl...@gmail.com Crawley WA 6009 http://staffhome.ecm.uwa.edu.au/~00043886/ Australiahttp://www.researcherid.com/rid/A-4995-2008 __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] xvfb? cron job updates R packages, fails on some requiring X11
G'day Paul, On Thu, 19 Jan 2017 10:52:16 -0600 Paul Johnson wrote: > In Centos 7 systems, I wrote a script that runs on the cron and I > notice some package updates and installs fail like this: > > [] > > I understand I need something like xvfb to simulate an X11 session, > but I don't understand how to make it work. Can one of you give me an > idiot's guide on what to do? I do not know about Centos (i.e. how to install the necessary software), but on my Ubuntu machines I have since years the following in my crontab: 44 4 * * * cd /opt/src ; /usr/bin/xvfb-run ./R-devel-Doit 44 5 * * * cd /opt/src ; /usr/bin/xvfb-run ./R-aop-Doit Those scripts update the R source code (development version and patched version of latest official release, respectively), compile and install from scratch, and update (if necessary) some contributed packages for these R versions. I do not have iplots among the packages that I need for these R versions, but I do have rgl and upgrades of the latter work (without the need of setting any environment variables). Hope this helps. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Feature or bug?
G'day Sigbert, long time no see :) How is Berlin these days? On Thu, 21 May 2015 11:45:26 +0200 Sigbert Klinke wrote: It is a feature. > if I run > > update <- function (newtime) { ginput <<- list(time=newtime)} > > server <- function (input) { > print(paste("Before", input$time)) > update(1) > print(paste("After:", input$time)) > } > > ginput <- list(time=0) > server(ginput) > > then I get as result > > [1] "Before 0" > [1] "After: 0" The first print command evaluates input and after this the function server has an object named "input" in its local environment. The second print command reuses this object and extracts the component time from it (which has not changed). The change of the global variable has no effect. > If I uncomment the first print > > update <- function (newtime) { ginput <<- list(time=newtime) } > > server <- function (input) { > #print(paste("Before", input$time)) > update(1) > print(paste("After:", input$time)) > } > > ginput <- list(time=0) > server(ginput) > > then I get > > [1] "After: 1" Because the global variable is changed before input is evaluated. R has lazy argument evaluation, arguments are only evaluated once they are needed. You are essentially getting bitten by R's lazy evaluation plus "pass by value" syntax. > Even when I use a side effect (by assign some new value to a global > variable) I would have expected the same behaviour in both cases. To get the behaviour that you expect, you would have to write your code along the following lines: R> update <- function (newtime) { ginput <<- list(time=newtime)} R> server <- function(input){ + inp <- as.name(deparse(substitute(input))) + print(paste("Before", eval(substitute(XXX$time, list(XXX=inp) + update(1) + print(paste("After:", eval(substitute(XXX$time, list(XXX=inp) + } R> ginput <- list(time=0) R> server(ginput) [1] "Before 0" [1] "After: 1" A cleaner way is perhaps to use environments, as these are passed by reference: R> update <- function(env, newtime) env$time <- newtime R> server <- function(input){ + print(paste("Before", input$time)) + update(input, 1) + print(paste("After:", input$time)) + } R> ginput <- new.env() R> ginput$time <- 0 R> server(ginput) [1] "Before 0" [1] "After: 1" HTH. Cheers, Berwin == Full address A/Prof Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009 e-mail: berwin.turl...@gmail.com Australiahttp://www.maths.uwa.edu.au/~berwin http://www.researcherid.com/rid/A-4995-2008 __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] cbind question, please
G'day Erin, On Thu, 23 Apr 2015 20:51:18 -0400 Erin Hodgess wrote: > Here is the big picture. I have a character vector with all of the > names of the variables in it. > > I want to "cbind" all of the variables to create a matrix. > > Doing 3 is straightforward, but many, not so much. So I guess you want something like: R> do.call(cbind, sapply(big.char, as.name)) dog cat tree [1,] 1 25 [2,] 2 36 [3,] 3 47 Which does not seem to have the problem of confusing the numeric vector `cat' with the function 'cat'. HTH. Cheers, Berwin __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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] Where did lost variables go
G'day David, On Mon, 30 Dec 2013 20:42:53 -0500 David Parkhurst wrote: Some wild guesses in the absence of a reproducible example. > I have several variables in a data frame that aren't listed by ls() > after I attach that data frame. ls() list the objects in the global environment. If you attach a data frame it is attached to the search path, typically after the global environment. Type 'search()' to see your search path. ls() list the global environment, the first entry in the list and called ".GlobalEnv". Your data frame should be listed as an object in that environment. Assuming the name of your data frame is 'foo', then there should be the name 'foo' somewhere in the list of names returned by 'search()'. Assuming 'foo' is listed in the second position, then 'ls(2)' should list all the objects found at that location of the search path, i.e. all the variables in your data frame. > Where did they go, See above. > and how can I stop the hidden ones from masking the local ones? Do you mean with "local ones" those in the global environment and by "hidden ones" those that you couldn't find? I.e. is there an object "bar" listed by 'ls()' but also an object "bar" listed by 'ls(2)' (i.e. your data frame 'foo' contained a variable with name 'bar')? Then it is the other way round, the local ones are hiding the hidden ones. For that reason attaching data frames has its dangers. It allows to easily access the variables in the data frame, but any changes to a variable creates a local copy. Thus, any change *will* not propagate back to the data frame! Hopefully the commands below will clarify further. Cheers, Berwin R> foo <- data.frame(bar=rnorm(2), fubar=runif(2)) R> ls() [1] "foo" R> attach(foo) R> search() [1] ".GlobalEnv""foo" "package:stats" [4] "package:graphics" "package:grDevices" "package:utils" [7] "package:datasets" "package:methods" "Autoloads" [10] "package:base" R> ls(2) [1] "bar" "fubar" R> bar [1] -0.07741633 1.05804653 R> fubar [1] 0.08516929 0.82718383 R> bar <- "what now" R> ls() [1] "bar" "foo" R> bar [1] "what now" R> ls(2) [1] "bar" "fubar" R> get("bar", pos=2) [1] -0.07741633 1.05804653 R> foo bar fubar 1 -0.07741633 0.08516929 2 1.05804653 0.82718383 R> detach(2) R> bar [1] "what now" R> fubar Error: object 'fubar' not found R> foo bar fubar 1 -0.07741633 0.08516929 2 1.05804653 0.82718383 R> attach(foo) The following object is masked _by_ .GlobalEnv: bar R> bar [1] "what now" R> fubar [1] 0.08516929 0.82718383 R> detach(2) R> bar [1] "what now" R> fubar Error: object 'fubar' not found R> foo bar fubar 1 -0.07741633 0.08516929 2 1.05804653 0.82718383 __ R-help@r-project.org 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] replicating C example from the Extensions Manual problem
G'day Erin, On Thu, 15 Mar 2012 01:03:59 -0500 Erin Hodgess wrote: > What is wrong, please? Missing #include #include In particular the latter is defining R_len_t. Guess that code snippet in the manual is only showing the code of the function, but not the complete file that needs to be compiled. :) HTH. Cheers, Berwin __ R-help@r-project.org 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] Argument validation within functions
G'day Johannes, On Tue, 06 Dec 2011 16:15:21 +0100 "Johannes Radinger" wrote: > Thank you, i didn't know that the !operator is > also working for is.numeric etc. > > Anyway I want to test if an argument is set in the > function call and if not a code is executed... So > far I tried: > > f <-function(a,b){ > if(!exists("b")) print("exists: b is not set") > if(is.null("b")) print("is.null : b is not set") > } > > f(a=1,b=2) > f(a=1) > f(b=2) > > I don't really know how to do it...e.g: for f(a=1) b is not set > so it also can't be NULL (thats why is.null is not working). I > just want to "test" if it is set with the function call not outside > the function or before etc. ?missing Cheers, Berwin == Full address A/Prof Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: berwin.turl...@gmail.com Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] expand.gird with constraints?
G'day Kathie, On Tue, 2 Aug 2011 19:56:00 -0700 (PDT) Kathie wrote: > Hi, R users, > > Here is an example. > > k <- c(1,2,3,4,5) > i <- c(0,1,3,2,1) > > if k=1, then j=0 from i > if k=2, then j=0, 1 from i > if k=3, then j=0, 1, 2, 3 from i > if k=4, then j=0, 1, 2 from i > if k=5, then j=0, 1 from i > > so i'd like to create a list like below. > > > list >k j > 1 1 0 > 2 2 0 > 3 2 1 > 4 3 0 > 5 3 1 > 6 3 2 > 7 3 3 > 8 4 0 > 9 4 1 > 10 4 2 > 11 5 0 > 12 5 1 > > I tried expand.grid, but I can't. > > Any suggestion will be greatly appreciated. One possibility is: R> k <- c(1,2,3,4,5) R> i <- c(0,1,3,2,1) R> tt <- c(1, 2, 4, 3, 2) R> data.frame(k=rep(k, tt), j=unlist(sapply(tt, function(ii) i[1:ii]))) k j 1 1 0 2 2 0 3 2 1 4 3 0 5 3 1 6 3 3 7 3 2 8 4 0 9 4 1 10 4 3 11 5 0 12 5 1 Not sure whether this is generalisable to your real problem... HTH. Cheers, Berwin == Full address A/Prof Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009 e-mail: berwin.turl...@gmail.com Australiahttp://www.maths.uwa.edu.au/~berwin http://www.researcherid.com/rid/A-4995-2008 __ R-help@r-project.org 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] Standards for delivery of GPL software in CRAN packages
G'day Gavin, On Mon, 27 Jun 2011 16:36:57 +0100 Gavin Simpson wrote: [...] > I don't recall the GPL mentioning anything requiring that the source > code be "helpful". The authors have most certainly fulfilled their > requirements under GPL, as has CRAN in distributing the package > sources. > > Or am I being obtuse and completely missing your point? Maybe :) While I agree with the first sentence cited, the GPL v2 says: The source code for a work means the preferred form of the work for making modifications to it. and GPL v3: The "source code" for a work means the preferred form of the work for making modifications to it. The OP perhaps suspects, though without offering any proof, that the authors of the package might not keep all the code in a single file but distributed over several files, perhaps with comments. If they use a script to remove all comments and collect all files into a single file before creating the .tar.gz file, then there would be a question whether the GPL is adhered too. Though, in the example provided, there would be no way to find out whether this is the case without contacting the authors. If memory serves correctly, it has definitely been suggested to me that putting into the src/ directory files from a tangled literate programming file without distributing this file too in the package might be questionable under the GPL. In this case, if requested, the author of such a package would presumably have to provide the file that contains the literate program from which the sources were tangled. As another example, as one could gather from occasional comments on r-devel (and r-help?), there must have been quite some discussions among the R core members whether it is compatible with the GPL to have a vignette in a package that depends on private style files (that are not distributed) and cannot be processed by anybody but the package author. Of course, since "R CMD check" now builds vignettes, and a package is only accepted on CRAN if it passes "R CMD check" on Kurt's machine, this discussion about how the GPL applies is rather academic. Cheers, Berwin __ R-help@r-project.org 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] %in% operator - NOT IN
G'day Dan, On Sun, 8 May 2011 05:06:27 -0400 Dan Abner wrote: > Hello everyone, > > I am attempting to use the %in% operator with the ! to produce a NOT > IN type of operation. Why does this not work? Suggestions? > > > data2[data1$char1 %in% c("string1","string2"),1]<-min(data1$x1) > > data2[data1$char1 ! %in% > > c("string1","string2"),1]<-max(data1$x1)+1000 > > Error: unexpected '!' in "data2[data1$char1 !" Try (untested) R> data2[!(data1$char1 %in% c("string1","string2")),1]<-max(data1$x1)+1000 HTH. Cheers, Berwin __ R-help@r-project.org 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] Using $ accessor in GAM formula
G'day Rolf, On Fri, 06 May 2011 09:58:50 +1200 Rolf Turner wrote: > but it's strange that the dodgey code throws an error with gam(dat1$y > ~ s(dat1$x)) but not with gam(dat2$cf ~ s(dat2$s)) > Something a bit subtle is going on; it would be nice to be able to > understand it. Well, R> traceback() 3: eval(expr, envir, enclos) 2: eval(inp, data, parent.frame()) 1: gam(dat$y ~ s(dat$x)) So the lines leading up to the problem seem to be the following from the gam() function: vars <- all.vars(gp$fake.formula[-2]) inp <- parse(text = paste("list(", paste(vars, collapse = ","), ")")) if (!is.list(data) && !is.data.frame(data)) data <- as.data.frame(data) Setting R> options(error=recover) running the code until the error occurs, and then examining the frame number for the gam() call shows that "inp" is "expression(list( dat1,x ))" in your first example and "expression(list( dat2,s ))" in your second example. In both examples, "data" is "list()" (not unsurprisingly). When, dl <- eval(inp, data, parent.frame()) is executed, it tries to eval "inp", in both cases "dat1" and "dat2" are found, obviously, in the parent frame. In your first example "x" is (typically) not found and an error is thrown, in your second example an object with name "s" is found in "package:mgcv" and the call to eval succeeds. "dl" becomes a list with two components, the first being, respectively, "dat1" or "dat2", and the second the body of the function "s". (To verify that, you should probably issue the command "debug(gam)" and step through those first few lines of the function until you reach the above command.) The corollary is that you can use the name of any object that R will find in the parent frame, if it is another data set, then that data set will become the second component of "inp". E.g.: R> dat=data.frame(min=1:100,cf=sin(1:100/50)+rnorm(100,0,.05)) R> gam(dat$cf ~ s(dat$min)) Family: gaussian Link function: identity Formula: dat$cf ~ s(dat$min) Estimated degrees of freedom: 3.8925 total = 4.892488 GCV score: 0.002704789 Or R> dat=data.frame(BOD=1:100,cf=sin(1:100/50)+rnorm(100,0,.05)) R> gam(dat$cf ~ s(dat$BOD)) Family: gaussian Link function: identity Formula: dat$cf ~ s(dat$BOD) Estimated degrees of freedom: 3.9393 total = 4.939297 GCV score: 0.002666985 > Just out of pure academic interest. :-) Hope your academic curiosity is now satisfied. :) HTH. Cheers, Berwin == Full address A/Prof Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: berwin.turl...@gmail.com Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] Strange R squared, possible error
G'day Gabor, On Thu, 17 Mar 2011 20:38:21 -0400 Gabor Grothendieck wrote: > > Or am I missing something? O.k., because the residuals don't add to zero, there may be a non-zero correlation between residuals and fitted values, which messes up the equation at the variance level. > Try it on an example to convince yourself: > > > fm <- lm(demand ~ Time, BOD) > > var(fitted(fm)) + var(resid(fm)) - var(BOD$demand) > [1] 3.552714e-15 > > > > fm0 <- lm(demand ~ Time - 1, BOD) > > var(fitted(fm0)) + var(resid(fm0)) - var(BOD$demand) > [1] 59.28969 But, and this is of course the geometry of least squares: R> sum(fitted(fm)^2) + sum(resid(fm)^2) - sum(BOD$demand^2) [1] 0 R> sum(fitted(fm0)^2) + sum(resid(fm0)^2) - sum(BOD$demand^2) [1] 2.273737e-13 and the reason why the formula changes if there is no (explicit) intercept term in the model. Cheers, Berwin __ R-help@r-project.org 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] Strange R squared, possible error
G'day Gabor, On Thu, 17 Mar 2011 11:36:38 -0400 Gabor Grothendieck wrote: > The idea is that if you have a positive quantity that can be broken > down into two nonnegative quantities: X = X1 + X2 then it makes sense > to ask what proportion X1 is of X. For example: 10 = 6 + 4 and 6 is > .6 of the total. > > Now, in the case of a model with an intercept its a mathematical fact > that the variance of the response equals the variance of the fitted > model plus the variance of the residuals. Thus it makes sense to ask > what fraction of the variance of the response is represented by the > variance of the fitted model (this fraction is R^2). > > But if there is no intercept then that mathematical fact breaks down. > That is, its no longer true that the variance of the response equals > the variance of the fitted model plus the variance of the residuals. [...] Do you have any reference to back this up, as I am somewhat surprised. I know that if there is no intercept in the model, then the residuals may not add to zero, but we are still doing least-squares, i.e. the fitted values are orthogonal to the residual vector and the variance of the response is the sum of the variance of the fitted model plus the variance of the residuals. Or am I missing something? Cheers, Berwin == Full address ==== A/Prof Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] identical values not so identical? newbie help please!
G'day Maja, On Thu, 10 Mar 2011 11:44:28 -0800 (PST) maiya wrote: > Aaah, it truly is wonderful, this technology! > I guess I'm going to have to override it a bit though.. > Along the lines of > > tae <- ifesle(all.equal(obs, exp) == TRUE, 0, sum(abs(obs - exp))) Please read the help page on all.equal(). As all.equal() tries to report the difference between the objects, if one exists, one should not compare the result from all.equal directly to TRUE. HTH. Cheers, Berwin == Full address ======== A/Prof Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] two apparent anomalies
On Sat, 22 Jan 2011 06:16:43 -0800 (PST) "analys...@hotmail.com" wrote: > (1) > > > a = c("a","b") > > mode(a) > [1] "character" > > b = c(1,2) > > mode(b) > [1] "numeric" > > c = data.frame(a,b) > > mode(c$a) > [1] "numeric" R> str(c) 'data.frame': 2 obs. of 2 variables: $ a: Factor w/ 2 levels "a","b": 1 2 $ b: num 1 2 Character vectors are turned into factors by default by data.frame(). OTOH: R> c = data.frame(a,b, stringsAsFactors=FALSE) R> mode(c$a) [1] "character" > (2) > > > > a = c("a","a","b","b","c") > > levels(as.factor(a)) > [1] "a" "b" "c" > > levels(as.factor(a[1:3])) > [1] "a" "b" > > a = as.factor(a) > > levels(a) > [1] "a" "b" "c" > > levels(a[1:3]) > [1] "a" "b" "c" Subsetting factors does not get rid of no-longer used levels by default. OTOH: R> levels(a[1:3, drop=TRUE]) [1] "a" "b" or R> levels(factor(a[1:3])) [1] "a" "b" HTH. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] expand.grid
G'day Nick, On Wed, 19 Jan 2011 09:43:56 +0100 "Nick Sabbe" wrote: > Given a dataframe > > dfr<-data.frame(c1=c("a", "b", NA, "a", "a"), c2=c("d", NA, "d", "e", > "e"), c3=c("g", "h", "i", "j", "k")) > > I would like to have a dataframe with all (unique) combinations of > all the factors present. Easy: R> expand.grid(lapply(dfr, levels)) c1 c2 c3 1 a d g 2 b d g 3 a e g 4 b e g 5 a d h 6 b d h 7 a e h 8 b e h 9 a d i 10 b d i 11 a e i 12 b e i 13 a d j 14 b d j 15 a e j 16 b e j 17 a d k 18 b d k 19 a e k 20 b e k > In fact, I would like a simple solution for these two cases: given > the three factor columns above, I would like both all _possible_ > combinations of the factor levels, and all _present_ combinations of > the factor levels (e.g. if I would do this for the first 4 rows of > dfr, it would contain no combinations with c3="k"). R> dfrpart <- lapply(dfr[1:4,], factor) R> expand.grid(lapply(dfrpart, levels)) c1 c2 c3 1 a d g 2 b d g 3 a e g 4 b e g 5 a d h 6 b d h 7 a e h 8 b e h 9 a d i 10 b d i 11 a e i 12 b e i 13 a d j 14 b d j 15 a e j 16 b e j > It would also be nice to be able to choose whether or not NA's are > included. R> expand.grid(lapply(dfrpart, function(x) c(levels(x), + if(any(is.na(x))) NA else NULL))) c1 c2 c3 1 ad g 2 bd g 3 d g 4 a e g 5 be g 6 e g 7 a g 8 b g 9 g 10ad h 11bd h HTH. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] Retrieving Factors with Levels Ordered
G'day H.T. On Sat, 1 Jan 2011 00:41:10 -0500 (EST) "H. T. Reynolds" wrote: > When I create a factor with labels in the order I want, write the > data as a text file, and then retrieve them, the factor levels are no > longer in the proper order. Not surprisingly. :) [..big snip..] > testdf <- read.table(file = 'Junkarea/test.txt') Did you look at the file Junkarea/test.txt, e.g. with a text editor? You will see that read.table() stores only the observed values but no information about what the mode of each variable in the data frame is. In particular, it doesn't store that a variable is a factor and definitely not the levels and their ordering. Actually, write.table() saves a factor by writing out the observed labels as character stings. Only because read.table() by default turns character data into factors (a behaviour that some useRs don't like and why the option stringsAsFactors exists) you end up with a factor again. > # But levels are no longer ordered But the help file of factor (see ?factor) states in the warning section: The levels of a factor are by default sorted, but the sort order may well depend on the locale at the time of creation, and should not be assumed to be ASCII. Thus, arguably, the levels are ordered, just not the order you want. :) > Clearly I am missing something obvious, but I can't see it. If I save > "feducord" and load it, the order of the levels is as it should be. > But I don't know why writing to a test file should change anything. > Any help would be greatly appreciated. If you save() and load() an object, then it is saved in binary format, thus much more information about it can be stored. Indeed, I would expect that a faithful internal representation of the object is stored in the binary format so that a save(), rm() and load() would restore exactly the same object. Saving objects in text formats is prone to loss of information. E.g., as you experience with write.table() and read.table() no information about ordering of levels is stored using these functions. If care is not taken when writing numbers to text files, the internal representation can change, e.g. : R> x <- data.frame(x=seq(from=0, to=1, length=11)) R> write.table(x, file="/tmp/junk1") R> y <- read.table("/tmp/junk1") R> identical(x,y) [1] FALSE R> x-y x 1 0.00e+00 2 0.00e+00 3 0.00e+00 4 5.551115e-17 5 0.00e+00 6 0.00e+00 7 1.110223e-16 8 1.110223e-16 9 0.00e+00 10 0.00e+00 11 0.00e+00 Having said all this, if you want to save your data in a text file with a representation that remembers the ordering of the factor levels, look at dput(): R> fac <- gl(2,4, labels=c("White", "Black")) R> fac [1] White White White White Black Black Black Black Levels: White Black R> write.table(fac, file="/tmp/junk") R> str(read.table("/tmp/junk")) 'data.frame': 8 obs. of 1 variable: $ x: Factor w/ 2 levels "Black","White": 2 2 2 2 1 1 1 1 R> dput(fac, file="/tmp/junk") R> str(dget(file="/tmp/junk")) Factor w/ 2 levels "White","Black": 1 1 1 1 2 2 2 2 It might just be that the text representation used by dput() is not particularly digestible for the human eye. :) > (You're right, I don't have anything better to do on New Year's eve.) New Year's eve? The first day of the new year is already nearly over! :) HTH. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] using ``<-'' in function argument
G'day Ivan, On Fri, 03 Dec 2010 10:54:58 +0100 Ivan Calandra wrote: > Arf, yes it makes sense now! Well, my original post said: "R has lazy evaluation" and "the assignment takes place when the function evaluates the argument" :) > So the idea here is: never use "<-" in function argument... Never say never. :) But if you want to perform a (global) assignment within an actual argument to a function, you should better be sure that the formal argument gets evaluated. In fact, this trick is probably a good way, short of studying the source, to find out whether an argument gets evaluated or not (or only occasionally). Cheers, Berwin __ R-help@r-project.org 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] using ``<-'' in function argument
On Thu, 2 Dec 2010 23:34:02 -0500 David Winsemius wrote: > [...] Erik is telling you that your use of ncol<-4 got evaluated to > 4 and that the name of the resulting object was ignored, howevert the > value of the operation was passed on to matrix which used positional > matching since "=" was not used. Sounds like a fair summary of what Erik said, but it is subtly wrong. R has lazy evaluation of its arguments. There is nothing that forces the assignment to be evaluated and to pass the result into the function. On the contrary, the assignment takes place when the function evaluates the argument. For example: R> rm(list=ls(all=TRUE)) R> ls() character(0) R> foo <- function(x, y){ + if (x > 3) cat(y, "\n") + x} R> foo(4, bar <- 5) 5 [1] 4 R> ls() [1] "bar" "foo" R> bar [1] 5 R> foo(2, fubar <- 6) [1] 2 R> fubar Error: object 'fubar' not found R> ls() [1] "bar" "foo" > Usually the problem facing newbies is that they want to save > keystrokes and so use "=" for assignment (also a potential pitfall > although not as likely to mess you up as the choice to use the > two-keystroke path for argument assignment). On the contrary, the opposite is also very likely. One of my favourite idioms is: plot(fm <- lm(y~x, data=some.data)) to (1) fit a model, (2) assign the fitted model to an object and (3) look immediately at diagnostic plots. Students came to me and said that the code in the lab sheet didn't work and they were getting strange error messages "about objects not being found". They reassured me that they had typed in exactly what was on the lab sheet. Of course, once I got to their computer and looked at their screen, it was clear that they had typed: plot(fm = lm(y~x, data=some.data)) Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] DBLEPR?
G'day John, On Tue, 16 Nov 2010 14:02:57 -0500 "Prof. John C Nash" wrote: > Are the xxxPR routines now deprecated (particularly for 64 bit > systems) or still OK to use? They are still OK to use, and I use them occasionally. > If OK, can anyone point to documentation and examples? Section 6.5.1 "Printing from Fortran" in the "Writing R Extensions" manual has documentation (but no examples). Luckily, Section 5.7 "Fortran I/O", the section that I always look at first, has a link to Section 6.5.1. :) Cheers, Berwin ========== Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] R version 2-12.0 - running as 32 or as 64 bit?
G'day Dimitri, On Fri, 29 Oct 2010 16:45:00 -0400 Dimitri Liakhovitski wrote: > Question: I installed R verison 2-12.0 on my Windows 7 (64 bit) PC. > When I was installing it, it did not ask me anything about 32 vs. 64 > bit. So, if I run R now - is it running as a 32-bit or a 64-bit? Well, when I did the same, I got two shortcuts installed on my desktop, one named R 2.12.0 and the other named R x64 2.12.0. If I had any doubts which version of R these shortcuts would start, then the fourth line of the start-up message would put them to rest. Missing that message, you can always issue the command > .Machine$sizeof.pointer if the answer is 4, you are running 32bit, if the answer is 8, then you are running 64 bit. HTH. Cheers, Berwin == Full address ======== Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] what´s wrong with this code?
G'day Jose, On Fri, 29 Oct 2010 11:25:05 +0200 José Manuel Gavilán Ruiz wrote: > Hello, I want to maximize a likelihood function expressed as an > integral that can not be symbolically evaluated. I expose my problem > in a reduced form. > > g<- function(x){ > integrand<-function(y) {exp(-x^2)*y} > g<-integrate(integrand,0,1) > } > h<-function(x) log((g(x))) > > g is an object of the class function, but g(2) is a integrate object, > I can print(g(2)) an get a result, but > if I try h(2) R says is a nonnumeric argument for a mathematical > function. My goal is to maximize h. Which can be done without R quite trivially. The result of g(x) is exp(-x^2)/2. So h(x) is -x^2-log(2), which is maximised at x=0. Cheers, Berwin __ R-help@r-project.org 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] what´s wrong with this code?
G'day Jose, On Fri, 29 Oct 2010 11:25:05 +0200 José Manuel Gavilán Ruiz wrote: > Hello, I want to maximize a likelihood function expressed as an > integral that can not be symbolically evaluated. I expose my problem > in a reduced form. > > g<- function(x){ > integrand<-function(y) {exp(-x^2)*y} > g<-integrate(integrand,0,1) > } > h<-function(x) log((g(x))) > > g is an object of the class function, but g(2) is a integrate object, > I can print(g(2)) an get a result, but > if I try h(2) R says is a nonnumeric argument for a mathematical > function. R> print(g(2)) 0.00915782 with absolute error < 1.0e-16 Indeed print(g(2)) gives an output, but it is obviously not just a numeric number but something formatted, presumably based on the class of the returned object from g(2) and the values that this object contains. (You may want to read up on R's way(s) to object oriented programming). So what does g() return? R> str(g(2)) List of 5 $ value : num 0.00916 $ abs.error : num 1.02e-16 $ subdivisions: int 1 $ message : chr "OK" $ call: language integrate(f = integrand, lower = 0, upper = 1) - attr(*, "class")= chr "integrate" The object is a list with several values, and class "integrate". > My goal is to maximize h. what´s wrong? I guess you want to pass only the component value to h(). R> g <- function(x){ + integrand <- function(y) {exp(-x^2)*y} + integrate(integrand, 0, 1)$value} R> g(2) [1] 0.00915782 R> h(g(2)) [1] -0.693231 HTH, Cheers, Berwin __ R-help@r-project.org 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] cube root of a negative number
G'day Bill, On Wed, 27 Oct 2010 10:34:27 +1100 wrote: [...] > It is no surprise that this does not work when working in the real > domain, except "by fluke" with something like > > > -4^(1/3) > [1] -1.587401 > > > > where the precedence of the operators is not what you might expect. > Now that could be considered a bug, since apparently > > > -4^(1/2) > [1] -2 > > which comes as rather a surprise! [...] Mate, you must have been using Excel (or bc) too much recently if this takes you by surprise. :) This is exactly the behaviour that I would expect as I was always taught that exponentiation was the operator with the highest priority. The discussion at http://en.wikipedia.org/wiki/Order_of_operations points out that "there exist differing conventions concerning the unary operator -" but gives only Excel and bc as examples of programs/languages who give the unary operator the higher precedence. Moreover, http://www.burns-stat.com/pages/Tutor/spreadsheet_addiction.html points out that in Excel "the precedence of unary minus is different than many other languages (including VBA which is often used with Excel)". Cheers, Berwin __ R-help@r-project.org 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] cube root of a negative number
G'day Gregory, On Tue, 26 Oct 2010 19:05:03 -0400 Gregory Ryslik wrote: > Hi, > > This might be me missing something painfully obvious but why does the > cube root of the following produce an NaN? > > > (-4)^(1/3) > [1] NaN 1/3 is not exactly representable as a binary number. My guess is that the number that is closest to 1/3 and representable cannot be used as the exponent for negative numbers, hence the NaN. Essentially, don't expect finite precision arithmetic to behave like infinite precision arithmetic, it just doesn't. The resources mentioned in FAQ 7.31 can probably shed more light on this issue. Cheers, Berwin == Full address ==== Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] R-Fortran question (multiple subroutines)
G'day all, On Mon, 25 Oct 2010 06:52:15 -0400 Duncan Murdoch wrote: > Remko Duursma wrote: [...] > > I.e, my code looks something like this: > > > > subroutine f(x,y,z) > > > > call g(x,y,z) > > > > end > > > > subroutine g(x,y,z) > > > > z = x*y > > > > end > > > > > > calling this from R shows that subroutine g is not called. The code > > compiled as executable works fine. > > There are no such limitations imposed by R. I'd suggest your > diagnosis of the problem is wrong. If you can't spot the problem, > please post a real example (simplified if possible, but not as much > as the one above). Actually, it turns out that this example is simplified enough. :) I put this snippet into a file, compiled it via "R CMD SHLIB", loaded it into R and then was very surprised about the result of ".Fortran("f", x=1.1, y=2.2, z=0.0)". Eventually it dawned to me, Remko needs to read the table in section 5.2 of "Writing R Extensions". The simple fix in this case is to put a "double precision x, y, z" at the beginning of each subroutine. The default precision in FORTRAN is not double precision, that's why the code seems to work when compiled as stand alone but not when called from R. In general, I would recommend to start each subroutine with "implicit none" (o.k., I know that this is not standard FORTRAN77 but most compiler support that construct; the alternative, that would confrom with the FORTRAN77 standard, is to declare everything to be implicitly of type character and then let the fun begin) and then to explicitly declare all variables to be of the appropriate type. HTH. Cheers, Berwin __ R-help@r-project.org 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] R-Fortran question (multiple subroutines)
G'day Remko, On Mon, 25 Oct 2010 15:33:30 +1100 Remko Duursma wrote: > apologies if this is somewhere in a manual, I have not been able to > find anything relevant. You probably have to set some appropriate flags and the information should be somewhere in the manual of your compiler. Though, "Writing R Exensions" will tell you now to change/specify flags for compilers. > I run Windows Vista. My condolences. :) > But is it possible to have more than one subroutine in my source file, > one depending on the other? Yes. > I.e, my code looks something like this: > > subroutine f(x,y,z) > > call g(x,y,z) > > end > > subroutine g(x,y,z) > > z = x*y > > end > > calling this from R shows that subroutine g is not called. What do you mean with "g is not called"? How did you assert this? If the code functions correctly, then g has to be called, either explicitly or implicitly (e.g. if it was inlined). Do you mean that when you compile the code and load the resulting DLL from R, you can only call subroutine f via .Fortran but not subroutine g? > The code compiled as executable works fine. What do you mean with "compiled as executable"? The snippet that you showed would not compile as an executable as there is no main program. So what do you mean with "works fine"? That you can call the subroutine g from the main program that you use when creating an executable? My guess, based on this snippet, is that your compiler notices that subroutine g is not really needed and replaces the call to g within f by the body of g and does not create a separate entry point for g in the compiled object; a process known as in-lining (IIRC) and typically used by compilers if high levels of optimisations are requested. The compiler does not know that you intend to call g later directly from R via .Fortran, all it sees is the code in the file and, if high optimisation is requested, it may feel free to rearrange the code to stream-line it (e.g. by in-lining). So possible solutions could be: 1) ask for a lower level of optimisation 2) tell the compiler not to do in-lining (flag --no-inline??) 3) put the routines into separate files, compile the files separately and then link all the resulting object files together into a DLL. AFAIK, optimising/inlining across separate files is a tricky issue so few, if any, compilers would do so. 4) Check whether there is an option to specify which exportable symbols have to be in the DLL in the end. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] creating 'all' sum contrasts
G'day Michael, On Fri, 15 Oct 2010 12:09:07 +0100 Michael Hopkins wrote: > OK, my last question didn't get any replies so I am going to try and > ask a different way. > > When I generate contrasts with contr.sum() for a 3 level categorical > variable I get the 2 orthogonal contrasts: > > > contr.sum( c(1,2,3) ) > [,1] [,2] > 110 > 201 > 3 -1 -1 These two contrasts are *not* orthogonal. > This provides the contrasts <1-3> and <2-3> as expected. But I also > want it to create <1-2> (i.e. <1-3> - <2-3>). So in general I want > all possible orthogonal contrasts - think of it as the contrasts for > all pairwise comparisons between the levels. You have to decide what you want. The contrasts for all pairwise comparaisons between the levels or all possible orthogonal contrasts? The latter is actually not well defined. For a factor with p levels, there would be p orthogonal contrasts, which are only identifiable up to rotation, hence infinitely many such sets. But there are p(p-1) pairwise comparisons. So unless p=2, yo have to decide what you want > Are there are any options for contrast() or other functions/libraries > that will allow me to do this automatically? Look at package multcomp, in particular functions glht and mcp, these might help. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] declaring GPL license
G'day Stacey, On Thu, 14 Oct 2010 12:36:20 -0400 Stacey Wood wrote: > Thanks everyone. Now I know that I should not include another copy > of the license, but where should I refer to the copies on > http://www.r-project.org/Licenses/? In the DESCRIPTION file? I used to mention that location in the License: line of the DESCRIPTION file, after GPL (>=2). But some versions ago, that line got standardised and now it is no longer possible to put additional text onto that line. I guess you could just add a line LicenceLocation: As far as I understand the documentation, you are free to add additional lines beyond those described in "Writing R Extensions" to the DESCRIPTION file. Some entries described in the manual are optional and for those that are not some restrictions may apply that have to be followed. Personally, I decided not to put any such reference into the DESCRIPTION file of my packages as I guess that people who worry about the lisence will either look into the source package (where I have a copy) or know how to find one. :) HTH. Cheers, Berwin __ R-help@r-project.org 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] declaring GPL license
G'day Marc, On Thu, 14 Oct 2010 10:46:39 -0500 Marc Schwartz wrote: > If you want (and you should), create and include a file called > "COPYING" in the 'inst' folder in the package, so that it gets copied > to the main package directory upon installation. [...] But that would be against the wishes of R core, from the "Writing R Extension" manual: Whereas you should feel free to include a license file in your source distribution, please do not arrange to @emph{install} yet another copy of the @acronym{GNU} @file{COPYING} or @file{COPYING.LIB} files but refer to the copies on @url{http://www.r-project.org/Licenses/} and included in the @R{} distribution (in directory @file{share/licenses}). :) Cheers, Berwin __ R-help@r-project.org 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] exponentiate elements of a whole matrix
On Wed, 13 Oct 2010 11:51:39 +0100 "Maas James Dr (MED)" wrote: > I've tried hard to find a way to exponentiate each element of a whole > matrix such that if I start with A > > A = [ 2 3 > 2 4] > > I can get back B > > B = [ 7.38 20.08 > 7.38 54.60] > > I've tried > > B <- exp(A) but no luck. What have you tried exactly? And with which version? This should work with all R versions that I am familiar with, e.g.: R> A <- matrix(c(2,2,3,4),2,2) R> A [,1] [,2] [1,]23 [2,]24 R> B <- exp(A) R> B [,1] [,2] [1,] 7.389056 20.08554 [2,] 7.389056 54.59815 Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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 to make list() return a list of *named* elements
On Mon, 4 Oct 2010 19:45:23 +0800 Berwin A Turlach wrote: Mmh, > R> my.return <- function (vector.of.variable.names) { > sapply(vector.of.variable.names, function(x) list(get(x))) > } make that: R> my.return <- function (vector.of.variable.names) { sapply(vector.of.variable.names, function(x) get(x)) } Cheers, Berwin __ R-help@r-project.org 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 to make list() return a list of *named* elements
G'day Hans, On Mon, 4 Oct 2010 11:28:15 +0200 Hans Ekbrand wrote: > On Thu, Sep 30, 2010 at 09:34:26AM -0300, Henrique Dallazuanna wrote: > > You should try: > > > > eapply(.GlobalEnv, I)[c('b', 'my.c')] > > Great! > > b <- c(22.4, 12.2, 10.9, 8.5, 9.2) > my.c <- sample.int(round(2*mean(b)), 4) > > my.return <- function (vector.of.variable.names) { > eapply(.GlobalEnv, I)[vector.of.variable.names] > } Well, if you are willing to create a vector with the variable names, then simpler solutions should be possible, i.e. solutions that only operate on the objects of interest and not on all objects in the global environment (which could be a lot depending on your style). For example: R> my.return <- function (vector.of.variable.names) { sapply(vector.of.variable.names, function(x) list(get(x))) } R> str(my.return(c("b","my.c"))) List of 2 $ b : num [1:5] 22.4 12.2 10.9 8.5 9.2 $ my.c: int [1:4] 7 5 23 4 Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] R Spracheinstellung
G'day Tobias, On Wed, 29 Sep 2010 14:01:10 + "Keller Tobias (kelleto0)" wrote: > Für unseren Statistik Unterricht brauchen wir das R-Programm. > Ich habe das Programm heruntergeladen, deutsch gewählt und > installiert. Das Problem ist nach 3mal neu installieren, dass das > Programm immer auf englisch kommt. > > Können Sie mir helfen? R for Windows FAQ 3.2: http://cran.ms.unimelb.edu.au/bin/windows/base/rw-FAQ.html#I-want-R-in-English_0021 HTH. Cheers, Berwin __ R-help@r-project.org 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] OT: What distribution is this?
G'day Rainer, On Sun, 26 Sep 2010 10:29:08 +0200 Rainer M Krug wrote: > I realized that I am simply interested in the pdf p(y) that y > *number* of entities (which ones is irrelevant) in N are are *not* > drawn after the sampling process has been completed. Even simpler (I > guess), in a first step, I would only need the mean number of > expected non-drawn entities in N (pMean). > > But what about the range in between? You can calculate those probilities. Your problem was sounding vaguely familiar to me yesterday but I could not put my finger onto its abstract counterpart. Now with this clarification I can. :) Your set up is exactly the one of lottery draws. In each draw n of N numbers are drawn without replacement. So your question is "what is the probability that I have not seen k of the N numbers after x draws?". This question can easily be answered by using some Markov Chain theory. Let Y_l be the number of entities that you have not seen after the l-th draw, taking possible values 0, 1, , N. Obviously, Y_0=N with probability 1, and Y_1=N-n with probability 1. Now, the probability that Y_l, l>=1, is equal to k is given by; 0 if k=N-n+1, N-n+2,..., N and sum_{i=0,...,n} P[Y_l=k | Y_{l-1} = k+i] P[Y_{l-1} = k+i] o/w. P[Y_l=k | Y_{l-1} = k+i] is the probability that the n entities sampled in the l-th draw consists of i entities of the k+i entities that were not seen by draw l-1 and n-i entities of the N-k-i entities that were already seen by draw l-1. This probability is choose(k+i, i)*choose(N-k-i, n-i) / choose(N, n) Note that this probability is independent of l, i.e., Y_0, Y_1, Y_2,... form a stationary Markov Chain. Thus, to answer your question, the strategy would be to calculate the transition matrix once and then start with either the state vector of Y_0 or of Y_1 (both of which are rather trivial as mentioned above) and calculate the state vector of Y_x by repeatedly multiplying the chosen initial state vector with the transition matrix for the appropriate number of times. The transition matrix is (N+1)x(N+1), so it can be rather large, but the number of non-zero entries in the transition matrix is at most (N+1)*(n+1), so depending on the relative magnitude of n and N, sparse matrix techniques might be helpful (or using a language such as C, FOTRAN, ... for the calculations). HTH. Cheers, Berwin __ R-help@r-project.org 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] OT: What distribution is this?
G'day Rainer, On Sat, 25 Sep 2010 16:24:17 +0200 Rainer M Krug wrote: > This is OT, but I need it for my simulation in R. > > I have a special case for sampling with replacement: instead of > sampling once and replacing it immediately, I sample n times, and > then replace all n items. > > > So: > > N entities > x samples with replacement > each sample consists of n sub-samples WITHOUT replacement, which are > all replaced before the next sample is drawn > > My question is: which distribution can I use to describe how often > each entity of the N has been sampled? Surely, unless I am missing something, any given entity would have (marginally) a binomial distribution: A sub-sample of size n either contains the entity or it does not. The probability that a sub-sample contains the entity is a function of N and n alone. x sub-samples are drawn (with replacement), so the number of times that an entity has been sampled is the number of sub-samples in which it appears. This is given by the binomial distribution with parameters x and p, where p is the probability determined in the previous paragraph. I guess the fun starts if you try to determine the joint distribution of two (or more) entities. HTH. Cheers, Berwin ====== Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] minor diagonal in R
On Tue, 7 Sep 2010 06:26:09 -0700 (PDT) "Peng, C" wrote: > > This this what you want? > > > A=matrix(1:16,ncol=4) > > A > [,1] [,2] [,3] [,4] > [1,]159 13 > [2,]26 10 14 > [3,]37 11 15 > [4,]48 12 16 > > diag(A[1:4,4:1]) > [1] 13 10 7 4 Or > A[cbind(1:4,4:1)] [1] 13 10 7 4 one character more to type, but could be more efficient for large A :) Cheers, Berwin ========== Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] Call to rgamma using .C causes R to hang
G'day Steve, On Tue, 20 Jul 2010 17:20:49 +0930 Steve Pederson wrote: I suggest to insert the following two lines (untested as I usually don't work on Windows): > # include > # include > void test1 (double *x, double *result) > { GetRNGstate(); > result[0] = rgamma(*x, 2.0); PutRNGstate(); > } > [...] > I'm running Windows XP & I followed the instructions at > http://cran.r-project.org/doc/manuals/R-admin.html#Windows-standalone > after building R from source. But you didn't follow the instructions of Section 6.3 in the "Writing R Extensions manual. :) http://cran.r-project.org/doc/manuals/R-exts.html#Random-numbers HTH. Cheers, Berwin ========== Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] Namibia becoming NA
G'day Ted, On Sun, 18 Jul 2010 09:25:09 +0100 (BST) (Ted Harding) wrote: > On 18-Jul-10 05:47:03, Suresh Singh wrote: > > I have a data file in which one of the columns is country code and > > NA is the > > code for Namibia. > > When I read the data file using read.csv, NA for Namibia is being > > treated as > > null or "NA" > > > > How can I prevent this from happening? > > > > I tried the following but it didn't work > > input <- read.csv("padded.csv",header = TRUE,as.is = c("code2")) > > > > thanks, > > Suresh > > I suppose this was bound to happen, and in my view it represent > a bit of a mess! With a test file temp.csv: > > Code,Country > DE,Germany > IT,Italy > NA,Namibia > FR,France Thanks for providing an example. > leads to exactly the same result. And I have tried every variation > I can think of of "as.is" and "colClasses", still with exactly the > same result! Did you think of trying some variations of "na.strings"? ;-) IMO, the simplest way of coding missing values in CSV files is to have two consecutive commas; not some code (whether NA, 99, 999, -1, ...) between them. > Conclusion: If an entry in a data file is intended to become the > character value "NA", there seems to be no way of reading it in > directly. This should not be so: it should be preventable! It is, through simple use of the "na.strings" argument: R> X <- read.csv("temp.csv", na.strings="") R> X Code Country 1 DE Germany 2 IT Italy 3 NA Namibia 4 FR France R> which(is.na(X)) integer(0) HTH. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] forcing a zero level in contr.sum
G'day Stephen, On Wed, 7 Jul 2010 16:41:18 -0400 "Bond, Stephen" wrote: > Please, do not post if you do not know the answer. People will see > this has answers and skip. > > I tried with > mat1=contrasts(fixw$snconv) > mat1=mat1[,-2] > summary(frm2sum <- glm(resp.frm ~ > C(snconv,contr=mat1)+mprime+mshape,data=fixw,family="quasibinomial")) > > the unwanted level is still there. Unbelievable. model.matrix() instead of summary() would reveal that the estimate corresponding to "the unwanted level" is probably estimating something quite different then what you think it estimates. ?C and look at the how.many argument, presumably you want C(snconv, contr=mat1, how.many=NCOL(mat1)) To see what is going on, consider: R> ( tmp <- gl(4,3) ) [1] 1 1 1 2 2 2 3 3 3 4 4 4 Levels: 1 2 3 4 R> options(contrasts=c("contr.sum", "contr.poly")) R> ( tt <- contrasts(tmp) ) [,1] [,2] [,3] 1100 2010 3001 4 -1 -1 -1 R> tt <- tt[,-2] R> contrasts(tmp, 2) <- tt R> tmp [1] 1 1 1 2 2 2 3 3 3 4 4 4 attr(,"contrasts") [,1] [,2] 110 200 301 4 -1 -1 Levels: 1 2 3 4 R> contrasts(tmp) <- tt R> tmp [1] 1 1 1 2 2 2 3 3 3 4 4 4 attr(,"contrasts") [,1] [,2] [,3] 110 0.2886751 200 -0.8660254 301 0.2886751 4 -1 -1 0.2886751 Levels: 1 2 3 4 HTH. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] rgl installation failure
G'day, On Sat, 5 Jun 2010 12:51:08 +0530 "khush " wrote: > I am trying to install rgl package under R and getting some errors > which is below. > > > install.packages("rgl") > Warning in install.packages("rgl") : > argument 'lib' is missing: using '/usr/lib/R/library' > trying URL 'http://cran.csie.ntu.edu.tw/src/contrib/rgl_0.91.tar.gz' > Content type 'application/x-gzip' length 1677498 bytes (1.6 Mb) > opened URL > == > downloaded 1.6 Mb > > * installing *source* package _rgl_ ... > checking for gcc... gcc -m32 -std=gnu99 > checking for C compiler default output file name... a.out > checking whether the C compiler works... yes > checking whether we are cross compiling... no > checking for suffix of executables... > checking for suffix of object files... o > checking whether we are using the GNU C compiler... yes > checking whether gcc -m32 -std=gnu99 accepts -g... yes > checking for gcc -m32 -std=gnu99 option to accept ISO C89... none > needed checking how to run the C preprocessor... gcc -m32 -std=gnu99 > -E checking for gcc... (cached) gcc -m32 -std=gnu99 > checking whether we are using the GNU C compiler... (cached) yes > checking whether gcc -m32 -std=gnu99 accepts -g... (cached) yes > checking for gcc -m32 -std=gnu99 option to accept ISO C89... (cached) > none needed > checking for libpng-config... no > checking libpng... checking for grep that handles long lines and -e... > /bin/grep > checking for egrep... /bin/grep -E > checking for ANSI C header files... yes > checking for sys/types.h... yes > checking for sys/stat.h... yes > checking for stdlib.h... yes > checking for string.h... yes > checking for memory.h... yes > checking for strings.h... yes > checking for inttypes.h... yes > checking for stdint.h... yes > checking for unistd.h... yes > checking png.h usability... no > checking png.h presence... no > checking for png.h... no > checking for png_read_update_info in -lpng... no You should heed indications as the above, i.e. that items checked for are not found. [...] > In file included from pixmap.cpp:14: > pngpixmap.h:3:17: error: png.h: No such file or directory And this compiler error seems to be directly linked to the issues identified by configure; and all further error come from the fact that png.h cannot be found. You seem to be using some flavour of Unix, presumably GNU/Linux; but from your posting it is impossible to say which one. On a Debian based system (e.g. Kubuntu 9.10 which I am using), I would issue the command (assuming that apt-file is installed and initialised): ber...@bossiaea:~$ apt-file find include/png.h libpng12-dev: /usr/include/png.h and then check whether the libpng12-dev package is installed. If not, I would install it and retry. Perhaps that would be all that is needed to compile rgl, but there might be other packages missing. If you have a system based on RedHat and Suse, you will have to figure out how to use the package management tools on those systems to find out which packages have to be installed additionally so that you are able to compile rgl. HTH. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] Gelman 2006 half-Cauchy distribution
G'day Chris, On Fri, 28 May 2010 08:29:30 -0500 Christopher David Desjardins wrote: > Hi, > I am trying to recreate the right graph on page 524 of Gelman's 2006 > paper "Prior distributions for variance parameters in hierarchical > models" in Bayesian Analysis, 3, 515-533. I am only interested, > however, in recreating the portion of the graph for the overlain > prior density for the half-Cauchy with scale 25 and not the posterior > distribution. However, when I try: > > curve(dcauchy, from=0, to=200, location=0, scale=25) Which version of R do you use? This command creates 12 warnings under R 2.11.0 on my linux machine. Reading up on the help page of curve() would make you realise that you cannot pass the location and scale parameter to dcauchy in the manner you try. I guess you want: R> prior <- function(x) 2*dcauchy(x,location=0, scale=25) R> curve(prior, from=0, to=200) or, more compactly, R> curve(2*dcauchy(x, location=0, scale=25), from=0, to=200) Cheers, Berwin ========== Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] pdf files in loops
On Wed, 31 Mar 2010 22:06:48 -0400 James Rome wrote: > print() did not help, and I get strange messages about mode(onefile) > not being changed: Either: R> pic <- histogram(...) R> print(pic) or R> print(histogram(...)) > I would actually like all the histograms in one pdf file too... Then move the pdf() and the dev.off() commands outside the loop. HTH. Cheers, Berwin __ R-help@r-project.org 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] pdf files in loops
G'day James, On Wed, 31 Mar 2010 21:44:31 -0400 James Rome wrote: > I need to make a bunch of PDF files of histograms. [...] > What am I doing wrong? http://cran.ms.unimelb.edu.au/doc/FAQ/R-FAQ.html#Why-do-lattice_002ftrellis-graphics-not-work_003f HTH. Cheers, Berwin == Full address ======== Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] From THE R BOOK -> Warning: In eval(expr, envir, enclos) : non-integer #successes in a binomial glm!
G'day all, On Tue, 30 Mar 2010 16:19:46 +0100 Corrado wrote: > David Winsemius wrote: > > A) It is not an error, only a warning. Wouldn't it seem reasonable > > to issue such a warning if you have data that violates the > > distributional assumptions? > I am not questioning the approach. I am only trying to understand why > a (rather expensive) source of documentation and the behaviour of a > function are not aligned. 1) Also expensive books have typos in them. 2) glm() is from a package that is part of R and the author of this book is AFAIK not a member of R core, hence has no control on whether his documentation and the behaviour of a function are aligned. a) If he were documenting a function that was part of a package he wrote as support for his book, as some authors do, there might be a reason to complain. But then 1) would still apply. b) Even books written by members of R core have occasionally misalignments between the behaviour of a function and the documentation contained in such books. This can be due to them documenting a function over whose implementation they do not have control (e.g. a function in a contributed package) or the fact that R is improving/changing from version to version while books are rather static. For these reasons it is always worthwhile to check the errata page for a book, if such exists. The source of the warning is due to the fact that you do not provide all necessary information about your response. If your response is binomial (with a mean depended on some explanatory variables), then each response consists of two numbers, the number of trials and the number of success. If you calculate the observed proportion of successes from these two numbers and feed this into glm as the response, you are omitting necessary information. In this case, you should provide the number of trials on which each proportion is based as prior weights. For example: R> x <- seq(from=-1,to=1,length=41) R> px <- exp(x)/(1+exp(x)) R> nn <- sample(8:12, 41, replace=TRUE) R> yy <- rbinom(41, size=nn, prob=px) R> y <- yy/nn R> glm(y~x, family=binomial, weights=nn) Call: glm(formula = y ~ x, family = binomial, weights = nn) Coefficients: (Intercept)x 0.2461.124 Degrees of Freedom: 40 Total (i.e. Null); 39 Residual Null Deviance: 91.49 Residual Deviance: 50.83AIC: 157.6 R> glm(y~x, family=binomial) Call: glm(formula = y ~ x, family = binomial) Coefficients: (Intercept)x 0.2143 1.1152 Degrees of Freedom: 40 Total (i.e. Null); 39 Residual Null Deviance: 9.256 Residual Deviance: 5.229AIC: 49.87 Warning message: In eval(expr, envir, enclos) : non-integer #successes in a binomial glm! HTH, Cheers, Berwin ====== Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] language!!!
G'day Katharina, On Thu, 18 Mar 2010 11:16:26 +0100 Katharina Gerstendorf wrote: > there_s something that really bothers me about R and after hours and > hours of internet research, I_m still stuck with the same problem: And during these hours of research you did not come across the R Windows FAQ available on CRAN and all mirrors? E.g.: http://cran.ms.unimelb.edu.au/bin/windows/base/rw-FAQ.html particularly question 3.2 of that FAQ: http://cran.ms.unimelb.edu.au/bin/windows/base/rw-FAQ.html#I-want-R-in-English_0021 > Thanks a lot! You are welcome. Cheers, Berwin == Full address ======== Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] OT Sorta: Odds Are, It's Wrong...
G'day all, On Mon, 15 Mar 2010 12:46:42 -0500 Marc Schwartz wrote: > In turn, that reminds me of Stephen Senn's writing in Dicing with > Death: Chance, Risk and Health: > > "We can predict nothing with certainty but we can predict how > uncertain our predictions will be, on average that is." And I am reminded Stephen Senn's citation of Maurice Kendal in the same book: `Life', the statistician Maurice Kendal once remarked, `would be simpler if the Bayesians followed the example or their master and published posthumously'. Cheers, Berwin __ R-help@r-project.org 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] capturing errors in Sweave
G'day Sundar, On Tue, 2 Mar 2010 01:03:54 -0800 Sundar Dorai-Raj wrote: > Thanks, Berwin. That works just great! You are welcome. I noticed by now that "cat(tmp)" is sufficient; the "tmp[1]" in "cat(tmp[1])" was a left over from earlier attempts to get the output to look correct. Cheers, Berwin __ R-help@r-project.org 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] capturing errors in Sweave
G'day Sundar, On Mon, 1 Mar 2010 23:46:55 -0800 Sundar Dorai-Raj wrote: > Thanks for the input, but I don't want "try" in the Sweave output. I > want the output to look just like it does in the console, as if an > uncaptured error really did occur. I don't think that you will get around using "try"; and you will have to work moderately hard to make the output appear as it does on the console. Probably somewhere along the lines: Sweave code start ++ <>= MySqrt <- function(x) { if (missing(x)) { stop("'x' is missing with no default") } if (!is.numeric(x)) { stop("'x' should only be numeric") } if (x < 0) { stop("'x' should be non-negative") } return(sqrt(x)) } @ <>= tmp <- try(MySqrt()) @ <>= MySqrt() @ <>= cat(tmp[1]) @ <>= tmp <- try(MySqrt("a")) @ <>= MySqrt("a") @ <>= cat(tmp[1]) @ <>= tmp <- try(MySqrt(-2)) @ <>= MySqrt(-2) @ <>= cat(tmp[1]) @ <<>>= MySqrt(4) @ +++ Sweave code end ++ Now what I would like to know is how to include easily warning messages in my Sweave output without having to try whether Jean Lobry's [1] hack still works. :) HTH. Cheers, Berwin [1] https://www.stat.math.ethz.ch/pipermail/r-help/2006-December/121975.html == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] advice/opinion on "<-" vs "=" in teaching R
G'day all, On Fri, 15 Jan 2010 09:08:44 - (GMT) (Ted Harding) wrote: > On 15-Jan-10 08:14:04, Barry Rowlingson wrote: > >> [...] > > I would regard modifying a variable within the parameters of a > > function call as pretty tasteless. What does: > > > > foo(x<-2,x) > > or > > foo(x,x<-3) > > > > do that couldn't be done clearer with two lines of code? Depending on the definition of foo() more than two lines might be necessary to clarify what this code is actually supposed to do. :) > > Remember: 'eschew obfuscation'. Absolutely. > Tasteless or not, the language allows it to be done; and therefore > discussion of distinctions between ways of doing it is relevant to > Erin's question! And a full discussion of these examples would include to warn about the side effects of lazy evaluation: R> rm(list=ls(all=TRUE)) R> foo <- function(x, y){ + if(x > 0 ) x else x + y} R> foo(2, x <- 10) [1] 2 R> x Error: object 'x' not found R> foo(-2, x <- 10) [1] 8 R> x [1] 10 Having said this, I often use "plot(fm <- lm(y ~ . , somedata))". And, yes, students have complained that the code I gave them was not working, that they got a strange error message and they promised that they had typed exactly what I had written on the lab sheet. On checking, they had, of course, typed "plot(fm = lm(y ~ . , somedata))" Cheers, Berwin __ R-help@r-project.org 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] Relase positive with log and zero of negative with 0
G'day Kevin, On Sun, 15 Nov 2009 7:18:18 -0800 wrote: > This is a very simple question but I couldn't form a site search > quesry that would return a reasonable result set. > > Say I have a vector: > > x <- c(0,2,3,4,5,-1,-2) > > I want to replace all of the values in 'x' with the log of x. > Naturally this runs into problems since some of the values are > negative or zero. So how can I replace all of the positive elements > of x with the log(x) and the rest with zero? If you do not mind a warning message: R> x <- c(0,2,3,4,5,-1,-2) R> x <- ifelse(x <= 0,0, log(x)) Warning message: In log(x) : NaNs produced R> x [1] 0.000 0.6931472 1.0986123 1.3862944 1.6094379 0.000 0.000 If you do mind, then: R> x <- c(0,2,3,4,5,-1,-2) R> ind <- x>0 R> x[!ind] <- 0 R> x[ind] <- log(x[ind]) R> x [1] 0.000 0.6931472 1.0986123 1.3862944 1.6094379 0.000 0.000 HTH. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] R and HDF5 Question
G'day Scott, On Fri, 13 Nov 2009 09:52:43 -0700 Scott MacDonald wrote: > I am trying to load an hdf5 file into R and running into some > problems. It's a while that I used hdf5 files and that package in R, but: > This builds fine. The library seems to load without issue, but no > data is returned when I try to load a file: > > > library(hdf5) > > hdf5load("test.h5") > > NULL Is NULL the return of the hdf5load command or are you typing it on the command line? Anyway, .hdf5 files can contain several objects, just as R's .rda file. load() will load an .rda file and put all objects in that file into the workspace. Likewise, hdf5load() loads an hdf5 file and puts all objects in that file into the workspace. > Yet, > > osx:data scott$ h5dump test.h5 HDF5 "test.h5" { GROUP > "/" { DATASET "dset" { DATATYPE H5T_STD_I32LE DATASPACE SIMPLE > { ( 31 ) / ( 31 ) } DATA { (0): 1, 2, 4, 8, 16, 32, 64, 128, 256, > 512, 1024, 2048, 4096, 8192, (14): 16384, 32768, 65536, 131072, > 262144, 524288, 1048576, 2097152, (22): 4194304, 8388608, 16777216, > 33554432, 67108864, 134217728, (28): 268435456, 536870912, > 1073741824 } } } } > > Any thoughts? Did you try an ls() after the hdf5load() command? If the hdf5load() command was successfull, an ls() should show you that an object with name "dset" is now in your workspace; if I read the output above correctly. HTH. Cheers, Berwin __ R-help@r-project.org 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] Unexpected behaviour for as.date()
G'day Isabella, On Tue, 10 Nov 2009 20:11:31 -0800 "Isabella Ghement" wrote: > I tried the solution you suggested and find that am having problems > getting R to extract the year from an object created by as.date(): Perhaps it would be best to first clarify what you really need. :) If you have to extract the year, why go via as.date? Why not just: R> strsplit("02-MAY-00", "-")[[1]][3] [1] "00" R> as.numeric(strsplit("02-MAY-00", "-")[[1]][3]) [1] 0 Or if you have a vector of character strings, each a date in that format: R> x <- c("02-MAY-00", "02-MAY-01") R> sapply(x, function(x) as.numeric(strsplit(x, "-")[[1]][3])) 02-MAY-00 02-MAY-01 0 1 HTH. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] Unexpected behaviour for as.date()
G'day Isabella, On Tue, 10 Nov 2009 18:40:11 -0800 "Isabella Ghement" wrote: > I am trying to use the function as.date() from the "dates" package As far as I can tell, there is no package called "dates", did you mean the package "date"? > in R 2.10.0 to convert a character date to a Julian date, as follows: > > > as.date("02-MAY-01", order="mdy") # convert May 2, 2001 to a Julian > > date > [1] 2May1 Are you sure that you are doing what your comments imply? Try: R> library(date) R> as.date("02-MAY-01", order="mdy") [1] 2May1 R> as.date("02-MAY-2001", order="mdy") [1] 2May2001 R> as.numeric(as.date("02-MAY-01", order="mdy")) [1] -21428 R> as.numeric(as.date("02-MAY-2001", order="mdy")) [1] 15097 > However, when trying to convert a character date from the year 2000 > to a Julian date, I get an instead of the desired Julian date: > > > as.date("02-MAY-00", order="mdy") # convert May 2, 2000 to a Julian > > date > [1] > > Not quite sure why R is unable to handle this type of date (perhaps it > thinks it comes from the year 1900?!). My guess it thinks it comes from the year 0. Not sure why it cannot handle such dates. But then, as far as I know, there is actually some discussion about whether the year 0 exist or whether we went straight from 1BC to 1AD.. > Is there a way to force R to convert character dates from the year > 2000 into Julian dates? Presumably you will need something like: R> as.date(sub("-00", "-2000", "02-MAY-00")) [1] 2May2000 HTH. Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] Constrained Optimization
On Wed, 4 Nov 2009 14:48:08 -0800 (PST) s t wrote: > I'm trying to do the following constrained optimization example. > Maximize x1*(1-x1) + x2*(1-x2) + x3*(1-x3) > s.t. x1 + x2 + x3 = 1 > x1 >= 0 and x1 <= 1 > x2 >= 0 and x2 <= 1 > x3 >= 0 and x3 <= 1 > which are the constraints. > I'm expecting the answer x1=x2=x3 = 1/3. This is a quadratic programming problem (mininimising/maximising a quadratic function under linear equality and inequality constraints). There are several R packages that solve such problems, e.g.: R> library(quadprog) R> Dmat <- diag(3) R> dvec <- rep(1,3)/3 R> Amat <- cbind(rep(1,3), diag(3), -diag(3)) R> bvec <- c(1, rep(0,3), rep(-1,3)) R> solve.QP(Dmat, dvec, Amat, bvec, meq=1) $solution [1] 0.333 0.333 0.333 $value [1] -0.167 $unconstrainted.solution [1] 0.333 0.333 0.333 $iterations [1] 2 0 $iact [1] 1 You may want to consult: http://cran.r-project.org/web/views/Optimization.html HTH. Cheers, Berwin ====== Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] acf gives correlations > 1
G'day Steve, On Mon, 14 Sep 2009 13:44:56 +0200 Steve Jones wrote: > Apologies for the missing data. It can be downloaded from here (22Kb): > http://www.squaregoldfish.co.uk/sekrett/series.csv Well, the Details section of acf's help page states: By default, no missing values are allowed. If the 'na.action' function passes through missing values (as 'na.pass' does), the covariances are computed from the complete cases. This means that the estimate computed may well not be a valid autocorrelation sequence, and may contain missing values. [...] And you have seem to have a massive amount of missing data: R> dat <- scan(url("http://www.squaregoldfish.co.uk/sekrett/series.csv";)) Read 6940 items R> mean(!is.na(dat)) [1] 0.02881844 And, not surprisingly, an even smaller proportion of consecutive, non-missing observations. R> mean(!is.na(dat[-1]) & !is.na(dat[-length(dat)])) [1] 0.006340971 You can find out which formulae are used exactly by acf by studying the code, but this might give you an idea about what is going on: R> ind <- !is.na(dat) R> (mu <- mean(dat[ind])) ## too lazy for mean(dat, na.rm=TRUE) [1] 373.5165 R> (sig2 <- var(dat[ind])) ## ditto [1] 463.4041 R> ind <- which(!is.na(dat[-1]) & !is.na(dat[-length(dat)])) R> sum( (dat[ind]-mu) * (dat[ind+1] - mu)) / length(ind) [1] 593.3041 R> sum( (dat[ind]-mu) * (dat[ind+1] - mu)) / length(ind) / sig2 [1] 1.280317 HTH Cheers, Berwin == Full address Berwin A Turlach Tel.: +61 (8) 6488 3338 (secr) School of Maths and Stats (M019)+61 (8) 6488 3383 (self) The University of Western Australia FAX : +61 (8) 6488 1028 35 Stirling Highway Crawley WA 6009e-mail: ber...@maths.uwa.edu.au Australiahttp://www.maths.uwa.edu.au/~berwin __ R-help@r-project.org 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] Constrained fits: y~a+b*x-c*x^2, with a,b,c >=0
G'day Alex, On Wed, 27 May 2009 11:51:39 +0200 Alex van der Spek wrote: > I wonder whether R has methods for constrained fitting of linear > models. > > I am trying fm<-lm(y~x+I(x^2), data=dat) which most of the time gives > indeed the coefficients of an inverted parabola. I know in advance > that it has to be an inverted parabola with the maximum constrained to > positive (or zero) values of x. > > The help pages for lm do not contain any info on constrained fitting. > > Does anyone know how to? Look at the package nnls on CRAN. According to your subject line, you are trying to solve what is known as a quadratic program, and there are at least two quadratic programming solvers (ipop in kernlab and solve.qp in quadprog) available for R. HTH. Cheers, Berwin === Full address = Berwin A TurlachTel.: +65 6516 4416 (secr) Dept of Statistics and Applied Probability+65 6516 6650 (self) Faculty of Science FAX : +65 6872 3919 National University of Singapore 6 Science Drive 2, Blk S16, Level 7 e-mail: sta...@nus.edu.sg Singapore 117546http://www.stat.nus.edu.sg/~statba __ R-help@r-project.org 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] Getting lm() to work with a matrix
G'day Luc, On Wed, 20 May 2009 09:58:41 -0400 Luc Villandre wrote: > MikSmith wrote: > > [...] > Indeed, functions like /lm()/ require the object fed to the /data/ > argument to be either [...] But the data argument is optional and does not need to be specified. > In your situation, I suggest you typecast your matrix into a data > frame using /as.data.frame()/. [...] My guess is that he is already working with a data frame and does not work with matrices, otherwise he should not have encountered problems: R> response <- matrix(rnorm(120), ncol=4) R> spectra.spec <- matrix(rnorm(900), ncol=30) R> spectra.lm <- lm(response[,3]~spectra.spec[,2:20]) R> spectra.lm Call: lm(formula = response[, 3] ~ spectra.spec[, 2:20]) Coefficients: (Intercept) spectra.spec[, 2:20]1 -0.48404 0.42503 spectra.spec[, 2:20]2 spectra.spec[, 2:20]3 -0.08955-0.27605 spectra.spec[, 2:20]4 spectra.spec[, 2:20]5 -0.16832-0.14107 spectra.spec[, 2:20]6 spectra.spec[, 2:20]7 -0.47009-0.23672 spectra.spec[, 2:20]8 spectra.spec[, 2:20]9 0.12920 0.23306 spectra.spec[, 2:20]10 spectra.spec[, 2:20]11 -0.28586 0.03579 spectra.spec[, 2:20]12 spectra.spec[, 2:20]13 0.10676-0.34407 spectra.spec[, 2:20]14 spectra.spec[, 2:20]15 0.20253-0.17259 spectra.spec[, 2:20]16 spectra.spec[, 2:20]17 0.19765 0.40705 spectra.spec[, 2:20]18 spectra.spec[, 2:20]19 -0.12448-0.17149 Cheers, Berwin __ R-help@r-project.org 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] Unintended loading of package:datasets
G'day David, On Sun, 10 May 2009 21:35:30 -0400 "David A Vavra" wrote: > Thanks. Perhaps something else is going on. There is a large time > period (about 20 sec.) after the message about loading the package. > More investigation, I suppose. What R version are you using? I do not remember ever getting a message that the package datasets is loaded, nor a message for any of the other packages that are loaded by default. (There once, long ago, may have been such a message, but if so, I do not remember this anymore.) Cheers, Berwin __ R-help@r-project.org 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] Unintended loading of package:datasets
G'day David, On Sun, 10 May 2009 17:17:38 -0400 "David A Vavra" wrote: > The dataset package is being loaded apparently by one of the packages > that I am using. The loading of the datasets takes a long time and I > would like to eliminate it. I thought the datasets were effectively > examples so don't understand why they would be required at all. > > 1) How can I determine what is causing the datasets to be loaded? help(options) and then read the part on "defaultPackages" > 2) How can I stop them from doing so? Create an appropriate .Rprofile. R-DownUnder had long time ago a discussion on what people put into their .Rprofile, and Bill Venables' .Rprofile seem to contain the following snippet: ### puts four more packages on to the default ### search path. I use them all the time local({ old <- getOption("defaultPackages") options(defaultPackages = c(old, "splines", "lattice", "mgcv", "MASS")) }) > I am using the following: > > Rpart, grDevices, graphics, stats, utils, methods, base > There is also an environment named 'Autoloads' Rpart?? Or rpart?? I know of a package with the latter name but not the former, and R is case sensitive. So you may try: local({ options(defaultPackages=c("rpart", "grDevices", "graphics", "stats", "utils", "methods") }) FWIW, note that for command XXX, help(XXX) will give you a help page that has usually some example code at the end, that code can be run by example(XXX). Typically, such example code is using data sets from the datasets package; so if you do not load it, the example() command might not work anymore for some functions. Don't complain if your R installation doesn't work anymore if you mess around with defaultPackages, in particular if you remove packages that are usually in the default of defaultPackages. :) And I agree with Rolf (Turner), it is hard to believe that the datasets package would produce a noticeable delay on start-up; for me it also loads instantaneously. I guess that your problem is more along David's (Winsemuis) guess. Cheers, Berwin __ R-help@r-project.org 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] Beyond double-precision?
G'day all, On Sat, 09 May 2009 08:01:40 -0700 spencerg wrote: > The harmonic mean is exp(mean(logs)). Therefore, log(harmonic > mean) = mean(logs). > > Does this make sense? I think you are talking here about the geometric mean and not the harmonic mean. :) The harmonic mean is a bit more complicated. If x_i are positive values, then the harmonic mean is H= n / (1/x_1 + 1/x_2 + ... + 1/x_n) so log(H) = log(n) - log( 1/x_1 + 1/x_2 + ... + 1/x_n) now log(1/x_i) = -log(x_i) so if log(x_i) is available, the logarithm of the individual terms are easily calculated. But we need to calculate the logarithm of a sum from the logarithms of the individual terms. At the C level R's API has the function logspace_add for such tasks, so it would be easy to do this at the C level. But one could also implement the equivalent of the C routine using R commands. The way to calculate log(x+y) from lx=log(x) and ly=log(y) according to logspace_add is: max(lx,ly) + log1p(exp(-abs(lx-ly))) So the following function may be helpful: logadd <- function(x){ logspace_add <- function(lx, ly) max(lx, ly) + log1p(exp(-abs(lx-ly))) len_x <- length(x) if(len_x > 1){ res <- logspace_add(x[1], x[2]) if( len_x > 2 ){ for(i in 3:len_x) res <- logspace_add(res, x[i]) } }else{ res <- x } res } R> set.seed(1) R> x <- runif(50) R> lx <- log(x) R> log(1/mean(1/x)) ## logarithm of harmonic mean [1] -1.600885 R> log(length(x)) - logadd(-lx) [1] -1.600885 Cheers, Berwin =========== Full address = Berwin A TurlachTel.: +65 6515 4416 (secr) Dept of Statistics and Applied Probability+65 6515 6650 (self) Faculty of Science FAX : +65 6872 3919 National University of Singapore 6 Science Drive 2, Blk S16, Level 7 e-mail: sta...@nus.edu.sg Singapore 117546http://www.stat.nus.edu.sg/~statba __ R-help@r-project.org 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] Splitting a vector into equal groups
G'day Utkarsh, On Mon, 04 May 2009 11:51:21 +0530 utkarshsinghal wrote: > I have vector of length 52, say, x=sample(30,52,replace=T). I want to > sort x and split into five *nearly equal groups*. What do you mean by *nearly equal groups*? The size of the groups should be nearly equal? The sum of the elements of the groups should be nearly equal? > Note that the observations are repeated in x so in case of a tie I > want both the observations to fall in same group. Then it becomes even more important to define what you mean with "nearly equal groups". As a start, you may consider: R> set.seed(1) R> x=sample(30,52,replace=T) R> xrle <- rle(sort(x)) R> xrle Run Length Encoding lengths: int [1:25] 2 1 2 2 3 1 1 1 5 1 ... values : int [1:25] 1 2 4 6 7 8 9 11 12 13 ... R> cumsum(xrle$lengths) [1] 2 3 5 7 10 11 12 13 18 19 24 25 26 28 29 32 35 38 [19] 43 45 46 48 49 51 52 and use this to determine our cut-offs. E.g., should the first group have 10, 11 or 12 elements in this case? The information in xrle should enable you to construct your five groups once you have decided on a grouping. HTH. Cheers, Berwin === Full address = Berwin A TurlachTel.: +65 6516 4416 (secr) Dept of Statistics and Applied Probability+65 6516 6650 (self) Faculty of Science FAX : +65 6872 3919 National University of Singapore 6 Science Drive 2, Blk S16, Level 7 e-mail: sta...@nus.edu.sg Singapore 117546http://www.stat.nus.edu.sg/~statba __ R-help@r-project.org 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] The .tex version of the manual in foo.Rcheck
G'day Uwe, On Wed, 29 Apr 2009 11:03:43 +0200 Uwe Ligges wrote: > you can also take a source package an run > > R CMD Rd2dvi --no-clean packageName > > on it and you will get a temporary directory with the TeX sources in > it. Which is fine for manual processing and when done once or twice, but somewhat less helpful for automatic processing in scripts since the name of the temporary directory is hard to predict. `R CMD check foo' copies already unconditionally the foo-manual.log file from the temporary directory to the foo.Rcheck directory; and the foo-manual.tex file is copied if an error occurs during processing. What is the problem with also copying the foo-manual.tex file unconditionally to foo.Rcheck? Cheers, Berwin __ R-help@r-project.org 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] The .tex version of the manual in foo.Rcheck
G'day Bendix, On Mon, 27 Apr 2009 22:52:02 +0200 "BXC (Bendix Carstensen)" wrote: > In version 2.8.1, running Rcmd check on the package foo would leave > the file foo-manual.tex in the folder foo.Rcheck. > > But as of 2.9.0 only foo-manual.pdf and foo-manual.log are there. > > Is this intentional? I do not know whether it is intentional, but it seems that parts of `check' were substantially re-written due to the new Rd parser. As far as I can understand the Perl code, the .tex file is only copied to the .Rcheck directory if an error occurs while the file is processed. But this can easily be changed by moving two lines in the script a bit higher, see the attached patch. Cheers, Berwin === Full address ===== Berwin A TurlachTel.: +65 6516 4416 (secr) Dept of Statistics and Applied Probability+65 6516 6650 (self) Faculty of Science FAX : +65 6872 3919 National University of Singapore 6 Science Drive 2, Blk S16, Level 7 e-mail: sta...@nus.edu.sg Singapore 117546http://www.stat.nus.edu.sg/~statba __ R-help@r-project.org 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] Unable to re-import a table that was just exported to a file
On Sun, Apr 26, 2009 at 7:38 PM, Nigel Birney wrote: > > Hi all, > > I am saving a program's output to a file to be read by another > algorithm. But somehow such a simple operation (the reading) fails. > I get: > > Error in read.table("a_corr_data.txt", sep = ",", col.names = T, > row.names = F) : > __more columns than column names > > > Here is the write statement: > > write.table(a_corr_data,"a_corr_data.txt",sep=",",col.names=T,row.names=F) > > Here is the read statement: > > a_corr_data <- > read.table("a_corr_data.txt",sep=",",col.names=T,row.names=F) > > Nothing happens in-between (these actions are just 10-30 secs > apart). I tried to export/import without col.names, also tried > different deliminators ("/t") but the same error pops up again and > again. I am already quite unhappy with this. ?read.table You might find that you are using col.names and row.names wrongly in the read.table command. Best wishes, Berwin === Full address = Berwin A TurlachTel.: +65 6516 4416 (secr) Dept of Statistics and Applied Probability+65 6516 6650 (self) Faculty of Science FAX : +65 6872 3919 National University of Singapore 6 Science Drive 2, Blk S16, Level 7 e-mail: sta...@nus.edu.sg Singapore 117546http://www.stat.nus.edu.sg/~statba __ R-help@r-project.org 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] Code of sample() in C
G'day Ranjan, On Sun, 5 Apr 2009 17:00:56 -0500 Ranjan Maitra wrote: > [...[ > I haven't tried this function out much myself, so YMMV. There is an obvious problem in this code since len is not decreased but n is when an index is sampled. The last line in the last for loop should be x[j] = x[--len]; Otherwise this algorithm will produce samples in which an index can be repeated. And the probabilities with which an index is sampled would be non-trivial to determine; they would definitely not be uniform. HTH. Cheers, Berwin > #include > #ifndef USING_RLIB > #define MATHLIB_STANDALONE 1 /*It is essential to have this before > the call to the Rmath's header file because this decides >the definitions to be set. */ > #endif > #include > > /* Note that use of this function involves a prior call to the Rmath > library to get the seeds in place. It is assumed that this is done in > the calling function. */ > > /* Equal probability sampling; without-replacement case */ > /* Adapted from the R function called SampleNoReplace */ > > /** > * Stores k out of n indices sampled at random without replacement > * in y. > */ > int srswor(int n, int k, int *y) > { > if (k > n) { > return 1; > } > else { > const double len = (double) n; > int i; > int* x = malloc(n * sizeof(int)); > if (!x) { > return 2; > } > > for (i = 0; i < n; ++i) x[i] = i; > > for (i = 0; i < k; ++i) { > const int j = (int)(len * runif(0.0, 1.0)); > y[i] = x[j]; > x[j] = x[--n]; > } > free(x); > } > return 0; > } === Full address = Berwin A TurlachTel.: +65 6516 4416 (secr) Dept of Statistics and Applied Probability+65 6516 6650 (self) Faculty of Science FAX : +65 6872 3919 National University of Singapore 6 Science Drive 2, Blk S16, Level 7 e-mail: sta...@nus.edu.sg Singapore 117546http://www.stat.nus.edu.sg/~statba __ R-help@r-project.org 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] Eclipse and StatET Howto (also added Subversion, Rtools)
G'day Peter, On Sun, 05 Apr 2009 11:26:40 +0200 Peter Dalgaard wrote: > Berwin A Turlach wrote: > > G'day Dirk, > > > > On Sat, 4 Apr 2009 20:27:22 -0500 > > Dirk Eddelbuettel wrote: > > > >> On 4 April 2009 at 14:37, Ken-JP wrote: > >>> Yes, I have x11-common installed, and dpkg -S /etc/X11/rgb.txt > >>> shows "not found" for me. This is on Ubuntu 8.10 amd64. > >> Same 8.10 for both amd64 and i386 where I checked -- both have the > >> file. It could be a leftover from an earlier install of mine, or an > >> accidental deletion at your end. > > > > My guess that it is the former. :) Some time ago I wiped Debian of > > my laptop and installed Kubuntu 8.10 fresh (i386 flavour); so no > > left overs from previous versions and on my laptop the results are: > > > > ber...@berwin5:~$ apt-file find /etc/X11/rgb.txt > > ber...@berwin5:~$ dpkg -S /etc/X11/rgb.txt > > dpkg: /etc/X11/rgb.txt not found. > > Hum.. Fedora 9 doesn't have it either. > > It does have /usr/share/X11/rgb.txt though, so please check whether > it has moved. Apparently it has: ber...@berwin5:~$ apt-file find rgb.txt [...] x11-common: /usr/share/X11/rgb.txt [...] However, ber...@berwin5:~$ ls -l /usr/share/X11/rgb.txt lrwxrwxrwx 1 root root 16 Jan 14 03:28 /usr/share/X11/rgb.txt -> /etc/X11/rgb.txt ber...@berwin5:~$ more /usr/share/X11/rgb.txt /usr/share/X11/rgb.txt: No such file or directory Cheers, Berwin __ R-help@r-project.org 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] Eclipse and StatET Howto (also added Subversion, Rtools)
G'day Dirk, On Sat, 4 Apr 2009 20:27:22 -0500 Dirk Eddelbuettel wrote: > On 4 April 2009 at 14:37, Ken-JP wrote: > > Yes, I have x11-common installed, and dpkg -S /etc/X11/rgb.txt > > shows "not found" for me. This is on Ubuntu 8.10 amd64. > > Same 8.10 for both amd64 and i386 where I checked -- both have the > file. It could be a leftover from an earlier install of mine, or an > accidental deletion at your end. My guess that it is the former. :) Some time ago I wiped Debian of my laptop and installed Kubuntu 8.10 fresh (i386 flavour); so no left overs from previous versions and on my laptop the results are: ber...@berwin5:~$ apt-file find /etc/X11/rgb.txt ber...@berwin5:~$ dpkg -S /etc/X11/rgb.txt dpkg: /etc/X11/rgb.txt not found. Cheers, Berwin __ R-help@r-project.org 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.