Re: [R] 'singular gradient matrix? when using nls() and how to make the program skip nls( ) and run on
On Wed, Sep 05, 2007 at 04:43:19PM -0700, Yuchen Luo wrote: > Dear friends. > > I use nls() and encounter the following puzzling problem: > > > > I have a function f(a,b,c,x), I have a data vector of x and a vectory y of > realized value of f. > > > > Case1 > > I tried to estimate c with (a=0.3, b=0.5) fixed: > > nls(y~f(a,b,c,x), control=list(maxiter = 10, minFactor=0.5 > ^2048),start=list(c=0.5)). > > The error message is: "number of iterations exceeded maximum of 10" > > > > Case2 > > I then think maybe the value of a and be are not reasonable. So, I let nls() > estimate (a,b,c) altogether: > > nls(y~f(a,b,c,x), control=list(maxiter = 10, minFactor=0.5 > ^2048),start=list(a=0.3,b=0.5,c=0.5)). > > The error message is: > > "singular gradient matrix at initial parameter estimates". > > > > This is what puzzles me, if the initial parameter of (a=0.3,b=0.5,c=0.5) can > create 'singular gradient matrix', then why doesn't this 'singular gradient > matrix' appear in Case1? > > > > I have tried to change the initial value of (a,b,c) around but the problem > persists. I am wondering if there is a way out. > > > > My another question is, I need to run 220 of nls() in my program with > different y and x. When one of the nls() encounter a problem, the whole > program stops. In my case, the 3rd nls() runs into a problem. I would > still need the program to run the remaining 217 nls( )! Is there a way to > make the program skip the problematic nls() and complete the ramaining > nls()'s? ?try > > > > Your help will be highly appreciated! > > Yuchen Luo > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] non-linear fitting (nls) and confidence limits
dear list members, I apologize in advance for posting a second time, but probably after one week chances are, the first try went down the sink.. my question concerns computation of confidence intervals in nonlinear fits with `nls' when weigthing the fit. the seemingly correct procedure does not work as expected, as is detailed in my original post below. any remarks appreciated. greetings joerg original post: -- for unweighted fits using `nls' I compute confidence intervals for the fitted model function by using: #--- se.fit <- sqrt(apply(rr$m$gradient(), 1, function(x) sum(vcov(rr)*outer(x,x luconf <- yfit + outer(se.fit, qnorm(c(probex, 1 - probex))) #--- where `rr' contains an `nls' object, `x' is the independent variable vector, `yfit' the corresponding model prediction (`fitted(rr)'), `se.fit' the corresponding standard error and `luconf' the lower and upper confidence limits at some level specified by `probex'. when using the same approach with weighted fits (even if all weights are equal (but not equal to one)), I noted that the confidence limits depend on the weights in an counter-intuitive, probably erroneous way: consider the following simple example (yes, I know that this actually is a linar model :-)): #--- probex <- 0.05 x <- 1:10 y <- rnorm(x, x, .8) w1 <- rep(1, 10) w2 <- w1; w2[7:10] <- 0.01 * w2[7:10] rr1 <- nls(y ~ a*x + b, data = list(x = x, y = y), start = list(a = 1, b = 0), weights = w1) rr2 <- nls(y ~ a*x + b, data = list(x = x, y = y), start = list(a = 1, b = 0), weights = w2) yfit1 <- fitted(rr1) yfit2 <- fitted(rr2) se.fit1 <- sqrt(apply(rr1$m$gradient(), 1, function(x) sum(vcov(rr1)*outer(x,x luconf1 <- yfit1 + outer(se.fit1, qnorm(c(probex, 1 - probex))) se.fit2 <- sqrt(apply(rr2$m$gradient(), 1, function(x) sum(vcov(rr2)*outer(x,x luconf2 <- yfit2 + outer(se.fit2, qnorm(c(probex, 1 - probex))) op <- par(mfrow = c(2,1)) matplot(x, cbind(y, yfit1,luconf1), type = 'blll', pch = 1, col = c(1,2,4,4), lty = 1) matplot(x, cbind(y, yfit2,luconf2), type = 'blll', pch = 1, col = c(1,2,4,4), lty = 1) par(op, no.readonly = T) #--- the second fit uses unequal weights where the last 4 points are next to unimportant for the fit. the confidence intervals in this case are fine up to point no. 6, but nonsense afterwards I believe. I have tracked this down to the fact that `m$gradient' does _not_ contain the actual gradients w.r.t. the parameters but rather the gradients multiplied by the sqrt of the weights. without trying to understand the inner workings of the `nls*' functions in detail my questions are: 1. is the fact that `m$gradient' in an `nls' object does contain the scaled gradients documented somewhere (I did'nt find it)? if yes, where is it, if no, can someone please give me a hint what the rationale behind this approach is? 2. am I right to understand, that the above approach to compute `se.fit' (essentially in this compact form proposed by p. daalgaard on r-help some time ago) is erroneous (for weighted nls) in computing the confidence limits and that I have to undo the multiplication by sqrt(weights) which is included in `m$gradient' (or compute the gradients of my model function myself, e.g. with `deriv')? thanks, joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] weighted nls and confidence intervals
for unweighted fits using `nls' I compute confidence intervals for the fitted model function by using: #--- se.fit <- sqrt(apply(rr$m$gradient(), 1, function(x) sum(vcov(rr)*outer(x,x luconf <- yfit + outer(se.fit, qnorm(c(probex, 1 - probex))) #--- where `rr' contains an `nls' object, `x' is the independent variable vector, `yfit' the corresponding model prediction (`fitted(rr)'), `se.fit' the corresponding standard error and `luconf' the lower and upper confidence limits at some level specified by `probex'. when using the same approach with weighted fits (even if all weights are equal (but not equal to one)), I noted that the confidence limits depend on the weights in an counter-intuitive, probably erroneous way: consider the following simple example (yes, I know that this actually is a linar model :-)): #--- probex <- 0.05 x <- 1:10 y <- rnorm(x, x, .8) w1 <- rep(1, 10) w2 <- w1; w2[7:10] <- 0.01 * w2[7:10] rr1 <- nls(y ~ a*x + b, data = list(x = x, y = y), start = list(a = 1, b = 0), weights = w1) rr2 <- nls(y ~ a*x + b, data = list(x = x, y = y), start = list(a = 1, b = 0), weights = w2) yfit1 <- fitted(rr1) yfit2 <- fitted(rr2) se.fit1 <- sqrt(apply(rr1$m$gradient(), 1, function(x) sum(vcov(rr1)*outer(x,x luconf1 <- yfit1 + outer(se.fit1, qnorm(c(probex, 1 - probex))) se.fit2 <- sqrt(apply(rr2$m$gradient(), 1, function(x) sum(vcov(rr2)*outer(x,x luconf2 <- yfit2 + outer(se.fit2, qnorm(c(probex, 1 - probex))) op <- par(mfrow = c(2,1)) matplot(x, cbind(y, yfit1,luconf1), type = 'blll', pch = 1, col = c(1,2,4,4), lty = 1) matplot(x, cbind(y, yfit2,luconf2), type = 'blll', pch = 1, col = c(1,2,4,4), lty = 1) par(op, no.readonly = T) #--- the second fit uses unequal weights where the last 4 points are next to unimportant for the fit. the confidence intervals in this case are fine up to point no. 6, but nonsense afterwards I believe. I have tracked this down to the fact that `m$gradient' does _not_ contain the actual gradients w.r.t. the parameters but rather the gradients multiplied by the sqrt of the weights. without trying to understand the inner workings of the `nls*' functions in detail my questions are: 1. is the fact that `m$gradient' in an `nls' object does contain the scaled gradients documented somewhere (I did'nt find it)? if yes, where is it, if no, can someone please give me a hint what the rationale behind this approach is? 2. am I right to understand, that the above approach to compute `se.fit' (essentially in this compact form proposed by p. daalgaard on r-help some time ago) is erroneous (for weighted nls) in computing the confidence limits and that I have to undo the multiplication by sqrt(weights) which is included in `m$gradient' (or compute the gradients of my model function myself, e.g. with `deriv')? thanks, joerg __ R-help@stat.math.ethz.ch 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 store the parameter estimated by nls( ) to a variable?
On Sun, Aug 12, 2007 at 10:50:59AM -0700, Yuchen Luo wrote: > Dear Professor Murdoch. > Thank you so much > > Best Wishes > Yuchen Luo > > On 8/12/07, Duncan Murdoch <[EMAIL PROTECTED]> wrote: > > > > Yuchen Luo wrote: > > > Dear Colleagues. > > > > > > I believe this should be a problem encountered by many: > > > > > > nls( ) is a very useful and efficient function to use if we are just to > > > display the estimated value on screen. What if we need R to store the > > > estimated parameter in a variable? > > > > > > For example: > > > > > > x=rnorm(10, mean=1000, sd=10) > > > > > > y=x^2+100+rnorm(10) > > > > > > a=nls(y~(x^2+para),control=list(maxiter = 1000, minFactor=0.5 > > > ^1024),start=list(para=0.0)) > > > > > > How to store the estimated value of para in this case, in a variable, > > say, > > > b? > > > > > > It is easy to display a and find all the information. How ever, I need > > to > > > fit a different set of x and y in every loop of my code and I need to > > store > > > the estimated values for further use. I have checked both the online > > manual > > > and several S-plus books but no example as such showed up. > > > > > > Your help will be highly appreciated! > > > > coef(a) will get what you want. coef() works for most modelling > > functions where it makes sense. > > > > Duncan Murdoch > > > further information might be extracted from the `summary' return value. e.g., if you need not only the parameters but also their error estimates, you could get them via summary(a)$parameters[, "Std. Error"] (question: is this the canonical way to do it or is their a better one?) __ R-help@stat.math.ethz.ch 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] Reducing the size of pdf graphics files produced with R
On Wed, May 23, 2007 at 10:13:22AM -0700, Waichler, Scott R wrote: > > > Scott > > Scott Waichler, Senior Research Scientist > Pacific Northwest National Laboratory > MSIN K9-36 > P.O. Box 999 > Richland, WA 99352USA > 509-372-4423 (voice) > 509-372-6089 (fax) > [EMAIL PROTECTED] > http://hydrology.pnl.gov > > --- > > > > -Original Message- > > From: Joerg van den Hoff [mailto:[EMAIL PROTECTED] > > Sent: Wednesday, May 23, 2007 9:25 AM > > To: Waichler, Scott R > > Cc: r-help@stat.math.ethz.ch; [EMAIL PROTECTED] > > Subject: Re: [R] Reducing the size of pdf graphics files > > produced with R > > > > On Wed, May 23, 2007 at 07:24:04AM -0700, Waichler, Scott R wrote: > > > > as you are using MacOS X, you'll have ghostscript > > installed anyway. > > > > so try in R `dev2bitmap' with `type =pdfwrite'. I believe `gs' > > > > _does_ include compression. a quick test showed at least > > a reduction > > > > by about a factor of 2 relative to `pdf()'. probably one > > can fiddle > > > > with the ghostscript settings (cf. e.g. `Ps2pdf.htm' in the > > > > ghostscipt > > > > docs: you > > > > can adjust the resolution for images in the pdf file) to improve > > > > this, so as a last resort you could indeed export the graphics as > > > > postscript and do the conversion to `pdf' by adjusting > > the `ps2pdf' > > > > switches. but even with the default settings the pdf produced via > > > > dev2bitmap/ghostscript is the better solution. apart from > > file size > > > > I by and then ran into problems when converting `pdf()' output to > > > > postscript later on, for instance. > > > > > > Can you give an example of dev2bitmap usage? I tried using it in > > `dev2bitmap(file = "rf.pdf", type = "pdfwrite")' copies the > > current device to the pdf-file `rf.pdf', i.e. you should have > > plotted something on the screen prior to using this (the > > manpage tells you so much...). no `dev.off' is necessary in this case. > > > > in order to "plot directly" into the pdffile, you can use > > `bitmap' instead of `dev2bitmap', i.e. > > > > use either: > > > > plot(1:10) > > dev2bitmap(file = "rf1.pdf", type = "pdfwrite") > > > > or: > > > > bitmap(file = "rf2.pdf", type = "pdfwrite") > > plot(1:10) > > dev.off() > > > > both should produce the desired output file (at least after > > including the correct `width' and `height' settings). > > I tried the second method, using width=8.5, height=11, but the margins > in the output were bad, as if it assumed a4 paper or something. I tried > inserting paper="special" and paper="letter" in the bitmap call, but > both caused errors. Also, the plot was in grayscale instead of color. > > Scott Waichler if you look at the source code of `bitmap' you'll see that it essentially constructs a `gs' command (the string `cmd') which is fed by the output of the R `postscript' device via a pipe (the manpage of `bitmap' tells the same). so I presume if you don't get colour or seem to have wrong page size this has nothing to do with R, but you need to check your ghostscript installation. only guessing, though. joerg __ R-help@stat.math.ethz.ch 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] Reducing the size of pdf graphics files produced with R
On Wed, May 23, 2007 at 07:24:04AM -0700, Waichler, Scott R wrote: > > as you are using MacOS X, you'll have ghostscript installed anyway. so > > try in R `dev2bitmap' with `type =pdfwrite'. I believe `gs' _does_ > > include compression. a quick test showed at least a reduction by about > > a factor of 2 relative to `pdf()'. probably one can fiddle with the > > ghostscript settings (cf. e.g. `Ps2pdf.htm' in the ghostscipt > > docs: you > > can adjust the resolution for images in the pdf file) to > > improve this, so as a last resort you could indeed export the graphics > > as postscript and do the conversion to `pdf' by adjusting the `ps2pdf' > > switches. but even with the default settings the pdf produced via > > dev2bitmap/ghostscript is the better solution. apart from file size I > > by and then ran into problems when converting `pdf()' output to > > postscript later on, for instance. > > Can you give an example of dev2bitmap usage? I tried using it in place > of a pdf() statement. An X11 window opened and my figure flew by, but I > didn't get the file output. I also used dev2bitmap after opening a > pdf() and just before the dev.off() statement, since the help says it > works on the "current device", but again no written output. What am I > doing wrong? > > I tried: > dev2bitmap(file = plotfile2, type="pdfwrite", width=8.5, height=11, > pointsize=12) > print(myplot()) > dev.off() > > and > > pdf(file = plotfile, paper="letter", width=8.5, height=11, > pointsize=12) > print(myplot()) > dev2bitmap(file = plotfile2, type="pdfwrite", width=8.5, height=11, > pointsize=12) > dev.off() > > Thanks, > Scott Waichler > scott.waichler _at_ pnl.gov `dev2bitmap(file = "rf.pdf", type = "pdfwrite")' copies the current device to the pdf-file `rf.pdf', i.e. you should have plotted something on the screen prior to using this (the manpage tells you so much...). no `dev.off' is necessary in this case. in order to "plot directly" into the pdffile, you can use `bitmap' instead of `dev2bitmap', i.e. use either: plot(1:10) dev2bitmap(file = "rf1.pdf", type = "pdfwrite") or: bitmap(file = "rf2.pdf", type = "pdfwrite") plot(1:10) dev.off() both should produce the desired output file (at least after including the correct `width' and `height' settings). joerg __ R-help@stat.math.ethz.ch 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] Reducing the size of pdf graphics files produced with R
On Tue, May 22, 2007 at 01:06:06PM -0400, Chabot Denis wrote: > Thank you Prof. Ripley. > > Believe me, I do not have the skills to contribute such a thing as a > stream compressor and I DO appreciate the work and usefulness of the > pdf device as it is. I do most of my plots with pdf device, the rest > with quartz (especially when I'm not sure I'll want to save a plot) > and (rarely) png when the pdf output is too large or for > compatibility with microsoft applications. > > I find the statement you took from the help page promising: I often > include these large plots into LaTeX, so I'll investigate what form > of compression pdftex can do. > > Sincerely, > > Denis > Le 07-05-22 à 12:47, Prof Brian Ripley a écrit : > > >> From the help page > > > > 'pdf' writes uncompressed PDF. It is primarily intended for > > producing PDF graphics for inclusion in other documents, and > > PDF-includers such as 'pdftex' are usually able to handle > > compression. > > > > If you are able to contribute a stream compressor, R will produce > > smaller plots. Otherwise it is unlikely to happen (and it any case > > would be a > > smaller contribution than that of the author of pdf(), who is quite > > happy with external compressors). > > > > Acrobat does other things (not all of which it tells you about), > > but compression is the main advantage. > > > > On Tue, 22 May 2007, Chabot Denis wrote: > > > >> Hi, > >> > >> Without trying to print 100 points (see >> finzi.psych.upenn.edu/R/Rhelp02a/archive/42105.html>), I often print > >> maps for which I do not want to loose too much of coastline detail, > >> and/or plots with 1000-5000 points (yes, some are on top of each > >> other, but using transparency (i.e. rgb colors with alpha > >> information) this actually comes through as useful information. > >> > >> But the files are large (not as large as in the thread above of > >> course, 800 KB to about 2 MB), especially when included in a LaTeX > >> document by the dozen. > >> > >> Acrobat (not the reader, the full program) has an option "reduce file > >> size". I don't know what it does, but it shrinks most of my plots to > >> about 30% or original size, and I cannot detect any loss of detail > >> even when zooming several times. But it is a pain to do this with > >> Acrobat when you generate many plots... And you need to buy Acrobat. > >> > >> Is this something the pdf device could do in a future version? I > >> tried the "million points" example from the thread above and the 55 > >> MB file was reduced to 6.9 MB, an even better shrinking I see on my > >> usual plots. > >> > >> > >> Denis Chabot > >> > >> __ > >> R-help@stat.math.ethz.ch mailing list > >> https://stat.ethz.ch/mailman/listinfo/r-help > >> PLEASE do read the posting guide http://www.R-project.org/posting- > >> guide.html > >> and provide commented, minimal, self-contained, reproducible code. > >> > > > > -- > > Brian D. Ripley, [EMAIL PROTECTED] > > Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ > > University of Oxford, Tel: +44 1865 272861 (self) > > 1 South Parks Road, +44 1865 272866 (PA) > > Oxford OX1 3TG, UKFax: +44 1865 272595 > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. as an attempt to suggest something more helpful than "do write the compressor yourself if you have a problem with pdf()": as you are using MacOS X, you'll have ghostscript installed anyway. so try in R `dev2bitmap' with `type =pdfwrite'. I believe `gs' _does_ include compression. a quick test showed at least a reduction by about a factor of 2 relative to `pdf()'. probably one can fiddle with the ghostscript settings (cf. e.g. `Ps2pdf.htm' in the ghostscipt docs: you can adjust the resolution for images in the pdf file) to improve this, so as a last resort you could indeed export the graphics as postscript and do the conversion to `pdf' by adjusting the `ps2pdf' switches. but even with the default settings the pdf produced via dev2bitmap/ghostscript is the better solution. apart from file size I by and then ran into problems when converting `pdf()' output to postscript later on, for instance. hth, joerg __ R-help@stat.math.ethz.ch 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] save intermediate result
On Fri, May 04, 2007 at 04:55:36PM -0400, Weiwei Shi wrote:
> sorry for my English, staying in US over 7 years makes me think I were
> ok, now :(
>
> anyway, here is an example and why I need that:
>
> suppose I have a function like the following:
>
> f1 <- function(){
> line1 <- f2() # assume this line takes very long time, like more than 30 min
> # then I need to save line1 as an object into workspace b/c
> # the following lines after this might have some bugs
> # currently I use
> save(line1, file="line1.robj")
> # but I am thinking if there is anothe way, like enviroment, instead
> of save it into a file
>
> # codes as followed might have some bugs
> # blabla...
> }
yes, that's what I thought you meant. so use
line1 <<- f2()
i.e. the assignment operator `<<-' instead of `<-'. after f2() is
finished, the object `line1' is available from the R prompt, i.e.
visible in the workspace (for details, cf. `help("<<-")').
>
> if some codes as followed have bugs, my f1 function cannot return
> anything, which means I have to re-run the time-consuming line again.
you probably should debug it using some dummy data instead of the f2()
call...
>
> thanks,
>
> Weiwei
>
> On 5/4/07, Joerg van den Hoff <[EMAIL PROTECTED]> wrote:
> >On Fri, May 04, 2007 at 03:45:10PM -0400, Weiwei Shi wrote:
> >> hi,
> >>
> >> is there a way to save a R object into workspace instead of into a
> >> file during a running of function?
> >>
> >
> >if I understand the question correctly you want the 'super-assignment'
> >operator `<<-' as in
> >
> >-cut---
> >R> f <- function(x) y <<- 2*x
> >R> f(1)
> >R> y
> >[1] 2
> >R> f(10)
> >R> y
> >[1] 20
> >-cut---
> >
> >i.e. y is assigned a value in the 'workspace'. be careful, though, with
> >this kind of assignment. why not return the object you are interested
> >in as a list component together with other results when the function
> >call finishes?
> >
> >
>
>
> --
> Weiwei Shi, Ph.D
> Research Scientist
> GeneGO, Inc.
>
> "Did you always know?"
> "No, I did not. But I believed..."
> ---Matrix III
__
R-help@stat.math.ethz.ch 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] save intermediate result
On Fri, May 04, 2007 at 03:45:10PM -0400, Weiwei Shi wrote: > hi, > > is there a way to save a R object into workspace instead of into a > file during a running of function? > if I understand the question correctly you want the 'super-assignment' operator `<<-' as in -cut--- R> f <- function(x) y <<- 2*x R> f(1) R> y [1] 2 R> f(10) R> y [1] 20 -cut--- i.e. y is assigned a value in the 'workspace'. be careful, though, with this kind of assignment. why not return the object you are interested in as a list component together with other results when the function call finishes? __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] code/documentation mismatch with backslashes in argument list
I have a function definition such as
f <- function (pattern = ".*\\.txt") {}
in the manpage this has to be documented as
f <- function (pattern = ".*.txt")
in order to get the correct display (with double backslash) in the R console
when issuing `?f', but this causes complains from `R CMD CHECK' (code
documentation mismatch), which is not nice.
question: is there (or will there be with future releases) a way to avoid these
warnings without loosing the correct display of the manpage?
regards,
joerg
__
R-help@stat.math.ethz.ch 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] nls.control( ) has no influence on nls( ) !
On Mon, Apr 16, 2007 at 09:03:27AM +0200, Martin Maechler wrote:
> > "Yuchen" == Yuchen Luo <[EMAIL PROTECTED]>
> > on Sun, 15 Apr 2007 12:18:23 -0700 writes:
>
> Yuchen> Dear Friends.
> Yuchen> I tried to use nls.control() to change the 'minFactor' in nls( ),
> but it
> Yuchen> does not seem to work.
>
> yes, you do not seem to use it correctly.
> No reason to jump to the conclusion that you give in the subject
> of this posting ... which is hilariously wrong!
> You don't really think that a quality software like R has had nls()
> and nls.control(), and nls.control() would never have worked for all
> those years and tens of thousands of users ???
> Please get to your senses, and first read the posting guide (*)
> - and then try again, so we can help you further!
>
> Regards,
> Martin Maechler, ETH Zurich
>
> Yuchen> I used nls( ) function and encountered error message "step factor
> Yuchen> 0.000488281 reduced below 'minFactor' of 0.000976563". I then
> tried the
> Yuchen> following:
>
> Yuchen> 1) Put "nls.control(minFactor = 1/(4096*128))" inside the
> brackets of nls,
> Yuchen> but the same error message shows up.
>
> Yuchen> 2) Put "nls.control(minFactor = 1/(4096*128))" as a separate
> command before
> Yuchen> the command that use nls( ) function, again, the same thing
> happens,
> Yuchen> although the R responds to the nls.control( ) function
> immediately:
you need to provide a list for the `control' argument of `nls', i.e. you need
to use something like:
nls(other_arguments, control = nls.control(minFactor = 1/4096))
nls.control() only helps in so far, as that it always returns the complete list
of three components needed for the `control' argument of `nls'. you can equally
well use an explicit list instead.
> Yuchen> -
> Yuchen> $maxiter
>
> Yuchen> [1] 50
>
>
>
> Yuchen> $tol
>
> Yuchen> [1] 1e-05
>
>
>
> Yuchen> $minFactor
>
> Yuchen> [1] 1.907349e-06
>
> Yuchen> --
>
>
> Yuchen> I am wondering how may I change the minFactor to a smaller value?
> The manual
> Yuchen> that comes with the R software about nls( ) is very sketchy ---
> the only
> Yuchen> relevant example I see is a separate command like 2).
>
> Yuchen> A more relevent question might be, is lower the 'minFactor' the
> only
> Yuchen> solution to the problem? What are the other options?
check your model?
check your data?
change start values ("wrong" minimum)?
use another algorithm (cf. nls manpage)?
>
> Yuchen> Best Wishes
> Yuchen> Yuchen Luo
a final remark (off-topic?), concerning the response of m. maechler (but I
notice this attitude continously on the list): this is supposed to be a help
list, neither a boot camp nor a thought-police school, right?.
concerning the question at hand, the one word response
`?nls'
would have been justifiable (it usually _is_ in the manpage), if not overly
helpful, probably. but the actual response is worded such that --
would it be directed at me -- I would judge it simply offensive: "get to your
senses" is not what one should be told, simply because ond did not grasp how to
use the `control' argument and/or because in the subject (contrary to the body)
of
the mail 'has no influence' is used instead of `seems to have no influence'.
I really do think that if someone deems a certain question really stupid, it
would be wise simply to keep 'radio silence' (i.e. ignore it), which probably
is the best corrective anyway. the 'how dare you' responses are annoying, to
say the least.
for heavens sake...
__
R-help@stat.math.ethz.ch 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] replacing all NA's in a dataframe with zeros...
On Thu, Mar 15, 2007 at 10:21:22AM +0100, Peter Dalgaard wrote:
> Gavin Simpson wrote:
> > On Wed, 2007-03-14 at 20:16 -0700, Steven McKinney wrote:
> >
> >> Since you can index a matrix or dataframe with
> >> a matrix of logicals, you can use is.na()
> >> to index all the NA locations and replace them
> >> all with 0 in one command.
> >>
> >>
> >
> > A quicker solution, that, IIRC, was posted to the list by Peter
> > Dalgaard several years ago is:
> >
> > sapply(mydata.df, function(x) {x[is.na(x)] <- 0; x}))
> >
> I hope your memory fails you, because it doesn't actually work.
>
> > sapply(test.df, function(x) {x[is.na(x)] <- 0; x})
> x1 x2 x3
> [1,] 0 1 1
> [2,] 2 2 0
> [3,] 3 3 0
> [4,] 0 4 4
>
> is a matrix, not a data frame.
>
> Instead:
>
> > test.df[] <- lapply(test.df, function(x) {x[is.na(x)] <- 0; x})
> > test.df
> x1 x2 x3
> 1 0 1 1
> 2 2 2 0
> 3 3 3 0
> 4 0 4 4
>
> Speedwise, sapply() is doing lapply() internally, and the assignment
> overhead should be small, so I'd expect similar timings.
just an idea:
given the order of magnitude difference (factor 17 or so) in runtime
between the "obvious" solution and the fast one: would'nt it be
possible/sensible
to modify the corresponding subsetting method ("[.data.frame") such that it
recognizes the case when it is called with an arbitrary index matrix (the
problem is not restricted to indexing with a logical matrix, I presume?) and
switch internally to the fast solution given above?
in my (admittedly limited) experience it seems that one of the not so nice
properties of R is that one encounters in quite a few situations exactly the
above
situation: unexpected massive differences in run time between different
solutions (I'm not
talking about "explicit loop penalty"). what concerns me most, are the very
basic scenarios (not complex algorithms): data frames vs. matrices, naming
vector components or not, subsetting, read.table vs. scan, etc. if their were a
concise HOW TO list for the cases "when speed matters", that would be helpful,
too.
I understand that part of the "uneven performance" is unavoidable and one must
expect the user to go to the trouble to understand the reasons, e.g. for
differences between handling purely numerical data in either matrices or data
frames. but a factor of 17 between the obvious approach and the wise one seems
a trap in which 99% of the people will step (probably never thinking that their
might be a faster approach).
joerg
__
R-help@stat.math.ethz.ch 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] Deconvolution of a spectrum
On Fri, Mar 09, 2007 at 01:25:24PM +0100, Lukasz Komsta wrote:
>
> Dear useRs,
>
> I have a curve which is a mixture of Gaussian curves (for example UV
> emission or absorption spectrum). Do you have any suggestions how to
> implement searching for optimal set of Gaussian peaks to fit the curve?
> I know that it is very complex problem, but maybe it is a possibility
> to do it? First supposement is to use a nls() with very large functions,
> and compare AIC value, but it is very difficult to suggest any starting
> points for algotirithm.
>
> Searching google I have found only a description of commercial software
> for doing such deconvolution (Origin, PeakFit) without any information
> about used algorithms. No ready-to-use function in any language.
>
> I have tried to use a Mclust workaround for this problem, by generating a
> large dataset for which the spectrum is a histogram and feed it into
> the Mclust. The results seem to be serious, but this is very ugly and
> imprecise method.
>
> Thanks for any help,
>
> Luke
>
I would try `nls'. we have used `nls' for fitting magnetic resonance spectra
consisting of =~ 10 gaussian peaks. this works OK, if the input data are
reasonable (not too noisy, peak amplitudes above noise level, peak distance
not unreasonably smaller than peak width, i.e peak overlap such that peaks are
still more or less identifiable visually).
of course you must invest effort in getting the start values (automatically or
manually) right. if your data are good, you might get good start values for the
positions (the means of the gaussians) with an approach that was floating around
the r-help list in 11/2005, which I adopted as follows:
peaks <- function (series, span = 3, what = c("max", "min"), do.pad = TRUE,
add.to.plot = FALSE, ...)
{
if ((span <- as.integer(span))%%2 != 1)
stop("'span' must be odd")
if (!is.numeric(series))
stop("`peaks' needs numeric input")
what <- match.arg(what)
if (is.null(dim(series)) || min(dim(series)) == 1) {
series <- as.numeric(series)
x <- seq(along = series)
y <- series
}
else if (nrow(series) == 2) {
x <- series[1, ]
y <- series[2, ]
}
else if (ncol(series) == 2) {
x <- series[, 1]
y <- series[, 2]
}
if (span == 1)
return(list(x = x, y = y, pos = rep(TRUE, length(y))),
span = span, what = what, do.pad = do.pad)
if (what == "min")
z <- embed(-y, span)
else z <- embed(y, span)
s <- span%/%2
s1 <- s + 1
v <- max.col(z, "first") == s1
if (do.pad) {
pad <- rep(FALSE, s)
v <- c(pad, v, pad)
idx <- v
}
else idx <- c(rep(FALSE, s), v)
val <- list(x = x[idx], y = y[idx], pos = v, span = span,
what = what, do.pad = do.pad)
if (add.to.plot == TRUE)
points(val, ...)
val
}
this looks for local maxima in the vector ("y-values") or 2-dim array
("x/y-matrix") `series'in a neighborhood of each point defined by `span'.
if you first plot your data and then call the above on the data with
'add.to.plot = TRUE', the results of the peak search are added to your plot (and
you can modify this plotting via the `...' argument).
maybe this works for your data to get the peak position estimates (and the
amplitudes in the next step) right. frequently the standard deviations
estimates can be set to some fixed value for any given experiment.
and of course distant parts of your spectrum won't have anything to do which
each other, so you can split up the fitting to help `nls' along a bit.
joerg
__
R-help@stat.math.ethz.ch 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] reading a text file with a stray carriage return
On Thu, Mar 08, 2007 at 10:11:06AM -, Ted Harding wrote:
> On 08-Mar-07 jim holtman wrote:
> > How do you define a carriage return in the middle of a line
> > if a carriage return is also used to delimit a line? One of the
> > things you can do is to use 'count.fields' to determine the
> > number of fields in each line. For those lines that are not the
> > right length, you could combine them together with a 'paste'
> > command when you write them out.
>
>...
> cat funnyfile.csv | awk 'BEGIN{FS=","}
> {nline = NR; nfields = NF}
> {if(nfields != 3){print nline " " nfields}'
I think this might be more 'awk like':
cat funnyfile.csv | awk 'BEGIN{FS=","} NF != 3 {print NR, NF}'
(separating the "pattern" from the "action", I mean. it's slightly
faster too -- but not much). and since everyone has the same problems,
here's my fieldcounter/reporter (the actual awk program (apart from reporting
the results) is a single line, cute..).
#
#!/bin/sh
#scan input files and output statistics of number of fields.
#suitable mainly to check 'rectangular' data files for deviations
#from correct number of columns.
FS='[ \t]+' #default FS
COMMENT='#'
v1=F
v2=c
valop=$v1:$v2:
while getopts $valop option; do
if [ "$option" = "$v1" ]; then
FS="$OPTARG"
elif [ "$option" = "$v2" ]; then
COMMENT="$OPTARG"
else
echo ""
echo "usage: fieldcount [-$v1 field_separator] [-$v2 comment_char] file
..."
echo "default FS is white space"
echo "default comment_char is #"
echo ""
exit
fi
done
shift `expr $OPTIND - 1`
for i do
echo ""
echo "file $i (using FS = '$FS'):"
awk '
BEGIN {
FS = "'"$FS"'"
}
!/^'"$COMMENT"'/ {fc[NF]++}
END {
for (nf in fc) {
if (fc[nf] > 1) s1 = "s"; else s1 = ""
if (nf > 1) s2 = "s"; else s2 = ""
print fc[nf], "record" s1, "with", nf, "field" s2
}
}
' $i |sort -rn
echo ""
done
#
here, you could use this a la (it uses white space as separator per default):
fieldcount -F"," funnyfile1.csv funnyfile2.csv ...
it's handy for huge files, since it outputs cumulative results reverse sorted
by frequency of occurence (the bad guys should be infrequent and
at the bottom, contrary to real world experience...).
one could of course tweak this to report only records deviating from some
expected field count (say 3) by adding a further command line switch and using
this in the awk script as is done with the other shell variables ($FS and
$COMMENT).
joerg
>
>
> > On 3/7/07, Walter R. Paczkowski <[EMAIL PROTECTED]> wrote:
> >>
> >>
> >> Hi,
> >> I'm hoping someone has a suggestion for handling a simple problem.
> >> A client gave me a comma separated value file (call it x.csv)
> >> that has an id and name and address for about 25,000 people
> >> (25,000 records).
> >> I used read.table to read it, but then discovered that there are
> >> stray carriage returns on several records. This plays havoc with
> >> read.table since it starts a new input line when it sees the
> >> carriage return.
> >> In short, the read is all wrong.
> >> I thought I could write a simple function to parse a line and
> >> write it back out, character by character. If a carriage
> >> return is found, it would simply be ignored on the writing
> >> back out part. But how do I identify a carriage return? What is
> >> the code or symbol? Is there any easier way to rid the file
> >> of carriage returns in the middle of the input lines?
> >> Any help is appreciated.
> >> Walt Paczkowski
>
__
R-help@stat.math.ethz.ch 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] Some problems with X11
On Wed, Mar 07, 2007 at 05:21:52PM -0700, Rafael Rosolem wrote:
> Hi,
>
> I am really new with R, so I don't know anything about it. I have
> written a script (attached) which tries to do really basic stuff (such
> as computing basic statistics and basic plots). When I try to plot a
> histogram and pairs, for example, I get the following message:
>
> > source("project.R")
> Loading required package: sp
>
> -
> Analysis of geostatistical data
> For an Introduction to geoR go to http://www.est.ufpr.br/geoR
> geoR version 1.6-13 (built on 2006/12/26) is now loaded
> -
>
> Error in title(main = main, sub = sub, xlab = xlab, ylab = ylab, ...) :
> X11 font at size 8 could not be loaded
> >
>
> I have seen some threads about this problem, though none of them really
> states what should be done in order to solve that problem (At least it
> was not clear for me!).
well, some X11 font is not installed/not found.
in `R', getOpton("X11fonts") tells you which font family `X11()' is using by
default (but you can alter this in the `X11' call if need be) and the size 8
subtype seems to be missing.
you may compare the above with the output from `xlsfonts' in the unix shell,
or use `xfontsel' to browse interactively through the installed fonts
to verify whether `R's complaining is justified.
probable solution: use some existing sufficiently complete (different sizes)
font-family in the `X11()' call or copy over the default fonts from the labtop
where it works.
hth,
joerg
>
> I am also running R on a Ubuntu Linux Edgy distro. The interesting thing
> is that I have this problem on my desktop only. I also have Ubuntu Edgy
> installed on my laptop and it is working just fine!
>
> Thank you
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] R CMD CHECK question
hi,
second try...
I ran into problems when checking one of my packages with R CMD CHECK:
I have two packages, the first (named `pkc') depending on the second one (named
`roiutils'). The source code and DESCRIPTION files describes the dependency as
it should be, I think ('Imports', `require'). but if I run
R CMD CHECK pkc
I get "significant warnings" related to missing links (refering to functions
from the second package) in the manpages of the first package as can be
seen below. despite the warnings, after installing the two packages the help
system works just fine including the cross-references.
my question:
why is it, that R CMD CHECK is complaining? can one selectively switch of this
warning? or how have I to specify the links in the manpages to tell CHECK that
everything is basically OK?
CUT
* checking for working latex ... OK
* using log directory '/Users/vdh/rfiles/Rlibrary/.check/pkc.Rcheck'
* using R version 2.4.0 (2006-10-03)
* checking for file 'pkc/DESCRIPTION' ... OK
* this is package 'pkc' version '1.1'
* checking package dependencies ... OK
* checking if this is a source package ... OK
* checking whether package 'pkc' can be installed ... WARNING
Found the following significant warnings:
missing link(s): readroi readroi readroi figure readroi conv3exmodel
readroi
missing link(s): figure readroi
* checking package directory ... OK
* checking for portable file names ... OK
* checking for sufficient/correct file permissions ... OK
* checking DESCRIPTION meta-information ... OK
* checking top-level files ... OK
* checking index information ... OK
* checking package subdirectories ... OK
* checking R files for syntax errors ... OK
* checking R files for non-ASCII characters ... OK
* checking whether the package can be loaded ... OK
* checking whether the package can be loaded with stated dependencies ... OK
* checking whether the name space can be loaded with stated dependencies ... OK
* checking S3 generic/method consistency ... OK
* checking replacement functions ... OK
* checking foreign function calls ... OK
* checking R code for possible problems ... OK
* checking Rd files ... WARNING
Rd files with unknown sections:
/Users/vdh/rfiles/Rlibrary/pkc/man/fitdemo.Rd: example
See the chapter 'Writing R documentation files' in manual 'Writing R
Extensions'.
* checking Rd cross-references ... WARNING
Missing link(s) in documentation object 'compfit.Rd':
readroi readroi readroi figure readroi conv3exmodel readroi
Missing link(s) in documentation object 'exp3fit.Rd':
figure readroi
CUT
any hints appreciated,
joerg
__
R-help@stat.math.ethz.ch 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 function to its name?
On Fri, Mar 02, 2007 at 03:42:46PM +0100, Ido M. Tamir wrote:
> Hi,
>
> I can get from a string to a function with this name:
>
> >f1 <- function(x){ mean(x) }
>
> >do.call("f1",list{1:4})
> >get("f1")
> etc...
>
> But how do I get from a function to its name?
>
> >funcVec <- c(f1,median)
>
> >funcVec
> [[1]]
> function(x){ mean(x) }
> > str(funcVec)
> List of 2
> $ :function (x)
> ..- attr(*, "source")= chr "function(x){ mean(x) }"
> $ :function (x, ...)
> > deparse(funcVec[1])
> [1] "list(function (x) " "{" "mean(x)"
> [4] "})"
>
for any symbol/name
deparse(substitute(f1))
yields the string representation, but this won't give you "f1" for
deparse(substitute(funcVec[1])).
but rather the string "funcVec[1]".
if you actually want to access funcVec components via strings denoting the
components, maybe you simply could use
funcVec <- list(f1 = f1, median = median)
and access these via
funcVec[["f1"]]
which would enable using a string variable holding the name:
dum = "f1"
funcVec[[dum]]
hth
joerg
__
R-help@stat.math.ethz.ch 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] Double-banger function names: preferences and suggestions
On Thu, Mar 01, 2007 at 10:31:08AM -0700, Jason Barnhart wrote: > Definitely not #2. Prefer #1 but #3 is ok as well. Definitely not #1. Prefer #2 but #3 is ok as well. there is nothing like unanimous agreement :-) but in earnest: I don't think #1 is good on the already mentioned grounds that it masks the difference between a method call and a simple function name. I don't think a function name such as plot.mydata would be reasonable: is this a method call for plotting objects of class `mdyata' or is it some other function call (probably plotting my data)? in browsing through source code this is at least an unnecessary nuisance having to check this point. > > Thanks for contributing and inquiring. > > > - Original Message - > From: "hadley wickham" <[EMAIL PROTECTED]> > To: <[EMAIL PROTECTED]> > Sent: Sunday, February 25, 2007 7:44 AM > Subject: [R] Double-banger function names: preferences and suggestions > > > > What do you prefer/recommend for double-banger function names: > > > > 1 scale.colour > > 2 scale_colour > > 3 scaleColour > > > > 1 is more R-like, but conflicts with S3. 2 is a modern version of > > number 1, but not many packages use it. Number 3 is more java-like. > > (I like number 2 best) > > > > Any suggestions? > > > > Thanks, > > > > Hadley > > > > __ > > R-help@stat.math.ethz.ch mailing list > > https://stat.ethz.ch/mailman/listinfo/r-help > > PLEASE do read the posting guide > > http://www.R-project.org/posting-guide.html > > and provide commented, minimal, self-contained, reproducible code. > > > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch 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] str() to extract components
On Tue, Feb 27, 2007 at 03:52:54PM -, Simon Pickett wrote: > Hi, > > I have been dabbling with str() to extract values from outputs such as > lmer etc and have found it very helpful sometimes. > > but only seem to manage to extract the values when the output is one > simple table, any more complicated and I'm stumped :-( > > take this example of the extracted coeficients from a lmer analysis... > > using str(coef(lmer(resp3~b$age+b$size+b$pcfat+(1|sex), data=b))) yields > > Formal class 'lmer.coef' [package "Matrix"] with 3 slots > ..@ .Data :List of 1 > .. ..$ :`data.frame': 2 obs. of 4 variables: > .. .. ..$ (Intercept): num [1:2] 1.07 1.13 > .. .. ..$ b$age : num [1:2] 0.00702 0.00702 > .. .. ..$ b$size : num [1:2] 0.0343 0.0343 > .. .. ..$ b$pcfat: num [1:2] 0.0451 0.0451 > ..@ varFac: list() > ..@ stdErr: num(0) > > how do I "get inside" the first table to get the value 1.07 for instance? > > Any help much appreciated. > may `unlist' would be enough? __ R-help@stat.math.ethz.ch 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 with an exact size
On Fri, Feb 23, 2007 at 04:24:54PM -0200, Alberto Monteiro wrote:
> Is it possible to create a pdf output file with an (as nearly as
> possible) exact size?
>
> For example, if I want to draw in an A4 paper (210 x 297 mm) a
> square of 100 x 100 mm, how can I do it?
>
> FWIW, about 6 months ago I learned here how to create an exact
> png image. For example, if I want a 500 x 500 black square in
> a 1000 x 1000 white png, occupying the center of the png, the
> procedure is this:
>
> png("image.png", width=1000, height=1000, bg="white")
> par(mar=c(0,0,0,0)) # reset margins
> plot(0, xlim=c(0, 999), ylim=c(0, 999), col="white")
> par(usr=c(0, 999, 0, 999))
> points(c(250, 250, 749, 749, 250), c(250, 749, 749, 250, 250),
> type="l", col="black")
> dev.off()
>
> However, I don't know how do this with a pdf monstr... oops... file.
>
> Alberto Monteiro
>
going via `bitmap' and using the `pdfwrite' type is probably the better
idea since in this case `ghostscript` is used for pdf generation which seems
over all to generate cleaner/better pdf-output in comparsion
to the R internal pdf-device (I believe there was recently some
discussion of this on the list?).
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] R CMD CHECK question
hi,
I have two private packages, the first (`pkc') depending on the second one
(`roiutils'). The source code and DESCRIPTION files describes the dependency
as it should be ('Imports', `require'), at least I think so.
now, running
R CMD CHECK pkc
yields the following output in which I have inserted my
questions (lines starting with -->):
* checking for working latex ... OK
* using log directory '/Users/vdh/rfiles/Rlibrary/.check/pkc.Rcheck'
* using R version 2.4.0 (2006-10-03)
* checking for file 'pkc/DESCRIPTION' ... OK
* this is package 'pkc' version '1.1'
* checking package dependencies ... OK
* checking if this is a source package ... OK
* checking whether package 'pkc' can be installed ... WARNING
Found the following significant warnings:
missing link(s): readroi readroi readroi figure readroi conv3exmodel
readroi
missing link(s): figure readroi
--> there _are_ links to the mentioned functions (from `roiutils') in the
--> manpages of `pkc'. after installing the libs, the help system works just
--> fine. why is it, that CHECK is complaining? can one selectively switch of
--> this warning? or how have I to specify the links to tell CHECK that
--> everything is OK?
* checking package directory ... OK
* checking for portable file names ... OK
* checking for sufficient/correct file permissions ... OK
* checking DESCRIPTION meta-information ... OK
* checking top-level files ... OK
* checking index information ... OK
* checking package subdirectories ... OK
* checking R files for syntax errors ... OK
* checking R files for non-ASCII characters ... OK
* checking whether the package can be loaded ... OK
* checking whether the package can be loaded with stated dependencies ... OK
* checking whether the name space can be loaded with stated dependencies ... OK
* checking S3 generic/method consistency ... OK
* checking replacement functions ... OK
* checking foreign function calls ... OK
* checking R code for possible problems ... OK
* checking Rd files ... WARNING
Rd files with unknown sections:
/Users/vdh/rfiles/Rlibrary/pkc/man/fitdemo.Rd: example
See the chapter 'Writing R documentation files' in manual 'Writing R
Extensions'.
* checking Rd cross-references ... WARNING
Missing link(s) in documentation object 'compfit.Rd':
readroi readroi readroi figure readroi conv3exmodel readroi
Missing link(s) in documentation object 'exp3fit.Rd':
figure readroi
--> this seems the same problem as above, right?
any hints appreciated,
joerg
__
R-help@stat.math.ethz.ch 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] Confindence interval for Levenberg-Marquardt fit
On Wed, Feb 21, 2007 at 09:41:29AM -0600, Douglas Bates wrote: > On 2/21/07, joerg van den hoff <[EMAIL PROTECTED]> wrote: > >On Wed, Feb 21, 2007 at 11:09:52AM +, Prof Brian Ripley wrote: > >> Well, the algorithm used does not affect the confidence interval > >(provided > >> it works correctly), but what is nls.ml (presumably in some package you > >> have not mentioned) and why would I want to use an old-fashioned > >> algorithm? > > >is'nt this a bit strong? in what respect do you consider > >levenberg-marquardt > >(going back to the 19-forties, I think) as old-fashioned (especially > >in comparsion to the `nls' standard gauss-newton approach (both gentlemen > >seem to have done their major work a bit longer ago :-))? > >AFAICS levenberg-marquardt is generally appreciated for it's rapid > >convergence achieved by a smooth transition from an > >inverse-hessian approach to steepest descent. my own experience > >with non-linear least squares minimization using this algorithm > >are positive as well, but > >I have not tried out the levenberg-marquardt > >implementation in package `minpack.lm' (which originates from netlib.org) > >and don't know if it's good. but in any case there sure are implementations > >around (e.g. in the CERN MINUIT library) which have proven to be > >of high quality. > > The nls function provides implementations of the Gauss-Newton > algorithm and the Golub-Pereyra algorithm and the NL2SOL algorithm. > Check the possible values of the optional argument "algorithm" to > nls(). I'm using `nls' rather frequently and I'm aware of this. maybe the formulation "standard gauss-newton approach" was misleading (please attribute this to the fact that I'm not a native speaker: I actually meant the default algorithm used and was not implying that `nls' can only use some odd "standard" procedure) > > If you will allow me to wax philosophical for a moment, I submit that > there are three stages to any iterative optimization: > - formulating initial values of the parameters > - determining what step to take from iteration i to iteration i+1 > - determining when to declare convergence > > For the most part optimization research focuses on the second step. > The first and third steps are also important. One of the advantages > of the implementation of the Gauss-Newton and Golub-Pereyra algorithms > in nls is that they use a convergence criterion (based only on the > current parameter values) not a termination criterion (based on the > changes in the parameter values or the changes in the objective > function or both). See chapter 2 of Bates and Watts, "Nonlinear > Regression Analysis and Its Applications" (Wiley, 1988) for a > discussion of the particular convergence criterion used. I did this when starting to use `nls' seriously. being a physicist, not a mathematician, my attitude is nevertheless a bit more pragmatic (or simple-minded, if you like): I fully appreciate that these things need a thorough and rigorous mathematical foundation but any convergence criterion suits me just fine as long as (1) convergence is actually achieved (i.e. the algorithm reports some results...) and (2) the final least squares sum is sufficiently near the achievable minimum value for that model (sufficiently near meaning that the parameter estimation errors are much larger than any further changes in the parameters when wandering around any longer near the current point in parameter space -- possibly finding a slightly "better mininmum") and (3) the results are obtained in sufficiently short time. It's not obivous to me (I did'nt encounter an example case up to now), for instance, that a good termination criterion is either more unreliable (yielding wrong answers way off the true minimum) or less efficient (e.g. driving the iteration count up) than the convergence criterion used in `nls' (against which I do not have any objections -- apart from the minor handicap that it can't converge on 'ideal' data which by and then comes in the way as the help list shows) an impression which I cannot prove by hard data (so it may be wrong) is, that the `nls' gauss-newton approach tends to be more sensitive to suboptimal start values (sending it towards outer space in the worst case...) than levenberg-marquardt driven approaches and definitely has on average a higher iteration count. > > Another point that is often lost in a discussion of various > optimization algorithms is that one does not compare algorithms - one > compares implementations of algorithms which, at the very least, must > involve all of the above steps. sure > > For many p
Re: [R] Confindence interval for Levenberg-Marquardt fit
On Wed, Feb 21, 2007 at 11:09:52AM +, Prof Brian Ripley wrote: > Well, the algorithm used does not affect the confidence interval (provided > it works correctly), but what is nls.ml (presumably in some package you > have not mentioned) and why would I want to use an old-fashioned > algorithm? is'nt this a bit strong? in what respect do you consider levenberg-marquardt (going back to the 19-forties, I think) as old-fashioned (especially in comparsion to the `nls' standard gauss-newton approach (both gentlemen seem to have done their major work a bit longer ago :-))? AFAICS levenberg-marquardt is generally appreciated for it's rapid convergence achieved by a smooth transition from an inverse-hessian approach to steepest descent. my own experience with non-linear least squares minimization using this algorithm are positive as well, but I have not tried out the levenberg-marquardt implementation in package `minpack.lm' (which originates from netlib.org) and don't know if it's good. but in any case there sure are implementations around (e.g. in the CERN MINUIT library) which have proven to be of high quality. `nls' sure is a _very_ valuable function, but not necessarily the "last word" with regards to the chosen algorithm(s). > > You could start nls at the solution you got from nls.ml and use confint() > on that. maybe one should look at profile.nls and confint.nls and see what information of the usual `nls' object is actually used for the confidence intervall computation and mimick this for the `nls.lm' output? at a (admittedly) quick glance it seems that only parameters, std.errs. and the fitted/residual values are needed which should all be provided by nls.lm as well. maybe one could even try to map the nls.lm results into a structure of class `nls' (although this would not be a clean solution, probably) in order to use `confint.nls'? > > On Wed, 21 Feb 2007, Michael Dondrup wrote: > > > Dear all, > > I would like to use the Levenberg-Marquardt algorithm for non-linear > > least-squares regression using function nls.lm. Can anybody help me to > > find a a way to compute confidence intervals on the fitted > > parameters as it is possible for nls (using confint.nls, which does not > > work for nls.lm)? > > > > Thank you for your help > > Michael > > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] "[[" gotcha
Robin Hankin wrote: > The following gotcha caught me off-guard just now. > > I have two matrices, a and b: > > > a <- matrix(1,3,3) > b <- matrix(1,1,1) > > (note that both "a" and "b" are matrices). > > I want them in a list: > > > B <- NULL > > B[[1]] <- a > > B[[2]] <- b > > B > [[1]] > [,1] [,2] [,3] > [1,]111 > [2,]111 > [3,]111 > > [[2]] > [,1] > [1,]1 > > > > > This is fine. > > But swapping "a" and "b" over does not behave as desired: > > > > B <- NULL > > B[[1]] <- b > > B[[2]] <- a > Error in B[[2]] <- a : more elements supplied than there are to replace > > > > > > The error is given because after B[[1]] <- a, the variable B is > just a scalar and > not a matrix (why is this?) > > What's the bulletproof method for assigning matrices to a list (whose > length is > not known at runtime)? > > not sure about "bulletproof", but: you should tell R that B is really intended to be a list in the first place: B <- list() the rest then works as you intended. whether the 'simplification' of your 1x1 matrix to a scalar in your example is canonical (and desirable) behaviour seems a question for some of the experts (it's a bit reminiscent of the `drop = TRUE' vs. `drop = FALSE' problem) joerg > > > > > > > -- > Robin Hankin > Uncertainty Analyst > National Oceanography Centre, Southampton > European Way, Southampton SO14 3ZH, UK > tel 023-8059-7743 > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] sapply problem
Sundar Dorai-Raj wrote:
>
>
> Joerg van den Hoff said the following on 12/14/2006 7:30 AM:
>> I have encountered the following problem: I need to extract from
>> a list of lists equally named compenents who happen to be 'one row'
>> data frames. a trivial example would be:
>>
>> a <- list(list(
>> df = data.frame(A = 1, B = 2, C = 3)), list(df = data.frame(A = 4,B =
>> 5,C = 6)))
>>
>> I want the extracted compenents to fill up a matrix or data frame row
>> by row.
>> the obvious thing to do seems:
>>
>> b <- sapply(a, "[[", "df")
>> b <- t(b)
>>
>> now `b' looks all right:
>>
>> b
>> class(b)
>>
>> but it turns out that all elements in this matrix are one element lists:
>>
>> class(b[1,1])
>>
>> which prevents any further standard processing of `b' (like
>> `colMeans', e.g.)
>>
>> question 1: is their a straightforward way to enforce that `b' contains
>> simple numbers as elements right from the start (instead of something
>> like
>> apply(b, 1:2, "class<-", "numeric") afterwards)?
>>
>
> Try this:
>
> a <- list(list(df = data.frame(A = 1, B = 2, C = 3)),
> list(df = data.frame(A = 4, B = 5, C = 6)))
> b <- do.call("rbind", sapply(a, "[", "df"))
> b
>
>
>> question 2: should not sapply do this further 'simplification' anyway
>> in a situation
>> like this (matrix elements turn out to be one-element lists)?
>>
>
> I think it does as it much as it knows how. I think you might believe
> that matrix elements can only contain numeric values. This is not a
> valid assumption. Take this example:
>
> > a <- list(1)
> > b <- list(2)
> > (m <- matrix(c(a, b), 2, 1))
> [,1]
> [1,] 1
> [2,] 2
> > class(m[1, 1])
> [1] "list"
> > class(m[2, 1])
> [1] "list"
>
> HTH,
>
> --sundar
>
>> regards
>>
>> joerg
>>
thanks for the help to everybody. the proposal of sundar seems the most concise
one (if it
is also efficient I'll see with larger data sets ..).
with regards to my question no. 2: sure I realize that I _can_ have a sensible
matrix
consisting of elements that are lists. with regards to `sapply' and the
intention (as I
understand it) to simplify as much as is sensible possible (`sapply' unlists
usually
everything it can), I think it would be desirable if `sapply' would detect this
admittedly
very special case: "a matrix of one-element lists whose sole list elements are
atomic of
length 1" (which was my case and is exactly your example above).
of course, I don't know whether such a special treatment (probably meaning
further
recursive tests in `sapply') would slow down `sapply' significantly in general
(which
would be an argument against such a change of `sapply') or maybe it is against
the actual
purpose of sapply -- I'm not sure.
anyway, thanks again,
joerg
__
R-help@stat.math.ethz.ch 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] Better way to change the name of a column in a dataframe?
Ben Fairbank wrote: > Hello R users -- > > > > If I have a dataframe such as the following, named "frame" with the > columns intended to be named col1 through col6, > > > >> frame > > col1 col2 cmlo3 col4 col5 col6 > > [1,]3 10 2657 > > [2,]68 4 1071 > > [3,]75 1318 > > [4,] 106 5492 > > > > and I want to correct or otherwise change the name of one of the > columns, I can do so with > > > >> dimnames(frame)[[2]][which(dimnames(frame)[[2]]=="cmlo3")] <- "col3" > > > > which renames the offending column: > > > >> frame > > col1 col2 col3 col4 col5 col6 > > [1,]3 102657 > > [2,]684 1071 > > [3,]751318 > > [4,] 1065492 > > > > This seems cumbersome and not very intuitive. How can one accomplish > this more simply? > > well I would simply use names(frame)[3] = 'col3' (supposing you know the column number of your offending column anyway). __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Delete all dimnames
Serguei Kaniovski wrote: > Hello, how can I get rid of all dimnames so that: > $amat > Var3 Var2 Var1 > 8 1111 1 1 1 1 0 0 0 0 0 0 0 > 7 1110 1 0 0 0 1 0 0 0 0 0 0 > 6 1101 0 1 0 0 0 1 0 0 0 0 0 > 5 1100 0 0 0 0 0 0 1 0 0 0 0 > 4 1011 0 0 1 0 0 0 0 1 0 0 0 > 3 1010 0 0 0 0 0 0 0 0 1 0 0 > 2 1001 0 0 0 0 0 0 0 0 0 1 0 > 1 1000 0 0 0 0 0 0 0 0 0 0 1 > > is displayed with [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] in rows, and the > same in columns. The matrix was generated using "apply" > > dimnames(mat) <- NULL __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] sapply problem
I have encountered the following problem: I need to extract from a list of lists equally named compenents who happen to be 'one row' data frames. a trivial example would be: a <- list(list( df = data.frame(A = 1, B = 2, C = 3)), list(df = data.frame(A = 4,B = 5,C = 6))) I want the extracted compenents to fill up a matrix or data frame row by row. the obvious thing to do seems: b <- sapply(a, "[[", "df") b <- t(b) now `b' looks all right: b class(b) but it turns out that all elements in this matrix are one element lists: class(b[1,1]) which prevents any further standard processing of `b' (like `colMeans', e.g.) question 1: is their a straightforward way to enforce that `b' contains simple numbers as elements right from the start (instead of something like apply(b, 1:2, "class<-", "numeric") afterwards)? question 2: should not sapply do this further 'simplification' anyway in a situation like this (matrix elements turn out to be one-element lists)? regards joerg __ R-help@stat.math.ethz.ch 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 generate the random numbers uniformly distributed on the unit disc?
S.Q. WEN wrote:
> Hi,
> I want to get random number which is uniformly distributed on the unit
> disc.
> How can I do that with R?
>
>
> Best wishes,
> WAN WAN
>
> [[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
the following function should put "N" uniformly distributed points into
a disc of radius "R" (and plot the points in the (phi,rho) and (xx, yy)
coordinate systems).
disc <- function (N = 10^4, R = 1) {
phi <- runif(N) * 2 * pi
rho <- R * sqrt(runif(N))
xx <- rho * cos(phi)
yy <- rho * sin(phi)
layout(1:2)
plot(phi, rho, pch = ".")
plot(xx, yy, pch = ".", asp = 1)
layout(1)
invisible(list(phi = phi, rho = rho))
}
the trick is to transform a uniform distribution along `rho' in such a
way that it gets "denser" towards the rim in such a way that the
'smearing out' of the points over the circumference is compensated.
regards
joerg
__
R-help@stat.math.ethz.ch 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] extremely slow recursion in R?
Thomas Lumley wrote:
> On Thu, 24 Aug 2006, Jason Liao wrote:
>
>> I recently coded a recursion algorithm in R and ir ran a few days
>> without returning any result. So I decided to try a simple case of
>> computing binomial coefficient using recusrive relationship
>>
>> choose(n,k) = choose(n-1, k)+choose(n-1,k-1)
>>
>> I implemented in R and Fortran 90 the same algorithm (code follows).
>> The R code finishes 31 minutes and the Fortran 90 program finishes in 6
>> seconds. So the Fortran program is 310 times faster. I thus wondered if
>> there is room for speeding up recursion in R. Thanks.
>>
>
> Recursive code that computes the same case many times can often be sped up
> by memoization, eg
>
> memo<-new.env(hash=TRUE)
> chewse<-function(n,k) {
> if (n==k) return(1)
> if(k==1) return(n)
>
> if(exists(paste(n,k),memo,inherits=FALSE))
> return(get(paste(n,k),memo))
> rval<-chewse(n-1,k)+chewse(n-1,k-1)
> assign(paste(n,k),rval,envir=memo)
> return(rval)
> }
>
> This idea was discussed in an early Programmers' Niche article by Bill
> Venables in R News.
>
> However, I'm surprised that you're surprised that compiled Fortran 90 is
> 310 times faster than interpreted R. That would be about what I would
> expect for code that isn't making use of vectorized functions in R.
>
>
> -thomas
>
maybe someone's interested:
I made the same observation of seemingly very slow recursion recently:
just for fun I used the (in)famously inefficient
fib <- function(n = 1) {
if (n < 2)
fn <- 1
else
fn <- fib(n - 1) + fib(n - 2)
fn
}
for calculating the fibonacci numbers and compared `fib(30)' (about
1.3e6 recursive function calls ...) to some other languages (times in sec):
language time
==
C 0.034 (compiled, using integers)
Ocaml 0.041 (compiled, using integers)
Ocaml 0.048 (interpreted, using integers)
C 0.059 (compiled, using floats)
Lua1.1
ruby 3.4
R 21
octave120
apart from octave (which seems to have a _real_ issue with recursive
function calls), R is by far the slowest in this list and still a factor
7-20 slower than the interpreter based Lua and ruby. the speed loss
compared to C or Ocaml is about a factor of 350-600 here (Ocaml keeps
its speed more or less in this simple example even in 'script mode',
which is remarkable, I think (and it usually loses only a factor of
about 7 or so in script mode compared to the compiled variant)
for the specialists the bad performance of R in this situation might not
be surprising, but I was not aware that recursive function calls are
seemingly as expensive as explicit loops (where the execution time ratio
of R to C again is of the same order, i.e. =~ 400).
of course, such comparsions don't make too much sense: the reason to use
R will definitely not be its speed (which, luckily, often does not
matter), but the combination of flexible language, the number of
available libraries and the good 2D graphics.
joerg
__
R-help@stat.math.ethz.ch 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] Intro to Programming R Book
Raphael Fraser wrote: > I am new to R and am looking for a book that can help in learning to > program in R. I have looked at the R website suggested books but I am > still not sure which book best suite my needs. I am interesting in > programming, data manipulation not statistics. Any suggestions? > > Raphael > > __ > R-help@stat.math.ethz.ch 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. "S Programming" by Venables and Ripley (Springer) seems the only(?) one around targeting the language, not it's applications. luckily, it's very good. for the rest (things specific to R, e.g. package development, namespaces etc.) I think one can only resort to the R manuals . __ R-help@stat.math.ethz.ch 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] nls convergence problem
Earl F. Glynn wrote: > "Berton Gunter" <[EMAIL PROTECTED]> wrote in message > news:[EMAIL PROTECTED] >>> Or, maybe there's something I don't understand about the >>> algorithm being used. >> Indeed! So before making such comments, why don't you try to learn about >> it? >> Doug Bates is a pretty smart guy, and I think you do him a disservice >> when >> you assume that he somehow overlooked something that he explicitly warned >> you about. I am fairly confident that if he could have made the problem go >> away, he would have. So I think your vent was a bit inconsiderate and >> perhaps even intemperate. The R Core folks have produced a minor miracle >> IMO, and we should all be careful before assuming that they have >> overlooked >> easily fixable problems. They're certainly not infallible -- but they're a >> lot less fallible than most of the rest of us when it comes to R. > > I meant no disrespect to Doug Bates or any of the R Core folks. I thought > what I wrote had a "neutral" tone and was respectful. I am sorry if anyone > was offended by my comments and suggestions. I am certainly thankful for > all the hard work that has gone into developing R. > > efg > well, just a feedback to that: of course the tone of your mail was by no means inadequate (at least douglas bates did not object...). the tendency on this list to rather harshly rebuke people for some kind of (real or imagined) misconduct and to 'defend' R against 'attacks' is counterproductive and unnecessary. it goes without saying that the people 'behind' R can not and will not (and are not expected to) change the code after each mail on the help list which raises some question. and concerning your `nls' question: sure, the noise requirement is a pitfall in the beginning, but afterwards it's irrelevant: you don't fit noise free data in real life (in the sense that real data never follow you model exactly). and, sure, the convergence decision could be altered (given enough time and knowledge). whether convergence failure on exact data is a bug or a feature is a matter of taste, probably. getting better access to the trace output and especially access to intermediate pre-convergence values of the model parameters (this would 'solve' your problem, too) would really be an improvement, in my mind (I think this is recognized by d. bates, but simply way down his 'to do' list :-(). joerg __ R-help@stat.math.ethz.ch 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] least square fit with non-negativity constraints for absorption spectra fitting
Xu, Xiuli (NIH/NHLBI) [E] wrote: > I would really appreciate it if someone can give suggestions on how to > do spectra fitting in R using ordinary least square fitting and > non-negativity constraints. The lm() function works well for ordinary > least square fitting, but how to specify non-negativity constraints? It > wouldn't make sense if the fitting coefficients coming out as negative > in absorption spectra deconvolution. > > Thanks. > > Xiuli > I'm not sure, but would presume that constraints could not be imposed on a linear least squares fit. maybe someone can correct me. if you move to `nls', i.e. non-linear least squares fitting, you should be able to transform your model function. say, you want some parameter `a' to stay positive. then you could e.g. substitute `a = exp(b)' in the model function and fit `b' without constraints in the "new" model and calculate `a' afterwards (which obviously is guaranteed now to be positive). note that error estimates would than have to be computed by gaussian error propagation from `b' to `a'. joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] colors on graph
Jim Lemon wrote: > COMTE Guillaume wrote: >> Hy all, >> >> >> >> I need to draw something in 2 dimension that has 3 dimension, the choice >> has been made to use colors to display the third dimension into the >> graph. >> >> >> >> Has someone done something like that, i can't figure out how to >> parametize the colors. >> > Hi Guilllaume, > > Have a look at color.scale in the plotrix package. > > Jim > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html ?image ?filled.contour and the "See Also"s therein __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] simple question about variables....
Stéphane Cruveiller wrote: > Dear R users, > > I have a simple question on variable manipulation. > Imagine I have an object "OBJ" that has "toto" as one of its variables. > I would like to understand why if I do > > > varname <- "toto" > > >OBJ$varname returns no results > > whereas > > > OBJ[varname]returns the column entitled > "toto" > > > Thanks for your help. > > Stéphane. > because if the value of `varname' is substituted in the expressions, in the first case that yields OBJ$"toto" and in the second OBJ["toto"] the latter is valid, the former is not (you'd need `OBJ$toto' there), read ` ?"$" ': "...Both '[[' and '$' select a single element of the list. The main difference is that '$' does not allow computed indices, whereas '[[' does. 'x$name' is equivalent to 'x[["name"]]'..." not, too, the difference between `[' (sublist) and `[[' (single element extraction) joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Optional variables in function?
jim holtman wrote:
> ?missing
>
> On 7/2/06, Jonathan Greenberg <[EMAIL PROTECTED]> wrote:
>> I'm a bit new to writing R functions and I was wondering what the "best
>> practice" for having optional variables in a function is, and how to test
>> for optional and non-optional variables? e.g. if I have the following
>> function:
>>
>> helpme <- function(a,b,c) {
>>
>>
>> }
>>
>> In this example, I want c to be an optional variable, but a and b to be
>> required. How do I:
>> 1) test to see if the user has inputted c
>> 2) break out of the function of the user has NOT inputted a or b.
if(missing(c))
stop("need c")
or
if(missing(c))
return()
depending on your problem it might be better to provide sensible default
values for a,b,c in the function definition
helpme<-function(a=A, b=B, c=C) {...}
instead of enforcing that every argument is specified?
joerg
>>
>> Thanks!
>>
>> --j
>>
>> --
>>
>> Jonathan A. Greenberg, PhD
>> NRC Research Associate
>> NASA Ames Research Center
>> MS 242-4
>> Moffett Field, CA 94035-1000
>> Phone: 415-794-5043
>> AIM: jgrn3007
>> MSN: [EMAIL PROTECTED]
>>
>> __
>> R-help@stat.math.ethz.ch mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide!
>> http://www.R-project.org/posting-guide.html
>>
>
>
>
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Help needed understanding eval,quote,expression
Prof Brian Ripley wrote:
> You are missing eval(parse(text=)). E.g.
>
>> x <- list(y=list(y1="hello",y2="world"),z=list(z1="foo",z2="bar"))
> (what do you mean by the $ at the start of these lines?)
>> eval(parse(text="x$y$y1"))
> [1] "hello"
>
> However, bear in mind
>
>> fortune("parse")
>
> If the answer is parse() you should usually rethink the question.
> -- Thomas Lumley
>R-help (February 2005)
>
> In your indicated example you could probably use substitute() as
> effectively.
>
>
> On Wed, 28 Jun 2006, [EMAIL PROTECTED] wrote:
>
>> I am trying to build up a quoted or character expression representing a
>> component in a list in order to reference it indirectly.
>> For instance, I have a list that has data I want to pull, and another list
>> that has character vectors and/or lists of characters containing the names
>> of the components in the first list.
>>
>>
>> It seems that the way to do this is as evaluating expressions, but I seem
>> to be missing something. The concept should be similar to the snippet
>> below:
>>
>>
>> For instance:
>>
>> $x = list(y=list(y1="hello",y2="world"),z=list(z1="foo",z2="bar"))
>> $y = quote(x$y$y1)
>> $eval(y)
>> [1] "hello"
>>
>>
>> but, I'm trying to accomplish this by building up y as a character and
>> then evaluating it, and having no success.
>>
>> $y1=paste("x$y$","y1",sep="")
>> $y1
>> [1] "x$y$y1"
>>
>>
>> How can I evaluate y1 as I did with y previously? or can I?
>>
>>
>> Much Thanks !
>>
>>
if I understand you correctly you can achieve your goal much easier than
with eval, parse, substitute and the like:
x <- list(y=list(y1="hello",y2="world"),z=list(z1="foo",z2="bar"))
s1 <- 'y'
s2 <- 'y1'
x[[s1]][[s2]]
i.e. using `[[' instead of `$' for list component extraction allows to
use characters for indexing (in other words: x$y == x[['y']])
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R on MAC OS X
Sara Mouro wrote:
>
>> Dear all,
>>
>> I have been usig R for some time, but now I have a MAC instead of a
>> PC, am I am having problems in reading files...
>>
>>
>> I have tried:
>> Data<-read.table("Users/SaraMM/PhD/Analises-LitterBags/Dados-
>> Litter.txt",head=T)
>>
>> but it said:
>> Error in file(file, "r") : unable to open connection
>> In addition: Warning message:
>> cannot open file 'Users/SaraMM/PhD/Analises-LitterBags/Dados-
>> Litter.txt', reason 'No such file or directory'
>>
>>
>> I guess that might be simple...
>> ... related to the fact that in PCs we do "C:/"
>> but in iMAC, the only "path" I found was that one: "Users/
>> SaraMM(...)"...
>>
>> Can someone please help me on this?
>>
>> Wishes,
>> Sara
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
the (absolute) path to your file has to start with `/', i.e. your file
should be
"/Users/SaraMM/PhD/Analises-LitterBags/Dados-Litter.txt"
^
|
the leading slash was missing in your read.table call (in the original
way it is a relative path name and would only work, if your current
directory within R were `/' itself)
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] About filled.contour function
Prodromos Zanis wrote: > Dear R-projects users > > I would like like to ask if there is any way to produce a multipanel plot > with the filled.contour function. In the help information of filled > contour it is said that this function is restricted to a full page > display. > > With kind regards > > Prodromos Zanis > > > > `filled.contour' sets up a 2-by-1 grid (for colormap and image) using the `layout' function, hence a prior setup of a multipanel layout would be (and is) destroyed by the `filled.contour' call. you might edit a copy of `filled contour' to modify the `layout' call (e.g. to set up a 2-by-2 grid where only the upper row is used by the plots generated in filled contour). I tried this very shortly before answering but it seems that only 'within' the filled.contour copy one has access to the other subplots, after return from the function plot focus is again in subplot no. one --no time to check this further. in any case along this line one should be able to enforce multiple subplots where 2 of them are used by the modified `filled.contour' joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] problem with "code/documentation mismatch"
I have a package with a division method for special objects of class
"RoidataList", i.e. the function is named `/.RoidataList'.
documentation for this is in the file "Divide.RoidataList".
R CMD CHECK complains with:
=cut==
* checking for code/documentation mismatches ... WARNING
Functions/methods with usage in documentation object
'Divide.RoidataList' but not in code:
/
* checking Rd \usage sections ... WARNING
Objects in \usage without \alias in documentation object
'Divide.RoidataList':
/
=cut==
the `usage' section in the Rd file reads
=cut==
\usage{
x/y
}
=cut==
which, of course is the desired way to use the function.
what am I doing wrong, i.e. how should I modify the Rd file?
maybe obvious, but not to me.
joerg van den hoff
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Basic package structure question
Prof Brian Ripley wrote: > On Fri, 23 Jun 2006, Joerg van den Hoff wrote: > >> just to confirm duncan murdochs remark: >> >> our Windows machines lack proper development environments (mainly >> missing perl is the problem for pure R-code packages, I believe?) and we >> bypass this (for pure R-code packages only, of course) by >> >> 1.) install the package on the unix machine into the desired R library >> 2.) zip the _installed_ package (not the source tree!) found in the R >> library directory >> 3.) transfer this archive to the Windows machine >> 4.) unzip directly into the desired library destination >> >> this procedure up to now always worked including properly installed >> manpages (text + html (and I hope this remains the case in the future...) > > From README.packages: > > If your package has no compiled code it is possible that zipping up the > installed package on Linux will produce an installable package on > Windows. (It has always worked for us, but failures have been reported.) > > so this is indeed already documented, and does not work for everyone. > albeit sparsely... AFAIKS it's not in the `R Extensions' manual at all. If(!) this approach could be made an 'official' workaround and explained in the manual, that would be good. I'd appreciate if someone could tell me: are the mentioned failures confirmed or are they "UFOs"? if so, are (the) reasons for (or circumstances of) failure known (I'm always afraid walking on thin ice when using this transfer strategy)? what does "produce an installable package on Windows" in the README text mean? I presume it does not mean that R CMD INSTALL (or the Windows equivalent) does work? if it really means "unzip the package on the Windows machine into the library directory", should'nt the text be altered? and I forgot to mention in my first mail: I use the described procedure for transfer from a non-Intel apple machine under MacOS (a FreeBSD descendant) to Windows (and even (unneccessarily, I know) to Sun/Solaris). so the strategy is not restricted to transfer from Linux -> Windows. and it is useful (if it is not 'accidental' that it works at all): in this way one can keep very easily in sync several local incarnations of a package across platforms (if network access to a common disk is not the way to go): I simply `rsync' the affected (local, R-code only) library directories. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Basic package structure question
Gabor Grothendieck wrote: > On 6/22/06, Duncan Murdoch <[EMAIL PROTECTED]> wrote: >> Jay Emerson wrote: >>> At the encouragement of many at UseR, I'm trying to build my first real >>> package. I have no C/Fortran code, just plain old R code, so it should be >>> rocket science. On a Linux box, I used package.skeleton() to create a basic >>> package containing just one "hello world" type of function. I edited the >>> DESCRIPTION file, changin the package name appropriately. I edited the >>> hello.Rd file. Upon running R CMD check hello, the only warning had to do >>> with the fact that src/ was empty (obviously I had no source in such a >>> simple package). I doubt this is a problem. >>> >>> I was able to install and use the package successfully on the Linux system >>> from the .tar.gz file, so far so good! Next, on to Windows, where the >>> problem arose: >>> >>> I created a zip file from inside the package directory: zip -r ../hello.zip >>> ./* >>> >>> >> Which package directory, the source or the installed copy? I think this >> might work in the installed copy, but would not work on the source. >> It's not the documented way to build a binary zip, though. >>> When I moved this to my Windows machine and tried to install the package, I >>> received the following error: >>> >>> >>>> utils:::menuInstallLocal() >>>> >>> Error in unpackPkg(pkgs[i], pkgnames[i], lib, installWithVers) : >>> malformed bundle DESCRIPTION file, no Contains field >>> >>> I only found one mention of this in my Google search, with no reply to the >>> thread. The Contains field appears to be used for bundles, but I'm trying >>> to create a package, not a bundle. This leads me to believe that a simple >>> zipping of the package directory structure is not the correct format for >>> Windows. >>> >>> Needless to say, there appears to be wide agreement that making packages >>> requires precision, but fundamentally a package should (as described in the >>> documentation) just be a collection of files and folders organized a certain >>> way. If someone could point me to documentation I may have missed that >>> explains this, I would be grateful. >>> >> I think the "organized in a certain way" part is actually important. >> Using R CMD install --build is the documented way to achieve this. It's >> not trivial to do this on Windows, because you need to set up a build >> environment first, but it's not horribly difficult. >> >> Duncan Murdoch >>> Regards, >>> >>> Jay > > One idea that occurred to me in reading this would be to have a server > that one can send a package to and get back a Windows build to > eliminate having to set up a development environment. Not sure if > this is feasible, particularly security aspects, but if it were it would > open up package building on Windows to a larger audience. > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html just to confirm duncan murdochs remark: our Windows machines lack proper development environments (mainly missing perl is the problem for pure R-code packages, I believe?) and we bypass this (for pure R-code packages only, of course) by 1.) install the package on the unix machine into the desired R library 2.) zip the _installed_ package (not the source tree!) found in the R library directory 3.) transfer this archive to the Windows machine 4.) unzip directly into the desired library destination this procedure up to now always worked including properly installed manpages (text + html (and I hope this remains the case in the future...) joerg van den hoff __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] $
Davis, Jacob B. wrote: > If object is user defined is: > > > > object$df.residual > > > > the same thing as > > > > df.residual(object) > > > > This is my first time to encounter the $ sign in R, I'm new. I'm > reviewing "summary.glm" and in most cases it looks as though the $ is > used to extract some characteristic/property of the object, but I'm not > positive. > > > > Thanks > > ~~ > > > > Jacob Davis > > Actuarial Analyst > > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html `$' is for list component extraction. this is really very basic and, for once (I usually don't like this answer...), looking at the very first pages of the very first manual, would be not so bad an idea: http://cran.r-project.org/ --> Manuals --> An Introduction to R __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Function hints
hadley wickham wrote:
>> what I really would love to see would be an improved help.search():
>> on r-devel I found a reference to the /concept tag in .Rd files and the
>> fact that it is rarely used (again: I was not aware of this :-( ...),
>> which might serve as keyword container suitable for improving
>> help.search() results. what about changing the syntax here to
>> something like
>> \concept {
>> keyword = score,
>> keyword = score
>> ...
>> }
>> where score would be restricted to a small range of values (say, 1-3 or
>> 1-5). if package maintainer then would choose a handful of sensible
>> keywords (and scores) for a package and its functions one could expect
>> improved search results. this might be a naive idea, but could a
>> sort-by-relevance in the help.search() output profit from this?
>
> This is not something I think you can solve automatically. Good
> keywording requries a lot of effort, and needs to be consistent to be
> useful. The only way to achieve consistency is to have only person
I was thinking of manual keywording (by the package authors, nobody
else!) as a means to give the search engine (help.search()) reasonable
information including a (subjective) relevance score for each keyword.
of course, the problem is the same as with every (especially permuted)
index: to find the best compromise betweeen indexing next to nothing and
indexing everything (the best probably meaning to index comprehensively
but not excessively with reasonable index terms) in the documents at hand.
sure, consistency could not be enforced but it's not consistent right
now, simply because the real \keyword tag is far to restrictive for
indexing purposes(only a handful of predefined allowed keywords) and
otherwise only the name/alias and title in the Rd files seem to be
searched (and here the author must be really aware that these fields are
at the moment the ones which should be forced to contain the relevant
'keywords' if the function is to be found by help.search -- this imposes
sometimes rather artificial constraints on the wording, especially if
you try to include some general keyword in the title of a very
specialized function).
looking at the example I gave
(help.search("fitting") etc.) it's quite clear that `nls' simply is not
found because 'fitting' does not occur in the title, but I trust, if
asked to provide, say, three keywords, one of them would contain "fit"
or "fitting". I mean, every scientific journal asks you to do just this:
provide some free-text keywords, which you think to be relevant for the
paper. there are no restrictions/directives, usually, but the purpose
(to categorize the paper a bit) is served quite well.
and maybe the \concept tag really is meant for something different, I'm
not sure. what I have in mind really is similar to providing index terms
(plus scores to guide `help.search' in sorting). to stay with the `nls'
example:
\concept {
non-linear fitting = 4
non-linear least-squares = 5
non-linear models = 3
parameter estimimation = 2
gauss-newton = 1
}
would probably achieve that `nls' usually is correctly found (if this
syntax were allowed). apart from the scores (which would be nice, I
think) my main point is that extensive use of \concept (or a new tag
`\index', for instance, if \concept's purpose is actually different --
I'm not sure) should be pushed to get better hits from help.search().
I personally have decided to start using the \concept tag in its present
form for our local .Rd files extensively to "inject" a sufficient number
of free-text relevant keywords into help.search()
joerg
> keywording (difficult/expensive), or exhaustively document the process
> of keywording and then require all package authors to read and use
> (impossible).
>
> Hadley
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Function hints
Jonathan Baron wrote:
> On 06/19/06 13:13, Duncan Murdoch wrote:
>>> `help.search' does not allow full text search in the manpages (I can
>>> imagine why (1000 hits...), but without such a thing google, for
>>> instance, would probably not be half as useful as it is, right?) and
>>> there is no "sorting by relevance" in the `help.search' output, I think.
>>> how this sorting could be achieved is a different question, of course.
>> You probably want RSiteSearch("keyword", restrict="functions") (or even
>> without the "restrict" part).
>
> Yes. The restrict part will speed things up quite a bit, if you
> want to restrict to functions.
>
> Or, alternatively, you could use Namazu (which I use to generate
> what RSiteSearch provides) to generate an index specific to your
> own installed functions and packages. The trick is to cd to the
> directory /usr/lib/R/library, or the equivalent, and then say
>
> mknmz -q */html
>
> which will pick up the html version of all the man pages
> (assuming you have generated them, and I have no idea whether
> this can be done on Windows). To update, say
>
> mknmz --update=. -q */html
>
> Then make a bookmark for the Namazu search page in your browser,
> as a local file. (I haven't given all the details. You have to
> install Namazu and follow the instructions.)
>
> Or, if you have a web server, you could let Google do it for
> you. But, I warn you, Google will fill up your web logs pretty
> fast if you don't exclude it with robots.txt. I don't let it
> search my R stuff.
>
> I think that Macs and various Linux versions also have other
> alternative built-in search capabilities, but I haven't tried
> them. Beagle is the new Linux search tool, but I don't know what
> it does.
>
> Jon
thanks for theses tips. I was not aware of the `RSiteSearch' function
(I did know of the existence of the web sites, though) and this helps,
but of course this is depdendent on web access (off-line labtop
usage...) and does not know of 'local' (non-CRAN) packages (and knows of
maybe "too many" contributed packages, which I might not want to
consider for one reason or the other)
thanks also for the hint on `Namazu'. maybe I do as adviced to get a
index which is aware of my local configuration and private packages.
(under MacOS there is a very good and fast full text search engine, but
it cannot be told to only search the R documentation, for instance, so
one gets lots of other hits as well.)
what I really would love to see would be an improved help.search():
on r-devel I found a reference to the /concept tag in .Rd files and the
fact that it is rarely used (again: I was not aware of this :-( ...),
which might serve as keyword container suitable for improving
help.search() results. what about changing the syntax here to something like
\concept {
keyword = score,
keyword = score
...
}
where score would be restricted to a small range of values (say, 1-3 or
1-5). if package maintainer then would choose a handful of sensible
keywords (and scores) for a package and its functions one could expect
improved search results. this might be a naive idea, but could a
sort-by-relevance in the help.search() output profit from this?
to make it short: I'm not happy with the output, for instance, of
help.search("fitting") #1
vs.
help.search("linear fitting")#2
vs.
help.search("non-linear fitting")#3
I somehow feel that `lm' and `nls' should both be found in the first
search and that they should be near the top of the lists when they are
found.
but `lm' is found only in #1 (near the bottom of the list) and `nls' not
at all (which is really bad). this is partly a problem, of course, of
inconsistent nomenclature in the manpages but also due to the fact that
help.search() only accepts single phrases as pattern (and maybe the
absense of "concept" keywords including a score?)
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Function hints
hadley wickham wrote:
> One of the recurring themes in the recent UserR conference was that
> many people find it difficult to find the functions they need for a
> particular task. Sandy Weisberg suggested a small idea he would like
> to see: a hints function that given an object, lists likely
> operations. I've done my best to implement this function using the
> tools currently available in R, and my code is included at the bottom
> of this email (I hope that I haven't just duplicated something already
> present in R). I think Sandy's idea is genuinely useful, even in the
> limited form provided by my implementation, and I have already
> discovered a few useful functions that I was unaware of.
>
> While developing and testing this function, I ran into a few problems
> which, I think, represent underlying problems with the current
> documentation system. These are typified by the results of running
> hints on a object produced by glm (having class c("glm", "lm")). I
> have outlined (very tersely) some possible solutions. Please note
> that while these solutions are largely technological, the problem is
> at heart sociological: writing documentation is no easier (and perhaps
> much harder) than writing a scientific publication, but the rewards
> are fewer.
>
> Problems:
>
> * Many functions share the same description (eg. head, tail).
> Solution: each rdoc file should only describe one method. Problem:
> Writing rdoc files is tedious, there is a lot of information
> duplicated between the code and the documenation (eg. the usage
> statement) and some functions share a lot of similar information.
> Solution: make it easier to write documentation (eg. documentation
> inline with code), and easier to include certain common descriptions
> in multiple methods (eg. new include command)
>
> * It is difficult to tell which functions are commonly
> used/important. Solution: break down by keywords. Problem: keywords
> are not useful at the moment. Solution: make better list of keywords
> available and encourage people to use it. Problem: people won't
> unless there is a strong incentive, plus good keywording requires
> considerable expertise (especially in bulding up list). This is
> probably insoluable unless one person systematically keywords all of
> the base packages.
>
> * Some functions aren't documented (eg. simulate.lm, formula.glm) -
> typically, these are methods where the documentation is in the
> generic. Solution: these methods should all be aliased to the generic
> (by default?), and R CMD check should be amended to check for this
> situation. You could also argue that this is a deficiency with my
> function, and easily fixed by automatically referring to the generic
> if the specific isn't documented.
>
> * It can't supply suggestions when there isn't an explicit method
> (ie. .default is used), this makes it pretty useless for basic
> vectors. This may not really be a problem, as all possible operations
> are probably too numerous to list.
>
> * Provides full name for function, when best practice is to use
> generic part only when calling function. However, getting precise
> documentation may requires that full name. I do the best I can
> (returning the generic if specific is alias to a documentation file
> with the same method name), but this reflects a deeper problem that
> the name you should use when calling a function may be different to
> the name you use to get documentation.
>
> * Can only display methods from currently loaded packages. This is a
> shortcoming of the methods function, but I suspect it is difficult to
> find S3 methods without loading a package.
>
> Relatively trivial problems:
>
> * Needs wide display to be effective. Could be dealt with by
> breaking description in a sensible manner (there may already by R code
> to do this. Please let me know if you know of any)
>
> * Doesn't currently include S4 methods. Solution: add some more code
> to wrap showMethods
>
> * Personally, I think sentence case is more aesthetically pleasing
> (and more flexible) than title case.
>
>
> Hadley
>
>
> hints <- function(x) {
> db <- eval(utils:::.hsearch_db())
> if (is.null(db)) {
> help.search("abcd!", rebuild=TRUE, agrep=FALSE)
> db <- eval(utils:::.hsearch_db())
> }
>
> base <- db$Base
> alias <- db$Aliases
> key <- db$Keywords
>
> m <- all.methods(class=class(x))
> m_id <- alias[match(m, alias[,1]), 2]
> keywords <- lapply(m_id, function(id) key[key[,2] %in% id, 1])
>
> f.names <- cbind(m, base[match(m_id, base[,3]), 4])
> f.names <- unlist(lapply(1:nrow(f.names), function(i) {
> if (is.na(f.names[i, 2])) return(f.names[i, 1])
> a <- methodsplit(f.names[i, 1])
> b <- methodsplit(f.names[i, 2])
>
> if (a[1] == b[1]) f.names[i, 2] else f.names[i, 1]
> }))
>
>
Re: [R] can I call user-created functions without source() ?
Duncan Murdoch wrote:
> On 6/19/2006 10:19 AM, Joerg van den Hoff wrote:
>> Duncan Murdoch wrote:
>>> Just a few comments below on alternative ways to do the same things:
>>>
>>> On 6/19/2006 8:19 AM, Joerg van den Hoff wrote:
>>>
>>>> for short term usage of some specialized functions I have added some
>>>> lines to the `.Rprofile' in my home(!) directory as follows
>>>> (probably there are smarter solutions, but at least it works):
>>>>
>>>> #source some temporary useful functions:
>>>> fl <- dir(path='~/rfiles/current',patt='.*\\.R$',full.names=TRUE)
>>>> for (i in fl) {cat(paste('source("',i,'")\n',sep="")); source(i)}
>>>> rm(i,fl)
>>>
>>> Another way to do this without worrying about overwriting some
>>> existing variables is
>>>
>>> local({
>>> fl <- ...
>>> for (i in fl) ...
>>> })
>>
>>>
>>> No need to remove fl and i at the end; they were created in a
>>> temporary environment, which was deleted at the end.
>>>
>> sure, that's better (just one more case, where I did'nt know of the
>> existence of a certain function). but what is the difference (with
>> regards to scope) of `i' or `fl' and the functions defined via
>> sourcing? are'nt both objects defined within `local'? why _are_ the
>> functions visible in the workspace? probably I again don't understand
>> the `eval'/environment intricacies.
>
> Whoops, sorry. Yes, you'd need to be careful to assign them into
> globalenv(). Those ...'s weren't so simple.
>
no, no, you _were_ right (I tested it), but I did'nt understand why it
works. I found the answer only now in the `source' manpage: there is an
argument `local = FALSE' which by default enforces assignment into
globalenv(). so the `...' remain unaltered :-). sorry for the confusion.
joerg
PS: and as a follow up on the recent 'setwd vs. cd' thread: as I noted
in passing `source' has a further argument `chdir' to change the
directory prior to sourcing, which is neither the usual name (using
`cd') nor the R-way (using `setwd'), but rather identical to the
corresponding native matlab/octave command (both of these packages have
for quite some time changed to accepting `cd' as well). of course, the
`source' arguments have nothing to do with the commands but consistency
would'nt harm either :-). not important, but a bit funny.
> Duncan Murdoch
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] can I call user-created functions without source() ?
Duncan Murdoch wrote:
> Just a few comments below on alternative ways to do the same things:
>
> On 6/19/2006 8:19 AM, Joerg van den Hoff wrote:
>
>> for short term usage of some specialized functions I have added some
>> lines to the `.Rprofile' in my home(!) directory as follows (probably
>> there are smarter solutions, but at least it works):
>>
>> #source some temporary useful functions:
>> fl <- dir(path='~/rfiles/current',patt='.*\\.R$',full.names=TRUE)
>> for (i in fl) {cat(paste('source("',i,'")\n',sep="")); source(i)}
>> rm(i,fl)
>
> Another way to do this without worrying about overwriting some existing
> variables is
>
> local({
> fl <- ...
> for (i in fl) ...
> })
>
> No need to remove fl and i at the end; they were created in a temporary
> environment, which was deleted at the end.
>
sure, that's better (just one more case, where I did'nt know of the
existence of a certain function). but what is the difference (with
regards to scope) of `i' or `fl' and the functions defined via sourcing?
are'nt both objects defined within `local'? why _are_ the functions
visible in the workspace? probably I again don't understand the
`eval'/environment intricacies.
>>
>> here, I have put all the temporary stuff in a single dedicated dir
>> `~/rfiles/current', but of course you can use several dirs in this
>> way. all files in this dir with names ending in `.R' are sourced on
>> startup of R. this roughly works like one of the directories on
>> MATLAB's search path: every function definition in this directory is
>> 'understood' by R (but everything is loaded into the workspace on
>> startup, no matter, whether you really need it in the end: no real
>> `load on demand').
>
> It's possible to have load on demand in R, and this is used in packages.
> It's probably not worth the trouble to use it unless you're using a
> package.
>
>
> one
>> important difference, though: this is only sensible for function
>> definitions, not scripts ('executable programms' (which would directly
>> be executed on R startup, otherwise).
>> and, contrary to matlab/octave, this is not dynamic: everything is
>> read in at startup, later modifications to the directories are not
>> recognized without explicitely sourcing the files again.
>
> There isn't really any reasonable way around this. I suppose some hook
> could be created to automatically read the file if the time stamp
> changes, but that's not really the R way of doing things: generally in
> R active things are in the workspace, not on disk. A good way to work
> is prepare things on disk, then when they are ready, explicitly import
> them into R.
>
>>
>> if you in addition you want to load definitions from the startup
>> directory where you launch R (your project dir), the above could be
>> modified to:
>>
>> #source some temporary useful functions from startup dir:
>> fl <- dir(path=getwd(),patt='.*\\.R$',full.names=TRUE)
>> for (i in fl) {cat(paste('source("',i,'")\n',sep="")); source(i)}
>> rm(i,fl)
>>
>> in this way you at least don't need a separate `.Rprofile' in each
>> project dir.
>
> Another alternative if you want something special in the project is to
> create a .Rprofile file there, and put source("~/.Rprofile") into it, so
> both the local changes and the general ones get loaded.
>
> Duncan Murdoch
>>
>>
>>
>> joerg
>>
>> __
>> R-help@stat.math.ethz.ch mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide!
>> http://www.R-project.org/posting-guide.html
>
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] can I call user-created functions without source() ?
Duncan Murdoch wrote:
> On 6/19/2006 6:25 AM, (Ted Harding) wrote:
>> On 19-Jun-06 Rob Campbell wrote:
>>> Hi,
>>>
>>> I have to R fairly recently from Matlab, where I have been used to
>>> organising my own custom functions into directories which are adding to
>>> the Matlab search path. This enables me to call them as I would one of
>>> Matlab's built-in functions. I have not found a way of doing the same
>>> thing in R. I have resorted to using source() on symlinks located in
>>> the
>>> current directory. Not very satisfactory.
>>>
>>> I looked at the R homepage and considered making a "package" of my
>>> commonly used functions so that I can call them in one go:
>>> library(myFuncions, lib.loc="/path/to/library") Perhaps this is the
>>> only
>>> solution but the docs on the web make the process seem rather
>>> complicated--I only have a few script files I want to call! Surely
>>> there's a straightforward solution?
>>>
>>> How have other people solved this problem? Perhaps someone has a simple
>>> "package skeleton" into which I can drop my scripts?
>>>
>>>
>>> Thanks,
>>>
>>> Rob
>> There are pros and cons to this, but on the whole I sympathise
>> with you (having pre-R been a heavy matlab/octave user myself).
>>
>> Unfortunately (from this perspective) R does not seem to have
>> an automatic "load-on-demand" facility similar to what happens
>> in matlab (i.e. you call a function by name, and R would search
>> for it in whatever the current search-path is, and load its
>> definition plus what else it depends on).
>>
>> I have a few definitions which I want in every R session, so
>> I have put these in my ".Rprofile". But this is loaded from
>> my home directory, not from the directory where I was when I
>> started R, so it is the same every time.
>
> Which version of R are you using? This is not the current documented
> behaviour. It looks in the current directory first, and only in your
> home directory if that fails.
>
> Duncan Murdoch
>
>
> Again, one of the
>> conveniences of the matlab/octave approach is that you can
>> have a different sub-directory for each project, so if you
>> start work in a particular one then you have access to any
>> special definitions for that project, and not to others.
>>
>> I'm no expert on this aspect of R, but I suspect that the way
>> start-up is organised in R does not fit well with the other
>> kind of approach. I stand to be corrected, of course ...
>>
>> And others may well have formulated their own neat work-rounds,
>> so we wait eagerly to hear about these!
>>
>> Best wishes,
>> Ted.
>>
using `package.skeleton' (as already mentioned) for the package dir
layout and `prompt' for generating templates of the needed manpages is
not so bad (once you get used to it), if the things you have in mind are
at least of some long(er) term value (for you): at least you are forced
(sort of) to document your software...
for short term usage of some specialized functions I have added some
lines to the `.Rprofile' in my home(!) directory as follows (probably
there are smarter solutions, but at least it works):
#source some temporary useful functions:
fl <- dir(path='~/rfiles/current',patt='.*\\.R$',full.names=TRUE)
for (i in fl) {cat(paste('source("',i,'")\n',sep="")); source(i)}
rm(i,fl)
here, I have put all the temporary stuff in a single dedicated dir
`~/rfiles/current', but of course you can use several dirs in this way.
all files in this dir with names ending in `.R' are sourced on startup
of R. this roughly works like one of the directories on MATLAB's search
path: every function definition in this directory is 'understood' by R
(but everything is loaded into the workspace on startup, no matter,
whether you really need it in the end: no real `load on demand'). one
important difference, though: this is only sensible for function
definitions, not scripts ('executable programms' (which would directly
be executed on R startup, otherwise).
and, contrary to matlab/octave, this is not dynamic: everything is read
in at startup, later modifications to the directories are not recognized
without explicitely sourcing the files again.
if you in addition you want to load definitions from the startup
directory where you launch R (your project dir), the above could be
modified to:
#source some temporary useful functions from startup dir:
fl <- dir(path=getwd(),patt='.*\\.R$',full.names=TRUE)
for (i in fl) {cat(paste('source("',i,'")\n',sep="")); source(i)}
rm(i,fl)
in this way you at least don't need a separate `.Rprofile' in each
project dir.
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Access and assign list sub-elements using a string suchas "l$a$b"
Petr Pikal wrote: > Hi > > yes you are correct, I remembered there is something with eval from > older posts but did not find a connection to parse from eval help > page. Shouldn't there be a link? Or even an example? would be a good thing to do (there only is a link from parse to eval). after all eval(parse(text = some_string)) is what many users (which might come from matlab/octave where this _would_ suffice) want when they try out eval(some_string)). and a standard eval(parse(text=...) example on the eval manpage (where most people probably would look, would be very good. > > Best regards > Petr > > > On 15 Jun 2006 at 17:21, Dimitris Rizopoulos wrote: > > From: "Dimitris Rizopoulos" <[EMAIL PROTECTED]> > To: "Petr Pikal" <[EMAIL PROTECTED]>, <[EMAIL PROTECTED]> > Copies to: > Subject: Re: [R] Access and assign list sub-elements using a > string suchas "l$a$b" > Date sent:Thu, 15 Jun 2006 17:21:26 +0200 > >> - Original Message - >> From: "Petr Pikal" <[EMAIL PROTECTED]> >> To: "Gregory Jefferis" <[EMAIL PROTECTED]> >> Cc: >> Sent: Thursday, June 15, 2006 4:56 PM >> Subject: Re: [R] Access and assign list sub-elements using a string >> suchas "l$a$b" >> >> >>> Hi >>> very, very close >>> >>> >>> On 15 Jun 2006 at 13:27, Gregory Jefferis wrote: >>> >>> Date sent: Thu, 15 Jun 2006 13:27:05 +0100 >>> From: Gregory Jefferis <[EMAIL PROTECTED]> >>> To: "[EMAIL PROTECTED]" >>> <[EMAIL PROTECTED]> >>> Forwarded to: >>> Forwarded by: Gregory Jefferis <[EMAIL PROTECTED]> >>> Date forwarded: Thu, 15 Jun 2006 14:54:13 +0100 >>> Subject:[R] Access and assign list sub-elements using a >>> string such as "l$a$b" >>> If I have a list I can set a sub-element as follows on the command line: people=list() people$tom$hair="brown" people But what if I have a string containing the name of the sub-element that I want to access? subel= "people$tom$hair" get(subel) # returns error assign(subel,"red") # silent but doesn't change list people >>> See what happens when >>> >>> people<-assign(subel, "red") >> but I think this is not what Greg wanted; the above just assigns "red" >> to object 'people' (i.e., check `str(assign(subel, "red"))'). If I >> understood correctly, the following could be of help: >> >> people <- list() >> people$tom$hair <- "brown" >> people >> # >> subel <- "people$tom$hair" >> eval(parse(text = subel)) >> eval(parse(text = paste(subel, "<- 'red'"))) >> people >> >> >> Best, >> Dimitris >> >> >> >> Dimitris Rizopoulos >> Ph.D. Student >> Biostatistical Centre >> School of Public Health >> Catholic University of Leuven >> >> Address: Kapucijnenvoer 35, Leuven, Belgium >> Tel: +32/(0)16/336899 >> Fax: +32/(0)16/337015 >> Web: http://med.kuleuven.be/biostat/ >> http://www.student.kuleuven.be/~m0390867/dimitris.htm >> >> >>> HTH >>> Petr >>> >>> The attempts above using assign/get won't do what I am trying to do [nor according to the help should they]. I would be very grateful for any suggestions. Many thanks, Greg. -- Gregory Jefferis, PhD and: Research Fellow Department of Zoology St John's College University of Cambridge Cambridge Downing Street CB2 1TP Cambridge, CB2 3EJ United Kingdom Lab Tel: +44 (0)1223 336683 Office: +44 (0)1223 339899 Lab Fax: +44 (0)1223 336676 http://www.zoo.cam.ac.uk/zoostaff/jefferis.html [EMAIL PROTECTED] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html >>> Petr Pikal >>> [EMAIL PROTECTED] >>> >>> __ >>> R-help@stat.math.ethz.ch mailing list >>> https://stat.ethz.ch/mailman/listinfo/r-help >>> PLEASE do read the posting guide! >>> http://www.R-project.org/posting-guide.html >>> >> >> Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm >> > > Petr Pikal > [EMAIL PROTECTED] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] plot two graphs with different length of data
Eric Hu wrote: > Hi I am trying to plot two data set in the same picture window without > overlapping with each other. I am using the format plot(x1,y1,x2,y2) > but get the following error message: > >> plot(as.numeric(r0[,4]),as.numeric(r0[,7]),as.numeric(r0[,4]),as.numeric(r0[,7][ind[,1]])) > Error in plot.window(xlim, ylim, log, asp, ...) : > invalid 'ylim' value > > Can anyone tell me what went wrong? Thanks. > > Eric > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html matplot(cbind(x1,x2), cbind(y1,y2)) might be what you want. (if x1, x2 and y1,y2 are of equal length. otherwise pad the short ones with NAs or use `matplot' with type =n (to get the scaling of the plot right), followed by `plot(x1,y1), lines(x2,y2)') __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] solving first-order differential equation
ZhanWu Dai wrote:
> I am an initial user of R. Could you give me some explanations or examples on
> how to solve the first order differential equations by the first-order
> Runge-Kutta method?
>
> Thank you very much
>
> Kind regards
>
>
>
>
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
not really an answer, but a remark:
if your ODE is of the form
dy
--- - k y = A f(x)
dx
(k, A const.) it might be a better idea to use the 'analytic' solution
instead of runge-kutta (faster, probably more accurate).
for instance, if the initial condition is
y(x=0) = 0
and you're looking only at x>0 the solution simply is
y(x) = A (x) {*} exp(-kx)
where {*} means the finite (continous) convolution extending from 0 to x:
y(x) = A integral from z=0 to z=x {f(z) exp(-k(x-z)) dz}
(which, of course, still has to be computed numerically in general.)
this closed-form solution can then
be used, for instance, to determine the unknown parameters (k, A) from a
least squares fit to measured f(x), y(x)
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Building packages in R - 'private' functions
Dan Rabosky wrote:
> Hello.
>
> I am creating an R package that I'd like to submit to CRAN (OS Windows
> XP). How do I distinguish among 'public' functions, e.g., those that are
> intended to be called by users of the package and for which I am providing
> documentation & examples, and 'private' functions, which are used
> internally by the 'public' functions, but for which I do not wish to
> provide documentation? The private functions are all coded in R (nothing
> in C or Fortran) and are essential to the operation of several public
> functions.
>
> I have been unable to find any documentation on this in the 'writing r
> extensions' manual', on previous posts to R-help, or through any other
> source. One possibility is to include the source code for the 'private'
> functions within the public functions. However, since multiple public
> functions utilize the same core set of 'private' functions, this seems
> unwieldy and redundant at best.
>
> If I simply include the source for the 'private' functions in the "R"
> directory (without corresponding *.Rd and *.html documentation in /man),
> then check the package with "R CMD check', it does appear to process the
> private functions (and successfully builds with R CMD build). However, I
> do receive a warning for including undocumented code objects. Is this the
> recommended approach and/or is there a better way to do this? One
> potential problem with this approach is that - should an error occur within
> a private function, it may be very difficult for the user to decipher the
> nature of the problem.
>
> Any suggestions will be greatly appreciated.
> ~Dan Rabosky
>
>
>
> Dan Rabosky
> Department of Ecology and Evolutionary Biology
> 237 Corson Hall
> Cornell University
> Ithaca, NY14853-2701 USA
> [EMAIL PROTECTED]
> web: http://www.birds.cornell.edu/evb/Graduates_Dan.htm
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
it's in the 'extensions' manual, including an example I believe:
1.
create a file `NAMESPACE' in the package top level dir (beside `R' and
`man') containing the single line
exportPattern("^[^\\.]")
2.
name all private functions with a leading `.' (more precisely: all
functions starting with a `.' are private in this setting).
of course, you can modify the pattern to suit another naming convention.
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] nls & fitting
Lorenzo Isella wrote: > Dear All, > I may look ridiculous, but I am puzzled at the behavior of the nls with > a fitting I am currently dealing with. > My data are: > >x N > 1 346.4102 145.428256 > 2 447.2136 169.530634 > 3 570.0877 144.081627 > 4 721.1103 106.363316 > 5 894.4272 130.390552 > 6 1264.9111 36.727069 > 7 1788.8544 52.848587 > 8 2449.4897 25.128742 > 9 3464.1016 7.531766 > 10 4472.1360 8.827367 > 11 6123.7244 6.600603 > 12 8660.2540 4.083339 > > I would like to fit N as a function of x according to a function > depending on 9 parameters (A1,A2,A3,mu1,mu2,mu3,myvar1,myvar2,myvar3), > namely > N ~ > (log(10)*A1/sqrt(2*pi)/log(myvar1)*exp(-((log(x/mu1))^2)/2/log(myvar1)/log(myvar1)) > > +log(10)*A2/sqrt(2*pi)/log(myvar2)*exp(-((log(x/mu2))^2)/2/log(myvar2)/log(myvar2)) > > +log(10)*A3/sqrt(2*pi)/log(myvar3)*exp(-((log(x/mu3))^2)/2/log(myvar3)/log(myvar3))) > > (i.e. N is to be seen as a sum of three "bells" whose parameters I need > to determine). > > > So I tried: > out<-nls(N ~ > (log(10)*A1/sqrt(2*pi)/log(myvar1)*exp(-((log(x/mu1))^2)/2/log(myvar1)/log(myvar1)) > > +log(10)*A2/sqrt(2*pi)/log(myvar2)*exp(-((log(x/mu2))^2)/2/log(myvar2)/log(myvar2)) > > +log(10)*A3/sqrt(2*pi)/log(myvar3)*exp(-((log(x/mu3))^2)/2/log(myvar3)/log(myvar3))) > ,start=list(A1 = 85, > A2=23,A3=4,mu1=430,mu2=1670,mu3=4900,myvar1=1.59,myvar2=1.5,myvar3=1.5 ) > ,algorithm = "port" > ,control=list(maxiter=2,tol=1) > ,lower=c(A1=0.1,A2=0.1,A3=0.1,mu1=0.1,mu2=0.1,mu3=0.1,myvar1=0.1,myvar2=0.1,myvar3=0.1) > ) > > getting the error message: > Error in nls(N ~ (log(10) * A1/sqrt(2 * pi)/log(myvar1) * > exp(-((log(x/mu1))^2)/2/log(myvar1)/log(myvar1)) + : > Convergence failure: singular convergence (7) > > > I tried to adjust tol & maxiter, but unsuccessfully. > If I try fitting N with only two "bells", then nls works: > > out<-nls(N ~ > (log(10)*A1/sqrt(2*pi)/log(myvar1)*exp(-((log(x/mu1))^2)/2/log(myvar1)/log(myvar1)) > > +log(10)*A2/sqrt(2*pi)/log(myvar2)*exp(-((log(x/mu2))^2)/2/log(myvar2)/log(myvar2)) > ) > ,start=list(A1 = 85, A2=23,mu1=430,mu2=1670,myvar1=1.59,myvar2=1.5 ) > ,algorithm = "port" > ,control=list(maxiter=2,tol=1) > ,lower=c(A1=0.1,A2=0.1,mu1=0.1,mu2=0.1,myvar1=0.1,myvar2=0.1) > ) > > out > Nonlinear regression model > model: N ~ (log(10) * A1/sqrt(2 * pi)/log(myvar1) * > exp(-((log(x/mu1))^2)/2/log(myvar1)/log(myvar1)) + log(10) * > A2/sqrt(2 * pi)/log(myvar2) * > exp(-((log(x/mu2))^2)/2/log(myvar2)/log(myvar2))) >data: parent.frame() > A1 A2mu1mu2 myvar1 myvar2 > 84.920085 40.889968 409.656404 933.081936 1.811560 2.389215 > residual sum-of-squares: 2394.876 > > Any idea about how to get nls working with the whole model? > I had better luck with the nls.lm package, but it does not allow to > introduce any constrain on my fitting parameters. > I was also suggested to try other packages like optim to do the same > fitting, but I am a bit unsure about how to set up the problem. > Any suggestions? BTW, I am working with R Version 2.2.1 > > Lorenzo > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html apart from the fact that fitting 9 parameters to 12 points quite genereally will not yield satisfactory results (at least estimates will have huge uncertainties), total failure in your case seems obviouus if you plot your data: it's not even obvious where the three peaks (means of the gaussians) should be: all below x=2000 or is there one peak at about x=4500 and one of the 'peaks' below x=2000 is spurious? if you can't decide, nls has problems. moreover: how should reliable estimates of the standard deviations (width) of the gaussian result if the peaks essentially consist of exactly one point? in short: I believe, you either need much more data points or you must put in substantial a priori information (e.g. either means or standard deviations of the gaussians). __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] problem with pdf() R 2.2.1, os 10.4.6
Marc Schwartz wrote: > On Fri, 2006-05-19 at 16:36 -0700, Betty Gilbert wrote: >> Hi, >> I'm trying to write a histogram to a pdf >> >> pdf() >> plot<-hist(c, xlim=c( 0.69, 0.84), ylim=c(0,100)) >> >> when I try to open the pdf I can't open it, there is always some >> error . Is there something I should add to make it run under this >> operation system? I had problems upgrading to 2.3 (problem >> downloading packages) so I'm not sure an upgrade will work out with >> me. I just want a publication quality histogram... >> >> thank you, >> betty > > Betty, > > You need to explicitly close the pdf device with: > > dev.off() > > once the plotting related code has completed. See the example in ?pdf > for more information. > > Otherwise, the result of hist() is not flushed to the disk file and the > file then properly closed. > > HTH, > > Marc Schwartz > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html and maybe you are better off with plot_something dev.print(pdf) with the standard settings of pdf() you get more or less a pdf file maintaining the current aspect ratio of the plot on the screen prior to the dev.print call __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] help
karim99.karim wrote: > Dear Sir, > > I’am a frensh student and i’am a new user of the R software. > > After using the command (x<-read.delim(“clipboard”) to read a spreadsheet of > Excel, I want to run the bds test and calculate the Lyapunov exponent. I have > charged the R software by the packages tseries and tseriesChaos. when i run > bds.test(x,m=2) Unfortunately the R software displays “error in > as.vector(x,mode= “double”) : the object (list) can not be automatically > converted on double” so what shall I do to run this two tests(lyapunov and > bds)? And what is my mistake? > > I thank you in advance and I’am waiting forward your e-mail > > This si my e-mail: [EMAIL PROTECTED] > > Accédez au courrier électronique de La Poste : www.laposte.net ; > 3615 LAPOSTENET (0,34 €/mn) ; tél : 08 92 68 13 50 (0,34€/mn) > > > > [[alternative HTML version deleted]] > > > > > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html seems that the data import failed. probably the exported excel data are not in a format accepted correctly by read.delim. check this. alternative (maybe): there is an additional package 'gdata' on CRAN which you can download and install. it contains a function read.xls which can read directly the excel binary format (with some limitations). __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Can't there be a cd command?
Duncan Murdoch wrote:
> On 5/16/2006 5:46 AM, Joerg van den Hoff wrote:
>> Manuel López-Ibáñez wrote:
>>> Jan T. Kim wrote:
>>>> That's an idea I like very much too -- much better than the currently
>>>> popular idea of "protecting" users from the "unfriendliness" of
>>>> programming, anyway...
>>>>
>>>
>>> It is just my opinion that the amount of mail in R-help speaks
>>> volumes about the current "friendliness" [1], or lack thereof, of R.
>>> Perhaps I am just the only one who thinks this way...
>>>
>>> [1] http://en.wikipedia.org/wiki/Usability
>>>
>>>
>>
>> I think you are 100% right: the r-help list says it all. needless to
>> say, R is a great achievment without any doubt, but claiming that it's
>> easy to use (beyond the most basic arithmetics) is really wishful
>> thinking.
>
> This is sloppy thinking. The volume of mail here shows that there are a
> lot of questions, perhaps because there are a lot of users.
well, as far as my english goes, 'sloppy' is a strong word (and apart
from mathematicians physicists (my background) probably are the people
who are most allergic to being accused of it :-)) and it's an overhasty
conclusion on your side, I'd say.
I want to make clear beforehand, that I do _not_ think this a very
important discussion, but rather an informal exchange of opinions, so
maybe this takes us all a bit to far, but anyway:
for one, I myself (and I think manuel, too) was not talking of the shear
volume of mails (this obviously would have to be 'calibrated' against
the total number of R users and the resulting quantity had to be
compared to other help-lists). at least my impression is, that there are
a number of reoccuring difficulties in the mail, which are rather
specific to R's design (whether this situation could or should be
altered: that would be a different topic). certainly, among these are
the subsetting/indexing issues, certainly lazy evaluation, certainly
anything related to environments, namespaces, computing on the language
(substitute, eval, ...).
> You're also misquoting Jan: he didn't say R was "easy to use", he said
> that the idea of urging people to program is better than trying to be
> too friendly and protecting them from it.
I did'nt quote anyone AFAIKS, so I can't have misquoted anyone (having
misinterpreted someone I cannot rule out). the 'easy to use' did not
refer to a special statement from this thread, but to the general
impression one can get from the list and contributed as well as official
manuals (I only checked now: the 'what is R?' section on the homepage
contains one 'ease', one 'easy', one 'easily' within a total of two or
three paragraphs...).
it is pointless to dwell on this to long: what is easy for you might be
difficult for me or vice versa, depending on the question to answer/
problem to solve. _if_ I take the freedom to interpret it as 'easy for
the pedestrians', the statements are simply not true ("easily extended
via packages"??).
with reference to the idea of urging people to programm: well, the idea
in itself is not objectionable, the question is how realistic the
approach is (i.e. how large will be the success rate of getting people
to programm, which otherwise would'nt have done _and_ is this rate
larger in R than in the other packages?).
>
>> I don't think programming R is easier than programming C, for example.
>
> I do both, and I think R programming is easier. It has a more sensible
> idea of scoping, it doesn't have the preprocessor doing bizarre
> transformations to the text, it doesn't have pointers writing to random
> memory locations, it can handle strings much more sensibly.
this is all very well, though I only partly agree, but this a very
technical assessment anyway and seems to indicate that a non-programmer
will not be much better off with R as a 'starting language' than with C
(since your criteria mostly will not be initially 'operational'). I
_bet_ this starting phase would be easier with MATLAB/octave (but I'm
not arguing for pushing beginners to MATLAB!).
>
> On the negative side, the vector orientation of R encourages people to
> come up with clever APL-style one-liners that are unreadable; the lack
> of type declarations, the weird handling of indexing, the strange object
> oriented programming models all make R programming hard.
yepp. and cascaded apply/lapply calls, I'd add.
>
> So R is not easy, but it's easier than C.
by a margin, maybe, though I have people in my group who definitely do
object (making especially a point o
Re: [R] Can't there be a cd command?
Manuel López-Ibáñez wrote:
> Jan T. Kim wrote:
>> That's an idea I like very much too -- much better than the currently
>> popular idea of "protecting" users from the "unfriendliness" of
>> programming, anyway...
>>
>
> It is just my opinion that the amount of mail in R-help speaks volumes
> about the current "friendliness" [1], or lack thereof, of R. Perhaps I
> am just the only one who thinks this way...
>
> [1] http://en.wikipedia.org/wiki/Usability
>
>
I think you are 100% right: the r-help list says it all. needless to
say, R is a great achievment without any doubt, but claiming that it's
easy to use (beyond the most basic arithmetics) is really wishful thinking.
I don't think programming R is easier than programming C, for example.
This is not to say that it takes the same time to solve the same problem
in both languages, since in R many, many things are already there
(either in the language (vectorized computations) or in the packages).
but the quantity 'number of new lines of working code per hour' should
be about the same.
I have used MATLAB/octave previously. in comparison to R, the MATLAB
language sure is not pretty, not consistent and less powerful, but
without any question easier. and of course, the interactive use is
facilitated by the absence of lots of paranthesis and quotes and the
way, unknown commands are resolved by looking them up in the search path.
but that's not a situation, one should complain about: R simply cannot
be everybody's darling.
what I always find a bit strange is the explicit claim/aim to blur the
distinction between users and programmers: one could postulate the same
aim (or property) for MATLAB, octave, IDL, or even any other language
providing an interactive interpreter (python, ruby, Tcl, C++/Cint/root,
...). in my experience the distinction between users (rather:
non-programmers) and programmers always is quite sharp.
again, in comparison to MATLAB (where I have some experience) I have
noted that the 'users' (a.k.a. non-programmers) have more difficulties
in using R interactively, mainly for the following reasons:
- difficulties in understanding the neccessety for all those paranthesis
(why can't I just say 'ls'?)
- subsetting ((x), [x], [[x]], $x, $"x", ['x'], ["x"], [["x"]], ... or
what!?)
- R's strong functional orientation: R encourages writing functions with
_lots_ of arguments which you must understand/know of (even if you
don't touch the defaults) and the user must construct the suitable
function call and launch it to get his/her results.
of course one can write interactive ('please enter ...') programms, but
usually does'nt since it 's much more convenient from the programmers
point of view to utilize the functional approach). in a way R's
functions behave like UNIX commands including lots of command line
parameters: they are great for experienced
users but the average guy neither understands
tar[EMAIL PROTECTED][block] [tarfile]
[exclude-file] {-I include-file | -C directory | file |
file} ...
nor
plot(x, y = NULL, type = "p", xlim = NULL, ylim = NULL,
log = "", main = NULL, sub = NULL, xlab = NULL, ylab = NULL,
ann = par("ann"), axes = TRUE, frame.plot = axes,
panel.first = NULL, panel.last = NULL,
col = par("col"), bg = NA, pch = par("pch"),
cex = 1, lty = par("lty"),
lwd = par("lwd"), asp = NA, ...)
easily.
so, to make a long story short: it would'nt harm avoiding the claim that
R is 'easy'. it's probably more of the "you either love it or hate it"
variety (so the unicycle metaphor occuring in this thread fits very
well, I believe).
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Can't there be a cd command?
Issac Trotts wrote:
> On 5/9/06, Prof Brian Ripley <[EMAIL PROTECTED]> wrote:
>> On Wed, 10 May 2006, Peter Dalgaard wrote:
>>
>>> "Issac Trotts" <[EMAIL PROTECTED]> writes:
>>>
I'm on Linux, so it didn't occur to me to look in the Windows FAQ
until it came up in the Google search. Why should Windows users be
the only ones who learn how to change directories?
>>> Well, on Linux et al., the normal mode of operation is to start R from
>>> the commmand line in the directory where you want to be. On Windows
>>> you start from a desktop icon or a menu entry, which implies that
>>> you're always in the same directory, and usually the wrong one.
>> Exactly. I don't think I have ever used setwd() on Linux.
>>
>> Also, I have never seen this asked before for a Unix-alike, so it seems
>> not to be a >F> problem to think everyone does too, and it isn't necessarily so.
>>
>> Frequently asked questions do make it to the FAQs: it is a defence
>> mechanism for the volunteers supporting R.
>>
>> help.search("directory") gets you there.
>
> OK, good to know. Thanks for helping.
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
probably this has used up sufficient disk space in the mail archive
already, but anyway a bit of support to the initial request: one of the
very first things I did, when I started using R, was adding
cd <- setwd
pwd <- getwd
to my .Rprofile. and when I encouraged colleagues to try R out, the 'top
three' questions where
"how can I get out? I tried ^C, exit(), quit(), but nothing works"
"how can I change the directory? `cd' does'nt work"
"where am I? `pwd' does'nt work" (`pwd' is something like a 'local
speciality', I presume...)
this is of course a marginal obstacle for new users, but avoidable all
the same. (and maybe most people (like myself) in the first place don't
think it important enough to post it -- maybe it _would_ be in the FAQ
otherwise...)
so, _if_ 'cd' would be recognized in future releases, it would'nt do any
harm, would it?
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] boxplot - labelling
Sachin J wrote: > Hi, > > How can I get the values of mean and median (not only points but values > too) on the boxplot. I am using boxplot function from graphics package. > Following is my data set > > > df > [1] 5 1 1 0 0 10 38 47 2 5 0 28 5 8 81 21 12 9 1 12 2 4 22 3 > > > mean.val <- sapply(df,mean) >> boxplot(df,las = 1,col = "light blue") >> points(seq(df), mean.val, pch = 19) > > I could get mean as dot symbol but i need values too? Also how to print the > x-axis labels vertically instead of horizontally? Is there any other function > to achieve these? > > Thanx in advance. > > Sachin > > > > > __ > > > > [[alternative HTML version deleted]] > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html use par(las = 3) bxp <- boxplot(df,las = 1,col = "light blue") to control label orientation and access the data displayed in the boxplot (see `?par', `?boxplot') use text(x, y, 'some_text_or_number') to add text to the boxplot. for x, use the 'number of the box' (i.e. `1' in the example above). for y use the mean or median or whatever you wanted). joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] symbols plot
before getting scolded for submitting a (non-)bug report: when using the 'symbols' function for plotting boxplot data (i.e. using 'boxplots' symbols), I noted that the x/y-position of the symbols is associated with the center of the box. while this is obviously natural for a usual plotting symbol (say a circle or a rectangle), it is probably not desired if one uses the 'boxplots' symbols: looking at such a plot, probably everyone will assume that y-position is that of the median within the box (in other words: that one can read off the values of the medians from the y-axis). the current behaviour is counter-intuitive, I believe, if the distributions are asymmetrical (and the median is not centered in it's box). (I even think, that such plots are misinterpreted easily: think what happens if the median lies very near one of the hinges in one box and is centered in another one within the same plot and the medians are actually the same) in short: I think the 'boxplots' should not be centered at the specified y-coordinates but rather drawn with a y-coordinate of y + bxh * (0.5 - bxm) where bxh and bxm are the second and fifths column of the 'boxplots' matrix. in this way, the median position is identical to the specified y-coordinate. at the very least, I think, the manpage should explicitely state that all symbols (including boxplots) are positioned with their geometrical center at the specified coordinates. right or wrong? regards, joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] bug in map('world') ?
hi,
did'nt see anything in the archive:
map('world',pro='rectangular',para=0)
yields a strange artifact (horizontal bar) extending over the whole map
at a certain latitude range (approx 65 deg. north), whereas
map('world',pro='rectangular',para=180)
(which should be the same) does not show the artifact.
the artifact shows up in other projections as well, e.g.
map('world',pro='bonne',para=45)
which seems to show that the problem might be in the data base of the
map (i.e. polygon definition)??
any ideas?
regards,
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] inconsistent behaviour of ifelse and if ... else
is the behaviour
val <- ifelse(TRUE, numeric(0), 123)
val #NA
intended or is it a bug, i.e. should an empty object be returned as
might be expected(also in comparsion to what an explicit
val <- {if(TRUE) numeric(0) else 123} yields)?
thanks,
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] nls and weighting
I posted this a week ago on r-devel but to no avail and hope this not considered cross-posting: ===cut=== hi everybody, which each release I hope that the section "weights: an optional numeric vector of (fixed) weights. When present, the objective function is weighted least squares. _not yet implemented_" in the help page of 'nls' is missing the last sentence. are their any plans to allow/include weighting in the upcoming releases? modifying the cost function to include the weights is probably not the problem, I presume. what is the reason for not including the weighting? are they related to the 'statistical' output (estimation of parameter uncertainties and significances?). I know of the "y ~ M" vs. "~ sqrt(W)*(y-M)" work around suggestion in MASS to include weighting. (I understand that resulting error estimates und confidence intervals from 'nls' are wrong in this case. right?) would'nt it be sensible to inlcude weighting in 'nls' at least on this level, i.e. weighted parameters provided now, correct error estimates and the like coming later? regards, joerg ===cut=== I post this here again hoping to learn about possible work arounds (apart from the MASS proposal) for the current situation with 'nls' (no weighting allowed). thanks in advance joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] problem with lapply(x, subset, ...) and variable select argument
Gabor Grothendieck wrote:
> The problem is that subset looks into its parent frame but in this
> case the parent frame is not the environment in tt but the environment
> in lapply since tt does not call subset directly but rather lapply does.
>
> Try this which is similar except we have added the line beginning
> with environment before the print statement.
>
> tt <- function (n) {
>x <- list(data.frame(a=1,b=2), data.frame(a=3,b=4))
>environment(lapply) <- environment()
>print(lapply(x, subset, select = n))
> }
>
> n <- "b"
> tt("a")
>
> What this does is create a new version of lapply whose
> parent is the environment in tt.
>
>
> On 10/10/05, joerg van den hoff <[EMAIL PROTECTED]> wrote:
>
>>I need to extract identically named columns from several data frames in
>>a list. the column name is a variable (i.e. not known in advance). the
>>whole thing occurs within a function body. I'd like to use lapply with a
>>variable 'select' argument.
>>
>>
>>example:
>>
>>tt <- function (n) {
>> x <- list(data.frame(a=1,b=2), data.frame(a=3,b=4))
>> for (xx in x) print(subset(xx, select = n)) ### works
>> print (lapply(x, subset, select = a)) ### works
>> print (lapply(x, subset, select = "a")) ### works
>> print (lapply(x, subset, select = n)) ### does not work as intended
>>}
>>n = "b"
>>tt("a") #works (but selects not the intended column)
>>rm(n)
>>tt("a") #no longer works in the lapply call including variable 'n'
>>
>>
>>question: how can I enforce evaluation of the variable n such that
>>the lapply call works? I suspect it has something to do with eval and
>>specifying the correct evaluation frame, but how?
>>
>>
>>many thanks
>>
>>joerg
>>
>>__
>>R-help@stat.math.ethz.ch mailing list
>>https://stat.ethz.ch/mailman/listinfo/r-help
>>PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
>>
>
>
many thanks to thomas and gabor for their help. both solutions solve my
problem perfectly.
but just as an attempt to improve my understanding of the inner workings
of R (similar problems are sure to come up ...) two more question:
1.
why does the call of the "[" function (thomas' solution) behave
different from "subset" in that the look up of the variable "n" works
without providing lapply with the current environment (which is nice)?
2.
using 'subset' in this context becomes more cumbersome, if sapply is
used. it seems that than I need
...
environment(sapply) <- environment(lapply) <- environment()
sapply(x, subset, select = n))
...
to get it working (and that means you must know, that sapply uses
lapply). or can I somehow avoid the additional explicit definition of
the lapply-environment?
again: many thanks
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] problem with lapply(x, subset, ...) and variable select argument
I need to extract identically named columns from several data frames in
a list. the column name is a variable (i.e. not known in advance). the
whole thing occurs within a function body. I'd like to use lapply with a
variable 'select' argument.
example:
tt <- function (n) {
x <- list(data.frame(a=1,b=2), data.frame(a=3,b=4))
for (xx in x) print(subset(xx, select = n)) ### works
print (lapply(x, subset, select = a)) ### works
print (lapply(x, subset, select = "a")) ### works
print (lapply(x, subset, select = n)) ### does not work as intended
}
n = "b"
tt("a") #works (but selects not the intended column)
rm(n)
tt("a") #no longer works in the lapply call including variable 'n'
question: how can I enforce evaluation of the variable n such that
the lapply call works? I suspect it has something to do with eval and
specifying the correct evaluation frame, but how?
many thanks
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R graphics
Sam Baxter wrote: > Hi > > I am trying to set up 16 graphs on one graphics page in R. I have used > the mfrow=c(4,4) command. However I get a lot of white space between > each graph. Does anyone know how I can reduce this? > > Thanks > > Sam > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html > ?layout as an alternative to par(mfrow) might be helpful anyway too large margins: ?par reduce value of "mar", for instance __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Michaelis-menten equation
I believe the following is correct: 1. first of all, as peter daalgaard already pointed out, your data Cp(t) are following a straight line very closely, i.e. 0.-order kinetics 2. for your diff. eq. this means that you are permanently in the range cp >> Km so that dCp/dt = - Vm/Vd = const. =: -b and, therefore, Cp = a - b*t 3. you can't get any reliable information concerning Km from the fit. the solution of the diff. eq (according to Maxima), namely t + const. = -(Km*log(Cp) + Cp)/(Vm/Vd) tells you the same: Km*log(cp) << Cp in your data. 4. in any case (even in the range Km ~= Cp) you can't determine Vm _and_ Vd separately according to your diff. eq., you only get the ratio b = Vm/Vd. this does make sense: what you are measuring is the decreasing plasma concentration, you don't have any information concerning the "relevant" volume fraction, i.e. the Vd, in you data. therefore any variation in the effective max. velocity can either be ascribed to a variation of Vm or to a modified Vd. in other words: you should view Vm* = Vm/Vd as your "effective" Vm. 5. x<-PKindex[,1] y<-PKindex[,2] res <- lm ( y~x ) yields a=8.561, b=Vm*= 0.279 summary(abs(residuals(r)/y)*100) Min. 1st Qu. MedianMean 3rd Qu.Max. 0.00863 0.09340 0.13700 0.26500 0.22900 1.37000 i.e. 0.265 percent deviation on average between data and fit. I believe that is the max. information you can get from your data ((the result for "b" is accidently is not so far away from the ratio of your results Vm = 10.04, Vd = 34.99 which actually must be completely unstable (double both parameters and nothing happens to the fit). 6. you state you simulated the data with km=4.69? using the above Vm and Vd, the resulting data are (not unexpectedly) quite different from those you used as input to the fit. maybe you made an error somewhere "upstream"? 7. in conclusion: don't try to fit Vm _and_ Vd separately, check whether your simulated data are correct, keep in mind that if km< Hi, > >We are doing a pharmaockinetic modeling. This model is > described as intravenous injection of a certain drug into > the body. Then the drug molecule will be eliminated (or decayed) > from the body. We here used a MM eq. to describe the elimination > process and the changes of the drug conc.. So the diff. eq. would > be: dCp/dt = -Vm/Vd * Cp/(Km+Cp). Vd is a volume of distribution. > We used lsoda to solve the diff. eq. first and fit the diff. eq. > with optim first (Nelder-Mead simplex) and followed by using nls > to take over the fitting process of optim. However, we can not > obtain the correct value for Km if we used the above model. The > correct Km can be obtained only when we modeled the diff eq. with > dCp/dt= -Vm/Vd * Cp/(Km/vd + Cp). Now we lost. The data were > from simulation with known Vm and Km. Any idea? Thanks. > > regards, > --- Chun-ying Lee > >>it is not clear to me what you are trying to do: >>you seem to have a time-concentration-curve in PKindex and you seem >>to set up a derivative of this time dependency according to some >>model in dCpdt. AFAIKS this scenario is not directly related to the >>situation described by the Michaelis-Menten-Equation which relates >>some "input" concentration with some "product" concentration. If Vm and >>Km are meant to be the canonical symbols, >>what is Vd, a volume of distribution? it is impossible to see (at least >>for me) what exactly you want to achieve. >> >>(and in any case, I would prefer "nls" for a least squares fit >>instead of 'optim'). >> >>joerg >> >>> >>> >>>__ >>>R-help@stat.math.ethz.ch mailing list >>>https://stat.ethz.ch/mailman/listinfo/r-help >>>PLEASE do read the posting guide! http://www.R-project.org/posting- > > guide.html > > > > > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Michaelis-menten equation
Chun-Ying Lee wrote:
> Dear R users:
>I encountered difficulties in michaelis-menten equation. I found
> that when I use right model definiens, I got wrong Km vlaue,
> and I got right Km value when i use wrong model definiens.
> The value of Vd and Vmax are correct in these two models.
>
> #-right model definiens
> PKindex<-data.frame(time=c(0,1,2,4,6,8,10,12,16,20,24),
>conc=c(8.57,8.30,8.01,7.44,6.88,6.32,5.76,5.20,4.08,2.98,1.89))
> mm.model <- function(time, y, parms) {
>dCpdt <- -(parms["Vm"]/parms["Vd"])*y[1]/(parms["Km"]+y[1])
>list(dCpdt)}
> Dose<-300
> modfun <- function(time,Vm,Km,Vd) {
>out <- lsoda(Dose/Vd,time,mm.model,parms=c(Vm=Vm,Km=Km,Vd=Vd),
> rtol=1e-8,atol=1e-8)
> out[,2] }
> objfun <- function(par) {
>out <- modfun(PKindex$time,par[1],par[2],par[3])
>sum((PKindex$conc-out)^2) }
> fit <- optim(c(10,1,80),objfun, method="Nelder-Mead)
> print(fit$par)
> [1] 10.0390733 0.1341544 34.9891829 #--Km=0.1341544,wrong value--
>
>
> #-wrong model definiens
> #-Km should not divided by Vd--
> PKindex<-data.frame(time=c(0,1,2,4,6,8,10,12,16,20,24),
>conc=c(8.57,8.30,8.01,7.44,6.88,6.32,5.76,5.20,4.08,2.98,1.89))
> mm.model <- function(time, y, parms) {
>dCpdt <- -(parms["Vm"]/parms["Vd"])*y[1]/(parms["Km"]/parms["Vd"]+y[1])
>list(dCpdt)}
> Dose<-300
> modfun <- function(time,Vm,Km,Vd) {
> out <- lsoda(Dose/Vd,time,mm.model,parms=c(Vm=Vm,Km=Km,Vd=Vd),
> rtol=1e-8,atol=1e-8)
>out[,2]
> }
> objfun <- function(par) {
> out <- modfun(PKindex$time,par[1],par[2],par[3])
> sum((PKindex$conc-out)^2)}
> fit <- optim(c(10,1,80),objfun, method="Nelder-Mead)
> print(fit$par)
> [1] 10.038821 4.690267 34.989239 #--Km=4.690267,right value--
>
> What did I do wrong, and how to fix it?
> Any suggestions would be greatly appreciated.
> Thanks in advance!!
>
>
>
it is not clear to me what you are trying to do:
you seem to have a time-concentration-curve in PKindex and you seem to
set up a derivative of this time dependency according
to some model in dCpdt. AFAIKS this scenario is not directly related to
the situation described by the Michaelis-Menten-Equation which relates
some "input" concentration with some "product" concentration. If Vm and
Km are meant to be the canonical symbols,
what is Vd, a volume of distribution? it is impossible to see (at least
for me) what exactly you want to achieve.
(and in any case, I would prefer "nls" for a least squares fit instead
of 'optim').
joerg
>
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] unexpected par('pin') behaviour
Martin Maechler wrote:
>>>>>>"joerg" == joerg van den hoff <[EMAIL PROTECTED]>
>>>>>>on Wed, 13 Jul 2005 16:00:58 +0200 writes:
>
>
> joerg> hi everybody,
> joerg> I noticed the following: in one of my scripts 'layout' is used to
> joerg> generate a (approx. square) grid of variable dimensions (depending
> on
> joerg> no. of input files). if the no. of subplots (grid cells) becomes
> joerg> moderately large (say > 9) I use a construct like
>
> joerg> ###layout grid computation and set up occurs here###
> joerg> ...
> joerg> opar <- par(no.readonly = T);
> joerg> on.exit(par(opar))
> joerg> par(mar=c(4.1, 4.1, 1.1, .1))
> joerg> ###plotting occurs here
> joerg> ...
>
> joerg> to reduce the figure margins to achieve a more
> joerg> compact display. apart from 'mar' no other par()
> joerg> setting is modified.
>
> yet another example of using par('no.readonly') when it's not
> needed and inefficient.
might be. but at least it is immune against modifying some more 'par'
settings in the course of modfications at some other place in the
programm. inefficiency: should be at the ppm level of total cpu-usage in
my case, :-). what's so bad with copying back and forth this moderately
large vector?
>
> Replacing the above by
>
> ###layout grid computation and set up occurs here###
> ...
> op <- par(mar=c(4.1, 4.1, 1.1, .1))
> on.exit(par(op))
> ###plotting occurs here
>
> will be much more efficient and even solve your problem with "pin".
>
right (solves the problem). I'll adopt this change for the time being.
thank you.
> But then, yes, there might be another par() problem hidden in
> your code / example,
> but unfortunately you have not specified reproducible code.
>
>
> joerg> this works fine until the total number of subplots becomes too
> large
> joerg> ("large" depending on the current size of the X11() graphics
> device
> joerg> window, e.g. 7 x 6 subplots for the default size fo x11()).
>
> joerg> I then get the error message (only _after_ all plots are correctly
> joerg> displayed, i.e. obviously during execution of the above on.exit()
> call)
>
> joerg> Error in par(opar) :
> joerg> invalid value specified for graphics parameter "pin"
>
> joerg> and par("pin") yields:
>
> joerg> [1] 0.34864 -0.21419
>
> you mean *after* all the plotting , not the "pin" values in
> 'opar', right?
yes
>
> joerg> which indeed is invalid (negative 2nd component).
>
> joerg> I'm aware of this note from ?par:
>
> joerg> The effect of restoring all the (settable) graphics parameters as
> joerg> in the examples is hard to predict if the device has been resized.
> joerg> Several of them are attempting to set the same things in different
> joerg> ways, and those last in the alphabet will win. In particular, the
> joerg> settings of 'mai', 'mar', 'pin', 'plt' and 'pty' interact, as do
> joerg> the outer margin settings, the figure layout and figure region
> joerg> size.
>
> {{which shows you the known but not widely known fact that
> traditional par() based graphics are ``flawed by design''
> and that's why there is the package "grid" for better
> designed graphics
... which seems to my simple mind a lot more complicated to come to
terms with than the graphics package. I understand that grid is more
powerful but the subset of functionality provided by 'graphics' seems
more difficult to use in 'grid'. wrong or right?
> }}
>
> joerg> but my problem occurs without any resizing of the
> joerg> x11() window prior to resetting par to par(opar).
>
> It still would be interesting to get a reproducible example
> here, as the posting guide asks for.
>
===cut
graphics.off()
f <- function(n=7, m=6) {
nm <- n*m
layout(matrix(1:(nm),n,m))
opar <- par(no.readonly = T)
on.exit(par(opar))
par(mar = c(4.1, 4.1, 1.1, 0.1))
for (i in 1:nm) plot(i, pch=(i-1)%%25+1)
layout(1)
}
f(5) #good
par('pin')
f() #bad (at least for x11() default size)
par('pin')
===cut
> Martin
thanks for bothering.
joerg
>
>
> joerg> any ideas, what is going on?
>
> jo
[R] unexpected par('pin') behaviour
hi everybody,
I noticed the following: in one of my scripts 'layout' is used to
generate a (approx. square) grid of variable dimensions (depending on
no. of input files). if the no. of subplots (grid cells) becomes
moderately large (say > 9) I use a construct like
###layout grid computation and set up occurs here###
...
opar <- par(no.readonly = T);
on.exit(par(opar))
par(mar=c(4.1, 4.1, 1.1, .1))
###plotting occurs here
...
to reduce the figure margins to achieve a more compact display. apart
from 'mar' no other par() setting is modified.
this works fine until the total number of subplots becomes too large
("large" depending on the current size of the X11() graphics device
window, e.g. 7 x 6 subplots for the default size fo x11()).
I then get the error message (only _after_ all plots are correctly
displayed, i.e. obviously during execution of the above on.exit() call)
Error in par(opar) :
invalid value specified for graphics parameter "pin"
and par("pin") yields:
[1] 0.34864 -0.21419
which indeed is invalid (negative 2nd component).
I'm aware of this note from ?par:
The effect of restoring all the (settable) graphics parameters as
in the examples is hard to predict if the device has been resized.
Several of them are attempting to set the same things in different
ways, and those last in the alphabet will win. In particular, the
settings of 'mai', 'mar', 'pin', 'plt' and 'pty' interact, as do
the outer margin settings, the figure layout and figure region
size.
but my problem occurs without any resizing of the x11() window prior to
resetting par to par(opar).
any ideas, what is going on?
platform powerpc-apple-darwin7.9.0
arch powerpc
os darwin7.9.0
system powerpc, darwin7.9.0
status Patched
major2
minor1.0
year 2005
month05
day 12
language R
regards,
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] nls(() and trace
hi everybody, is there a canonical way to get hold of the "trace=TRUE" output from nls, i.e. to copy it to a R variable (or at least to an external log file)? I have only found the possibility to "fix(nlsModel)" (and than the correct copy of that: namespace function ...) within the R-session by modifying the trace() definition within nlsModel. not really good for everyday use ... regards, joerg ps: if this comes to douglas bates' attention: would'nt it be helpful to add the trace output as a further component(possibly optitional, depending on the trace flag) to the nls-object returned by nls()? __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] x11 and pseudo-color
for some reason the following message seems not to have reached the list in the first try, at least I can't find it. my apologies if this is my fault: we are running R under Solaris with SunRay Terminals, which are set to 8 bit color to comply with some other software. In this configuration, X11() opens with colortype=true, i.e., it is not recognized that actually the display is only 8 bit. This leads to error messages (advising to use 'pseudo.cube'). question 1: is X11() supposed to recognize the actual color capabilities? i.e. is this a bug? question 2: is there a possibility to query the color capabilities from within R in order to being able to open the X11() displays always (for true color as well as for 8 bit) with the correct colortype setting from within a function? regards, joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] X11() and pseudo color
we are running R under Solaris with SunRay Terminals, which are set to 8 bit color to comply with some other software. In this configuration, X11() opens with colortype=true, i.e., it is not recognized that actually the display is only 8 bit. This leads to error messages (advising to use 'pseudo.cube'). question 1: is X11() supposed to recognize the actual color capabilities? i.e. is this a bug? question 2: is there a possibility to query the color capabilities from within R in order to being able to open the X11() displays always (for true color as well as for 8 bit) with the correct colortype setting from within a function? regards, joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] values of bars in barplot
luc tardieu wrote: Hi, I couldn't find how to have the values written on the top of each bar in a barplot. When using hist(), it is possible to use labels=T, but this option does not seem to exist for barplot(). Is there a trick I could use to do that ? Thanks to all Luc __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html something like: x <- 1:3 #ensure some free space above the bar and remember midpoints: bmp <- barplot(x, ylim = c(0, 1.1*max(x))) text(bmp, x + .03*max(x), x) should do. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] How to break an axis?
Bo Peng wrote:
Dear list,
I need to plot four almost horizontal lines with y-values around
1,3,4, 400. If I plot them directly, the first three lines will be
indiscernible so I am thinking of breaking y-axis into two parts, one
with range (0,5), another (395,400). Is there an easy way to do this?
I can think of two ways:
1. use two plots and draw axes manually. The plot margins, are however
difficult to adjust.
2. use one plot, adjust y-values of the lines and draw y-axis
manually. But, how would I break y-axis and add separation symbols
*on* yaxis? (By separation symbol, I mean something like
--//--
Many thanks in davance.
Bo
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
maybe something like
matplot(1:10, rep(1,10)%o%c(1,3,4), col=1:3, ylim = c(0,20), type='b')
par(new=T)
matplot(1:10, rep(400,10),axes=F,ann=F, col=4, ylim = c(0,400),type='b')
axis(4)
legend(par('usr')[2], par('usr')[4], bg='white', xjust=1, c('left
axis', 'left axis', 'left axis', 'right axis'),col=1:4,
pch=as.character(1:4))
solves your problem (double y-axis instead of splitting the axis).
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] double backslashes in usage section of .Rd files
I have written a package, where a function definition includes a regexp
pattern including double backslashes, such as
myfunction <- function (pattern = ".*\\.txt$")
when I R CMD CHECK the corresponding .Rd file, I get warnings
(code/documentation mismatch), if I enforce two backslashes in the
documentation print out by
\usage { myfunction (pattern = ".*.txt$") }
have I to live with this or is their a way to avoid the warnings (apart
from being satisfied with a wrong manpage ...)?
regards
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] y-label on right hand side of plot
Is there a better way than:
par(mar=c(6,6,6,6))
plot(1:10,yaxt="n",ylab="")
axis(4)
text(12,5.5,'y-label',xpd=T,srt=90)
to get the y-ticks _and_ the y-label to the rhs of the plot? I did not
find anything in the 'par', 'plot', 'axis' and 'title' manpages to
solve the problem. (the above is ugly, because one needs to hardcode the
text position or needs to calculate it 'manually' from par('usr'). it
would be much nicer, if there were a flag to 'title' controlling were
the labels occur).
thanks
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] OS X proxy question
Drew Balazs wrote:
All,
I'm currently using R 2.0.1 on a Powerbook G4 with OS X 10.3.8. So far
the only way I've found to set my proxy is by doing
Sys.putenv("http_proxy"=":) everytime I
start up R.
This works fine, but I'd like to find a solution that doesnt require
manual input everytime I start up R. Does any one have a better way of
accomplishing this?
Thanks,
Drew Balazs
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
if I understand your problem, you need a way to get some environment
variables right which usually are set in .profile or .cshrc.
to this end you need a directory `.MacOSX' in your home dir and there a
file `environment.plist'. this can best be edited with the uitility
`Property Editor' which is in `/Developer/Applications/Utilities'.
see also http://developer.apple.com/qa/qa2001/qa1067.html
hope this helps
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Installation problem MacOS X
Hector L. Ayala-del-Rio wrote: R gurus I have tried to install the R 2.0.1 binary for OS X and although the installation was successful I can get the application going. When I double click the icon R tries to load (R window shows briefly) and it quits immediately. This behavior was described in this list before and nobody found the answer to the problem. If you try to load the x11 version by typing "R" at the command line it loads up with no problem. This means that the app is partially working and there are no permissions issue. The most interesting thing is if I log to a different account (Dummy) and I double click the application it loads with no problem. This makes me think that there has to be some type of user specific file or directory that is causing the gui to quit. Any suggestions on what file(s) could be affecting R? Thanks Hector ** Héctor L. Ayala-del-Río, Ph.D. Assistant Professor Department of Biology University of Puerto Rico at Humacao CUH postal station 100 road 908 Humacao, PR 00791 Ph: 787-850- x 9001 Fax: 787-850-9439 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html you have moved the R app already to it's final destination (instead of starting it from the mounted disk image)? I would try the "Console" utility (in /Applications/Utilities) to have a look at console.log and system.log (or inspect the logs directly: /private/var/log/system.log and /Library/Logs/Console/vdh/console.log) immediately _after_ an attempt to start the R app. probably it's a permission problem all the same. regards joerg ps: there is a R mailing list exclusively for OS X users: https://stat.ethz.ch/mailman/listinfo/r-sig-mac __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] plot function
Cuichang Zhao wrote: Hello, how can use change the plot function to change the range of axises. I want my graph from a certain range [a, b], instead of from the min to max of of datas? if i want draw a line instead of dots, should i use both plot and lines function. for example: plot(x, y); lines(x, y); things seem not working if i only use lines(x, y). Thank you so much. Cuichang Zhao Jan 25, 2005 - [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html the first "see also" entry in the help text for "plot" is "plot.default". well, do that: try "?plot.default". time estimate for finding it yourself: one minute. you really should think about using the online help. if everybody would use your approach, the list would probably get a thousand 'hits' per hour of the quality: "is there a way to speed up matrix addition? the for loops run too slow." __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] more question
Cuichang Zhao wrote: Hello, thank you very much for your help in last email. it is very helpful. right now, i have more questions about my project, 1. solve can i remove the NA from a vectors: for exmample, if my vector is: v <- (NA, 1, 2, 3, 4, 5) which should read v <- c(NA, 1, 2, 3, 4, 5): how can I remove the NA from vector v v <- v[!is.na(v)] 2. how I can get input from the user? try "?readline" 3. can I read more than one data files in one program? and how i can write something to a file ?read.table, ?write.table, ?scan, ?readLines, ?write 4. how i can display/output something on the screen? ?print, ?x11, ?plot, ?Devices thank you so much Best wish C-Ming Jan 24, 2005 one remark: I personally think it's a bad idea to refer people too quickly (with implied reproach) to the manuals (if the manuals are good, _anything_ is written down _somewhere_ but it might take substantial time to find it). but at the very least your questions 2.-4. (and the solution to no. 1 is not so hard to find either) are so elementary, that you probably have not spent a single minute reading the documentation (online help or e.g. Introduction to R). and that is in my view not a healthy ratio of "time I invested myself to solve the problem" to "time others invest to solve my problem". no offense meant, but really ... joerg - [[alternative HTML version deleted]] __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] background color for plotting symbols in 'matplot'
thanks for the response.
Uwe Ligges wrote:
joerg van den hoff wrote:
something like
matplot2(matrix(1:6,3,2),matrix(7:12,3,2),pch=21,bg=c(2,3),type='b')
Where can we find "matplot2"?
oops. that should have been 'matplot' (not 'matplot2'), of course.
does not yield the expected (at least by me) result: only the points
on the first line get (successively) background colors for the
plotting symbols, the second line gets no background color at all for
its plotting symbols.
I think the natural behaviour should be two curves which (for the
example given above) symbol-background colors 2 and 3, respectively
(as would be obtained by a suitable manual combination of 'plot' and
'lines'). the modification of the matplot code to achieve this
behaviour is obvious as far as I can see (adding 'bg' to the explicit
arguments of matplot and handling similar to 'lty', 'cex' and the like
inside the function including transfer to 'plot' and 'lines' argument
list).
is the present behaviour a bug of 'matplot' or is it for some reason
intended behaviour?
The real point is that you might want to mark by rows *or* by columns,
so it's not that easy to specify a sensible default behaviour, at least
one has to think about it.
I'm aware of this: any specific behaviour could be the 'best' default
for someone. in terms of consistency, I would argue that matplot plots
"columns of x against columns of y", so these columns should be
addressed. that is how 'lty' and 'pch' and 'cex' do it. the present
behaviour of 'bg' ('bg' interpreted only for "column 1 of x against
column 1 of y") is not sensible.
If you want to implement it for all possible arguments, the well known
problem of huge number of arguments springs to mind as well.
that is indeed a problem, but I think mainly when reading the help
pages, which then are cluttered with many not often used graphic parameters.
Since you say the "modification [...] is obvious": I think R-core
welcomes your contribution.
well, I'm not a fluent R programmer. I'm not sure if the simple minded
modification of 'matplot' would be welcome by R-core. rather, I attach
here the modified code 'matplot2' (sic!), if someone wants to use it. a
'diff' vs. the original versions shows easily the few modified lines.
joerg
Uwe Ligges
regards,
joerg
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
#clone of the standard 'matplot' version augmented by using 'bg' as an
#additional explicit argument and modifications of the code leading to
#a bevaviour of 'bg' similar to 'lty', 'pch', 'cex' et cetera (columnwise
#recycling of 'bg' entries).
matplot2 <- function (x, y, type = "p", lty = 1:5, lwd = 1, pch = NULL, col =
1:6,
cex = NULL, xlab = NULL, ylab = NULL, xlim = NULL, ylim = NULL, bg=NULL,
..., add = FALSE, verbose = getOption("verbose"))
{
paste.ch <- function(chv) paste("\"", chv, "\"", sep = "",
collapse = " ")
str2vec <- function(string) {
if (nchar(string)[1] > 1)
strsplit(string[1], NULL)[[1]]
else string
}
xlabel <- if (!missing(x))
deparse(substitute(x))
ylabel <- if (!missing(y))
deparse(substitute(y))
if (missing(x)) {
if (missing(y))
stop("Must specify at least one of `x' and `y'")
else x <- 1:NROW(y)
}
else if (missing(y)) {
y <- x
ylabel <- xlabel
x <- 1:NROW(y)
xlabel <- ""
}
kx <- ncol(x <- as.matrix(x))
ky <- ncol(y <- as.matrix(y))
n <- nrow(x)
if (n != nrow(y))
stop("`x' and `y' must have same number of rows")
if (kx > 1 && ky > 1 && kx != ky)
stop("`x' and `y' must have only 1 or the same number of columns")
if (kx == 1)
x <- matrix(x, nrow = n, ncol = ky)
if (ky == 1)
y <- matrix(y, nrow = n, ncol = kx)
k <- max(kx, ky)
type <- str2vec(type)
if (is.null(pch))
pch <- c(paste(c(1:9, 0)), letters)[1:k]
else if (is.character(pch))
pch <- str2vec(pch)
if (verbose)
cat("matplot: doing ", k, " plots with ", paste(" col= (",
paste.ch(col), ")", sep = ""), paste(" pch= (", paste.ch(pch),
")", sep
[R] background color for plotting symbols in 'matplot'
something like matplot2(matrix(1:6,3,2),matrix(7:12,3,2),pch=21,bg=c(2,3),type='b') does not yield the expected (at least by me) result: only the points on the first line get (successively) background colors for the plotting symbols, the second line gets no background color at all for its plotting symbols. I think the natural behaviour should be two curves which (for the example given above) symbol-background colors 2 and 3, respectively (as would be obtained by a suitable manual combination of 'plot' and 'lines'). the modification of the matplot code to achieve this behaviour is obvious as far as I can see (adding 'bg' to the explicit arguments of matplot and handling similar to 'lty', 'cex' and the like inside the function including transfer to 'plot' and 'lines' argument list). is the present behaviour a bug of 'matplot' or is it for some reason intended behaviour? regards, joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R does not support UTF-8 (was german umlaut problem under MacOS)
Brian D Ripley wrote:
You wrote your mail in UTF-8. R does not support UTF-8, and that is both
documented and announced on startup in such a locale (at least on OSes
with standard-conforming implementations):
thanks for clarifying this point.
nevertheless:
1. the mail was (on purpose) sent in utf-8 to transport correctly the
output from the R command window (i.e. the GUI provided with the macOS
port). it is _this_ GUI (sorry for not explaining this correctly in the
first place) where the problem occurs. I'm not using (knowingly at
least) utf-8.
when starting the same binary from the command line in a terminal (where
I generally use ISO Latin 1 encoding) it is perfectly possible to get
the special characters into variables and into plots.
2. the OS is macos 10.3, i.e. essentially FreeBSD derivative and
hopefully conforms to the standardsbu R on startup in the GUI gives only:
cut=
R : Copyright 2004, The R Foundation for Statistical Computing
Version 2.0.1 (2004-11-15), ISBN 3-900051-07-0
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for a HTML browser interface to help.
Type 'q()' to quit R.
R>
cut=
i.e. no announcement whatsoever concerning missing utf-8 support,
despite the fact that following input is interpreted in such a way.
so, probably this is more a question to the maintainers of the macOS
port:_where_ did R (when startet with the GUI) get the notion that it
should interpret keyboard input as utf-8? can I change this (it's not
in the preferences, for instance)?
gannet% env LANG=en_GB.utf8 R
R : Copyright 2004, The R Foundation for Statistical Computing
Version 2.0.1 (2004-11-15), ISBN 3-900051-07-0
...
WARNING: UTF-8 locales are not currently supported
Solution: do not use an unsupported locale.
On Wed, 15 Dec 2004, joerg van den hoff wrote:
I did not find this in the archive (hope it isn't there...):
the current release of R (2.0.1) for MacOS (10.3.6) seems not to handle
german special characters like 'Ã' correctly:
I get two characters (Atilde quarter) here.
> f <- 'Ã'
can be entered at the prompt, but echoing the variable yields
You mean printing the contents, I presume.
yes ("shell speak").
[1] "\303\274" (I think the unicode of the character)
and inserting, for instance
text(1,2,f)
in some plot seems to insert two characters (âÂ) (probably an
interpretation of the first and second group of the unicode?).
I believe, this is a R problem or is there a simple configuration switch?
thanks
joerg
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
regards,
joerg
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] german umlaut problem under MacOS
I did not find this in the archive (hope it isn't there...): the current release of R (2.0.1) for MacOS (10.3.6) seems not to handle german special characters like 'Ã' correctly: > f <- 'Ã' can be entered at the prompt, but echoing the variable yields [1] "\303\274" (I think the unicode of the character) and inserting, for instance text(1,2,f) in some plot seems to insert two characters (âÂ) (probably an interpretation of the first and second group of the unicode?). I believe, this is a R problem or is there a simple configuration switch? thanks joerg __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Modifying Code in .rda based packages (e.g. lme4)
Prof Brian Ripley wrote:
No standard R package is supplied as a .rda, including not lme4. You
must be looking at a binary installation, and you would do best to
reinstall from the sources. You could use
R --vanilla
load("/all.rda")
fix(GLMM)
save(ls(all=T), file="/all.rda", compress = TRUE)
q()
but we would not recommend it. Indeed, we do not recommend your altering
functions in other people's packages. Why not just make a copy of GLMM
with another name and alter that?
On Fri, 11 Jun 2004, Dieter Menne wrote:
assume I want to make a minor local change in a package that is supplied as
.rda. For example, I want to get rid of the non-verbose-protected
"Iteration" message in GLMM/lme4.
Probably I have to load / change / save the package, but could someone help
me to get the syntax right?
I think, saving need to be done with
save(list=ls(all=T), file="/all.rda", compress = TRUE)
otherwise R complains about
Object "ls(all = T)" not found
(the '...' argument comes first in the 'save' argument list)
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] nls and R scoping rules
Prof Brian Ripley wrote:
On Thu, 10 Jun 2004, Prof Brian Ripley wrote:
Around R 1.2.x the notion was introduced that variables should be looked
for in the environment of a formula. Functions using model.frame got
converted to do that, but nls did not. I guess that the best way forward
is to ensure that nls (and nlsModel) does search the environment of the
formula for functions.
It transpires that is rather easy to achieve. At the top of nlsModel and
nlsModel.plinear use
env <- new.env(parent=environment(form))
instead of
env <- new.env()
This then automatically searches for objects used in the formula
Notes
1) nls is in a namespace, so you need to fix the copy in the namespace or
fix the source code and rebuild.
2) This will add to the baggage needed when you save() a nls fit in a
workspace. I think that is inevitable as we cannot just identify the
funtions used in the formula (as they might call other functions in the
local environment), and it is necessary to capture the objects needed to
evaluate the formula for the predict() method to be reliable.
We could call all.names() to get the names of functions used directly in
the formula, but that seems not good enough (previous para).
Question to the cognescenti: is the price in 2) too great for this to be
done for 1.9.1? I will put the change in R-devel for now.
thank's a lot for the two responses. that seems exactly what helps me out.
remaining question: how can I edit (I mean from within R, not by messing
around with the source code directly) an invisible function object such
as "nlsModel.plinear"?
help.search('invisible'), help('function'), help('edit') at least did
not tell me.
joerg
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] nls and R scoping rules
I apologize for posting this in essence the second time (no light at the
end of the tunnel yet..):
is there a way to enforce that "nls" takes both, the data *and* the
model definition from the parent environment? the following fragment
shows the problem.
# cut here==
wrapper <- function (choose=0)
{
x <- seq(0,2*pi,len=100)
y <- sin(1.5*x);
y <- rnorm(y,y,.1*max(y))
if (choose==0) {
rm(fifu,pos=1)
fifu <- function(w,x) {sin(w*x)}
}
else
assign('fifu',function(w,x) {sin(w*x)},.GlobalEnv)
res <- nls(y ~ fifu(w,x),start=list(w=1))
res
}
# cut here==
if called as "wrapper(1)" this runs fine because the fitting function
"fifu" is assigned in the GlobalEnv.
if called as "wrapper(0)", "fifu" is defined only locally and "nls"
(actually, "nlsModel", I think) does not know what I'm talking about.
I understand, the problem is that the scoping rules are such that "nls"
does not resolve 'fifu' in the parent environment, but rather in the
GlobalEnv. (this is different for the data, which *are* taken from the
parent environment of the nls-call).
I tried some variants of using "eval" but without starting to modify
"nls" itself there seems no way (up to now, anyway).
The solution to "assign" 'fifu' directly into the GlobalEnv does work,
of course, but leads to the undesirable effect of accumulating objects
in the workspace which are not needed there (and might overwrite
existing ones).
in response to my first post, I got the hint that for "lm" the situation
is different: it handles the above situation as desired (i.e. accepts
local model definition). in so far one might even argue that this
behaviour of "lm" and "nls" leads to an inconsistent behaviour of R in
quite similar situations (e.g. replacing at some point a linar model by
a nonlinear model in some large project is not achieved by simply
replacing "lm" by "nls" somewhere deep down in the source code).
regards,
joerg
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Is there an R-version of rayplot
maybe this q&d try helps?
#=cut herer=
vectorplot <- function (field) {
#input is a (N x 4 array) of N vectors:
# field[,1:2] - x/y position of vectors
# field[,3:4] - x/y componnent of vectors
# plotted are the 2-D vectors attached to the specified positions
if (is.null(dim(field))||dim(field)[2] != 4) stop("N x 4 array expected")
loc <- field[,1:2]
val <- field[,3:4]
alpha <- rbind(loc[,1],loc[,1]+val[,1])
omega <- rbind(loc[,2],loc[,2]+val[,2])
matplot(alpha,omega,type='l',lty=1,col='black') #the vector lines
points(loc) #the start points
points(loc+val,pch=20) #the end points
}
#example:
loc=matrix(rnorm(100),50,2)
field <- cbind(loc,loc)
vectorplot(field)
#=cut herer=
there are no nice arrow heads, of course, but one could construct some...
for now, vectors start with open circles and end with filled circles.
joerg
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] scoping rules
is there a good way to get the following fragment to work when calling
it as "wrapper(1)" ?
# cut here==
wrapper <- function (choose=0)
{
x <- seq(0,2*pi,len=100)
y <- sin(1.5*x);
y <- rnorm(y,y,.1*max(y))
if (choose==0) {
rm(fifu,pos=1)
fifu <- function(w,x) {sin(w*x)}
}
else
assign('fifu',function(w,x) {sin(w*x)},.GlobalEnv)
res <- nls(y ~ fifu(w,x),start=list(w=1))
res
}
# cut here==
I understand, the problem is that the scoping rules are such that "nls"
does not resolve 'fifu' in the parent environment, but rather in the
GlobalEnv. (this is different for the data, which *are* taken from the
parent environment of the nls-call).
The solution to "assign" 'fifu' directly into the GlobalEnv (which
happens when calling 'wrapper(1)") does obviously work but leads to the
undesirable effect of accumulating objects in the workspace which are
not needed there (and might overwrite existing ones).
so: is there a way to enforce that "nls" takes the model definition from
the parent environment together with the data?
joerg
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Confidence intervals for predicted values in nls
Cristina Silva wrote:
Dear all
I have tried to estimate the confidence intervals for predicted values of a
nonlinear model fitted with nls. The function predict gives the predicted
values and the lower and upper limits of the prediction, when the class of
the object is lm or glm. When the object is derived from nls, the function
predict (or predict.nls) gives only the predicted values. The se.fit and
interval aguments are just ignored.
Could anybody tell me how to estimate the confidence intervals for the
predicted values (not the model parameters), using an object of class nls?
Regards,
Cristina
--
Cristina Silva
IPIMAR - Departamento de Recursos Marinhos
Av. de Brasília, 1449-006 Lisboa
Portugal
Tel.: 351 21 3027096
Fax: 351 21 3015948
[EMAIL PROTECTED]
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
maybe this example helps:
==cut here===
#define a model formula (a and b are the parameters, "f" is "x"):
frml <- k1 ~ f*(1-a*exp(-b/f))
#simulate some data:
a0 <- .6
b0 <- 1.2
f <- seq(0.01,4,length=20)
k1true<- f*(1-a0*exp(-b0/f))
#include some noise
amp <- .1
k1 <- rnorm(k1true,k1true,amp*k1true)
#fit:
fifu <- deriv(frml,c("a","b"),function(a,b,x){})
rr<-nls(k1~fifu(a,b,f),start=list(a=.5,b=1))
#the derivatives and variance/covariance matrix:
#(derivs could be extracted from fifu, too)
dk1.da <- D(frml[[3]],'a')
dk1.db <- D(frml[[3]],'b')
covar <- vcov(rr)
#gaussian error propagation:
a <- coef(rr)['a']
b <- coef(rr)['b']
vark1 <-
eval(dk1.da)^2*covar[1,1]+
eval(dk1.db)^2*covar[2,2]+
2*eval(dk1.da)*eval(dk1.db)*covar[1,2]
errk1 <- sqrt(vark1)
lower.bound <- fitted(rr)-2*errk1 #two sigma !
upper.bound <- fitted(rr)+2*errk1 #dito
plot(f,k1,pch=1)
ff <- outer(c(1,1),f)
kk <- outer(c(1,1),k1)*c(1-amp,1+amp)
matlines(ff,kk,lty=3,col=1)
matlines(f,cbind(k1true,fitted(rr),lower.bound,upper.bound),col=c(1,2,3,3),lty=c(1,1,2,2))
xylim <- par('usr')
xpos <- .1*(xylim[2]-xylim[1])+xylim[1]
ypos <- xylim[4] - .1*(xylim[4]-xylim[3])
legend(xpos,ypos,
c(
'data',
'true',
'fit',
'confidence'
),
pch=c(1,-1,-1,-1),
lty=c(0,1,1,2),
col=c(1,1,2,3)
)
==cut here===
if you put this in a file and source it a few times from within R you'll
get an impression how often the fit deviates from the 'true' curve more
than expected from
the shown confidence limits.
I believe this approach is 'nearly' valid as long as gaussian error
probagation is valid (i.e. only to first order in covar and therefore
for not too large std. errors, am I right?).
to my simple physicist's mind this should suffice to get 'reasonable'
(probably, in strict sense, not completely correct?) confidence
intervals for the fit/the prediction.
If somebody objects, please let me know!
joerg
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
