Re: [R] filling a list faster
another solution:
l <- vector(mode = "list", length = 10^5)
system.time(for(i in (1:10^5)) l[[i]] <- c(i,i+1,i))
On my system this version is even slightly faster than the matrix version ...
Best,
Matthias
- original message
Subject: Re: [R] filling a list faster
Sent: Fri, 13 Jul 2007
From: Philippe Grosjean<[EMAIL PROTECTED]>
> If all the data coming from your iterations are numeric (as in your toy
> example), why not to use a matrix with one row per iteration? Also, do
> preallocate the matrix and do not add row or column names before the end
> of the calculation. Something like:
>
> > m <- matrix(rep(NA, 3*10^5), ncol = 3)
> > system.time(for(i in (1:10^5)) m[i, ] <- c(i,i+1,i))
> user system elapsed
>1.362 0.033 1.424
>
> That is, about 1.5sec on my Intel Duo Core 2.33Mhz MacBook Pro, compared to:
>
> > l <- list("1"<-c(1,2,3))
> > system.time(for(i in (1:10^5)) l[[length(l)+1]] <- c(i,i+1,i))
> user system elapsed
> 191.629 49.110 248.454
>
> ... more than 4 minutes for your code.
>
> By the way, what is your "very fast machine", that is actually four
> times faster than mine (gr!)?
>
> Best,
>
> Philippe Grosjean
>
> ..<
> ) ) ) ) )
> ( ( ( ( (Prof. Philippe Grosjean
> ) ) ) ) )
> ( ( ( ( (Numerical Ecology of Aquatic Systems
> ) ) ) ) ) Mons-Hainaut University, Belgium
> ( ( ( ( (
> ..
>
> Balazs Torma wrote:
> > hello,
> >
> > first I create a list:
> >
> > l <- list("1"<-c(1,2,3))
> >
> > then I run the following cycle, it takes over a minute(!) to
> > complete on a very fast mashine:
> >
> > for(i in (1:10^5)) l[[length(l)+1]] <- c(i,i+1,i)
> >
> > How can I fill a list faster? (This is just a demo test, the elements
> > of the list are calculated iteratively in an algorithm)
> >
> > Are there any packages and documents on how to use more advanced and
> > fast data structures like linked-lists, hash-tables or trees for
> > example?
> >
> > Thank you,
> > Balazs Torma
> >
> > __
> > R-help@stat.math.ethz.ch 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.
>
--- original message end
__
R-help@stat.math.ethz.ch 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] Generic distributions
Hello Alberto, hello Greg,
in distr you can do:
library(distr)
N <- Norm(mean = 1, sd = 2)
p(N)(0.5)
r(N)(100)
!!! not: p(N, 0.5) or r(N, 100) !!!
A detailed description of package "distr" is given in package "distrDoc".
library(distrDoc)
vignette("distr")
hth
Matthias
- original message
Subject: Re: [R] Generic distributions
Sent: Tue, 06 Mar 2007
From: Greg Snow<[EMAIL PROTECTED]>
> I think the distr package does this. There are also packages that link
> to winbugs if that is what you really want to do.
>
> --
> Gregory (Greg) L. Snow Ph.D.
> Statistical Data Center
> Intermountain Healthcare
> [EMAIL PROTECTED]
> (801) 408-8111
>
>
>
> > -Original Message-
> > From: [EMAIL PROTECTED]
> > [mailto:[EMAIL PROTECTED] On Behalf Of
> > Alberto Monteiro
> > Sent: Tuesday, March 06, 2007 1:38 PM
> > To: r-help@stat.math.ethz.ch
> > Subject: [R] Generic distributions
> >
> > Is there any class that generalizes distributions?
> >
> > For example, I could say
> > x <- generic_distribution("normal", list(mean=1, sigma=0.5))
> > and then use it like rgeneric_distribution(100, x) to get a
> > sample of 100, or pgeneric_distribution(0.5, x) to get the
> > pdf at (x = 0.5).
> >
> > In the openbugs/winbugs package, that uses a language that
> > looks like R/S, we can do things like x ~ dnorm(mu, tau),
> > forget that x is a normal with mean mu and variance 1/tau,
> > and then treat it generically.
> >
> > Alberto Monteiro
> >
> > PS: this is noise... but due to spam invasion, anything that
> > increases the nonspam/spam ratio should be welcome :-)
> >
> > __
> > R-help@stat.math.ethz.ch 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.
>
--- original message end
__
R-help@stat.math.ethz.ch 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] RE gap statistic in cluster analysis
there is an implementation in package SLmisc and also in the bioconductor package SAGx. hth Matthias GRAHAM LEASK schrieb: > Has anyone implemented Tibrishani's gap statistic in R or S plus? If so I > would greatly appreciate the relevant script file. > > Any help will be much appreciated > > > Kind regards > > > Dr Graham Leask > Economics and Strategy Group > Aston Business School > Aston University > Aston Triangle > Birmingham > B4 7ET > > Tel: Direct line 0121 204 3150 > email [EMAIL PROTECTED] > [[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. > -- Dr. rer. nat. Matthias Kohl E-Mail: [EMAIL PROTECTED] Home: www.stamats.de __ R-help@stat.math.ethz.ch 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] Integration + Normal Distribution + Directory Browsing Processing Questions
Hi Nils,
I would say, pnorm is faster and has a higher precision.
Best,
Matthias
- original message
Subject: Re: [R] Integration + Normal Distribution + Directory Browsing
Processing Questions
Sent: Mon, 22 Jan 2007
From: Nils Hoeller<[EMAIL PROTECTED]>
> Thank you,
>
> both work fine. Why is pnorm to prefer?
>
> Nils
>
>
> Matthias Kohl schrieb:
> > Hi,
> >
> > why don't you use pnorm?
> > E.g.,
> >
> > pnorm(1, mean = 0.1, sd = 1.2) - pnorm(0, mean = 0.1, sd = 1.2)
> >
> > Matthias
> >
> > - original message
> >
> > Subject: Re: [R] Integration + Normal Distribution + Directory Browsing
> Processing Questions
> > Sent: Sun, 21 Jan 2007
> > From: Dimitris Rizopoulos<[EMAIL PROTECTED]>
> >
> >
> >> you can use the `...' argument of integrate, e.g.,
> >>
> >> integrate(dnorm, 0, 1)
> >> integrate(dnorm, 0, 1, mean = 0.1)
> >> integrate(dnorm, 0, 1, mean = 0.1, sd = 1.2)
> >>
> >> look at ?integrate for more info.
> >>
> >> I hope it helps.
> >>
> >> Best,
> >> Dimitris
> >>
> >>
> >> Dimitris Rizopoulos
> >> Ph.D. Student
> >> Biostatistical Centre
> >> School of Public Health
> >> Catholic University of Leuven
> >>
> >> Address: Kapucijnenvoer 35, Leuven, Belgium
> >> Tel: +32/(0)16/336899
> >> Fax: +32/(0)16/337015
> >> Web: http://med.kuleuven.be/biostat/
> >> http://www.student.kuleuven.be/~m0390867/dimitris.htm
> >>
> >>
> >> Quoting Nils Hoeller <[EMAIL PROTECTED]>:
> >>
> >>
> >>> Hi everyone,
> >>>
> >>> I am new to R, but it's really great and helped me a lot!
> >>>
> >>> But now I have 2 questions. It would be great, if someone can help me:
> >>>
> >>> 1. I want to integrate a normal distribution, given a median and sd.
> >>> The integrate function works great BUT the first argument has to be a
> >>> function
> >>>
> >>> so I do integrate(dnorm,0,1) and it works with standard m. and sd.
> >>>
> >>> But I have the m and sd given.
> >>>
> >>> So for fixed m and sd I work around with a new function mynorm
> >>>
> >>> mynorm <- function(n) {
> >>> ret <- dnorm(n,0.6,0.15)
> >>> ret
> >>> }
> >>>
> >>> for example.
> >>>
> >>> BUT what can I do for dynamic m and sd?
> >>> I want something like integrate(dnorm(,0.6,0.15),0,1), with the first
> >>> dnorm parameter open for the
> >>> integration but fixed m and sd.
> >>>
> >>> I hope you can help me.
> >>>
> >>> 2. I am working with textfiles with rows of measure data.
> >>> read.table("file") works fine.
> >>>
> >>> Now I want R to read.table all files within a given directory and
> >>> process them one by the other.
> >>>
> >>> for(all files in dir xy) {
> >>> x <- read.table(nextfile)
> >>> process x
> >>> }
> >>>
> >>> Is that possible with R? I hope so. Can anyone give me a link to
> >>>
> >> examples.
> >>
> >>> Thanks for your help
> >>>
> >>> Nils
> >>>
> >>> __
> >>> R-help@stat.math.ethz.ch mailing list
> >>> https://stat.ethz.ch/mailman/listinfo/r-help
> >>> PLEASE do read the posting guide
> >>>
> >> http://www.R-project.org/posting-guide.html
> >>
> >>> and provide commented, minimal, self-contained, reproducible code.
> >>>
> >>>
> >>>
> >>
> >> Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm
> >>
> >> __
> >> R-help@stat.math.ethz.ch mailing list
> >> https://stat.ethz.ch/mailman/listinfo/r-help
> >> PLEASE do read the posting guide
> >> http://www.R-project.org/posting-guide.html
> >> and provide commented, minimal, self-contained, reproducible code.
> >>
> >>
> >
> > --- original message end
> >
> >
> > --
> > Dr. rer. nat. Matthias Kohl
> > [EMAIL PROTECTED]
> > www.stamats.de
> >
> >
>
> __
> R-help@stat.math.ethz.ch 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.
>
--- original message end
--
Dr. rer. nat. Matthias Kohl
[EMAIL PROTECTED]
www.stamats.de
__
R-help@stat.math.ethz.ch 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] Integration + Normal Distribution + Directory Browsing Processing Questions
Hi,
why don't you use pnorm?
E.g.,
pnorm(1, mean = 0.1, sd = 1.2) - pnorm(0, mean = 0.1, sd = 1.2)
Matthias
- original message
Subject: Re: [R] Integration + Normal Distribution + Directory Browsing
Processing Questions
Sent: Sun, 21 Jan 2007
From: Dimitris Rizopoulos<[EMAIL PROTECTED]>
> you can use the `...' argument of integrate, e.g.,
>
> integrate(dnorm, 0, 1)
> integrate(dnorm, 0, 1, mean = 0.1)
> integrate(dnorm, 0, 1, mean = 0.1, sd = 1.2)
>
> look at ?integrate for more info.
>
> I hope it helps.
>
> Best,
> Dimitris
>
>
> Dimitris Rizopoulos
> Ph.D. Student
> Biostatistical Centre
> School of Public Health
> Catholic University of Leuven
>
> Address: Kapucijnenvoer 35, Leuven, Belgium
> Tel: +32/(0)16/336899
> Fax: +32/(0)16/337015
> Web: http://med.kuleuven.be/biostat/
> http://www.student.kuleuven.be/~m0390867/dimitris.htm
>
>
> Quoting Nils Hoeller <[EMAIL PROTECTED]>:
>
> > Hi everyone,
> >
> > I am new to R, but it's really great and helped me a lot!
> >
> > But now I have 2 questions. It would be great, if someone can help me:
> >
> > 1. I want to integrate a normal distribution, given a median and sd.
> > The integrate function works great BUT the first argument has to be a
> > function
> >
> > so I do integrate(dnorm,0,1) and it works with standard m. and sd.
> >
> > But I have the m and sd given.
> >
> > So for fixed m and sd I work around with a new function mynorm
> >
> > mynorm <- function(n) {
> > ret <- dnorm(n,0.6,0.15)
> > ret
> > }
> >
> > for example.
> >
> > BUT what can I do for dynamic m and sd?
> > I want something like integrate(dnorm(,0.6,0.15),0,1), with the first
> > dnorm parameter open for the
> > integration but fixed m and sd.
> >
> > I hope you can help me.
> >
> > 2. I am working with textfiles with rows of measure data.
> > read.table("file") works fine.
> >
> > Now I want R to read.table all files within a given directory and
> > process them one by the other.
> >
> > for(all files in dir xy) {
> > x <- read.table(nextfile)
> > process x
> > }
> >
> > Is that possible with R? I hope so. Can anyone give me a link to
> examples.
> >
> > Thanks for your help
> >
> > Nils
> >
> > __
> > R-help@stat.math.ethz.ch mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
> >
> >
>
>
>
> Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
--- original message end
--
Dr. rer. nat. Matthias Kohl
[EMAIL PROTECTED]
www.stamats.de
__
R-help@stat.math.ethz.ch 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] package dependency tree
Hello,
http://bioconductor.org/packages/1.9/bioc/html/pkgDepTools.html
resp.
http://bioconductor.org/packages/2.0/bioc/html/pkgDepTools.html
may help you.
Best regards,
Matthias
Gabor Grothendieck schrieb:
> Try this, noting that available.packages() returns a matrix whose columns
> include "Depends" and "Suggests" and whose rownames are the package
> names:
>
>
>> AP <- available.packages()
>> rownames(AP)[grep("quantreg", AP[, "Depends"])]
>>
> [1] "cobs""emplik" "lss" "rankreg" "rqmcmb2"
>
>> rownames(AP)[grep("quantreg", AP[, "Suggests"])]
>>
> [1] "apTreeshape" "diveMove""dyn" "ggplot" "gsubfn"
> [6] "np"
>
> On 1/2/07, roger koenker <[EMAIL PROTECTED]> wrote:
>
>> Is there a painless way to find the names of all packages on CRAN
>> that "Depend" on a specified package?
>>
>>
>> url:www.econ.uiuc.edu/~rogerRoger Koenker
>> email[EMAIL PROTECTED]Department of Economics
>> vox: 217-333-4558University of Illinois
>> fax: 217-244-6678Champaign, IL 61820
>>
>> __
>> R-help@stat.math.ethz.ch 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.
>
--
Dr. rer. nat. Matthias Kohl
E-Mail: [EMAIL PROTECTED]
Home: www.stamats.de
__
R-help@stat.math.ethz.ch 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] Questions regarding "integrate" function
Hello,
your integrand needs to be a function which accepts a numeric vector as
first argument and returns a vector of the same length (see ?integrate).
Your function does not fulfill this requirement. Hence, you have to
rewrite your function or use sapply, apply or friends; something like
newintegrand <- function(x) sapply(x, integrand)
By the way you use a very old R version and I would recommend to update
to R 2.4.0.
hth
Matthias
Le Wang schrieb:
>Hi there. Thanks for your time in advance.
>
>I am using R 2.2.0 and OS: Windows XP.
>
>My final goal is to calculate 1/2*integral of
>(f1(x)^1/2-f2(x)^(1/2))^2dx (Latex codes:
>$\frac{1}{2}\int^{{\infty}}_{\infty}
>(\sqrt{f_1(x)}-\sqrt{f_2(x)})^2dx $.) where f1(x) and f2(x) are two
>marginal densities.
>
>My problem:
>
>I have the following R codes using "adapt" package. Although "adapt"
>function is mainly designed for more than 2 dimensions, the manual
>says it will also call up "integrate" if the number of dimension
>equals one. I feed in the data x1 and x2 and bandwidths h1 and h2.
>These codes worked well when my final goal was to take double
>integrals.
>
> integrand <- function(x) {
> # x input is evaluation point for x1 and x2, a 2x1 vector
> x1.eval <- x[1]
> x2.eval <- x[2]
> # n is the number of observations
> n <- length(x1)
> # x1 and x2 are the vectors read from data.dat
> # Compute the marginal densities
> f.x1 <- sum(dnorm((x1.eval-x1)/h1))/(n*h1)
> f.x2 <- sum(dnorm((x2.eval-x2)/h2))/(n*h2)
> # Return the integrand #
> return((sqrt(f.x1)-sqrt(f.x2))**2)
>
> }
>
> estimate<-0.5*adapt(1, lo=lo.default, up=up.default,
> minpts=minpts.default, maxpts=maxpts.default,
> functn=integrand, eps=eps.default, x1, x2,h1,h2)$value
>
>
>But when I used it for one-dimension, it failed. Some of my
>colleagues suggested getting rid of "x2.eval" in the "integrand"
>because it is only one integral. But after I changed it, it still
>didn't work. R gave the error msg: "evaluation of function gave a
>result of wrong length"
>
>I am not a frequent R user..although I looked up the mailing list
>for a while and there were few postings asking similar questions, I
>can't still figure out why my codes won't work. Any help will be
>appreciated.
>
>Le
>-
>~~
>Le Wang, Ph.D
>Population Center
>University of Minnesota
>
>______
>R-help@stat.math.ethz.ch 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.
>
>
--
Dr. rer. nat. Matthias Kohl
E-Mail: [EMAIL PROTECTED]
Home: www.stamats.de
__
R-help@stat.math.ethz.ch 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] Computing time for matrix addition or subtraction
not concerning your subject line, but function "crossprod" may be useful, too Matthias - original message Subject: Re: [R] Computing time for matrix addition or subtraction Sent: Mon, 13 Nov 2006 From: Prof Brian Ripley<[EMAIL PROTECTED]> > On Sun, 12 Nov 2006, YONGWAN CHUN wrote: > > > I wonder by chance if there is a way to reduce computing time for matrix > > addition or subtraction. With a lot of iterations, it would be helpful > > to reduce a little amount time. > > Yes, by making use of an optimized BLAS: see the R-admin manual. On my > (dual CPU) system this reduced your example from 36 to 6 seconds. > > BTW, it is the calculation of PP that is taking the most of time, not as > in your subject line. > > > Simple example is as below > > > > n <- 2000 > > P <- matrix(rnorm(n*n),n,n) > > PP <- P %*% P > > M <- diag(n) - P > > R <- M + t(M) - diag(n) + PP > > > > I would like to reduce time in calculating R. > > -- > Brian D. Ripley, [EMAIL PROTECTED] > Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ > University of Oxford, Tel: +44 1865 272861 (self) > 1 South Parks Road, +44 1865 272866 (PA) > Oxford OX1 3TG, UKFax: +44 1865 272595 > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > --- original message end __ R-help@stat.math.ethz.ch 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] median of gamma distribution
Hi, to compute the median (or expectation, var, sd, IQR, mad, ...) you can also use package "distrEx". library(distrEx) (G <- Gammad()) median(G) Matthias - original Nachricht Betreff: Re: [R] median of gamma distribution Gesendet: Fri, 30. Jun 2006 Von: [EMAIL PROTECTED] > On 30-Jun-06 Philip He wrote: > > Doese anyone know a R function to find the median of a gamma > > distribution? > > qgamma will do it. Test: > > > -log(0.5) > [1] 0.6931472 > > qgamma(0.5,1) > [1] 0.6931472 > > Ted. > > > > E-Mail: (Ted Harding) <[EMAIL PROTECTED]> > Fax-to-email: +44 (0)870 094 0861 > Date: 30-Jun-06 Time: 16:53:16 > -- XFMail -- > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! > http://www.R-project.org/posting-guide.html > --- original Nachricht Ende -- Dr. rer. nat. Matthias Kohl [EMAIL PROTECTED] www.stamats.de __ R-help@stat.math.ethz.ch 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] Generic code for simulating from a distribution.
Hi,
have a look at the packages distr and distrSim which are on CRAN.
hth
Matthias
- original Nachricht
Betreff: [R] Generic code for simulating from a distribution.
Gesendet: Mo, 10. Apr 2006
Von: [EMAIL PROTECTED]
> Hello all,
>
> I have the code below to simulate samples of certain size from a
> particular distribution (here,beta distribution) and compute some
> statistics for the samples.
>
> betasim2<-function(nsim,n,alpha,beta)
> {
> sim<-matrix(rbeta(nsim*n,alpha,beta),ncol=n)
> xmean<-apply(sim,1,mean)
> xvar<-apply(sim,1,var)
> xmedian<-apply(sim,1,median)
> simset<-data.frame(sampleno=seq1:nsim),means=xmean,vars=xvar,medians=xmedian
> )
> return(simset)
> }
>
> I can write a similar coding for any distribution individually.
> Now, I would like to have a generic code, say if I specify the
> distribution with the parameters and simulation and sample size I would
> like to have the simulations done for the mentioned distribution and the
> statistics performed.
>
> I would appreciate any help in doing so?
>
> Thanks for your time.
> Mathangi
>
>
> Mathangi Gopalakrishnan
> Graduate student
> Dept of Mathematics and Statistics
> University of Maryland, Baltimore County
> Baltimore, MD
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide!
> http://www.R-project.org/posting-guide.html
>
--- original Nachricht Ende
__
R-help@stat.math.ethz.ch 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] discrete probability distribution
Hi, have a look at http://mathworld.wolfram.com/GeometricDistribution.html respectively http://mathworld.wolfram.com/NegativeBinomialDistribution.html with r = 1. In R have a look at ?rnbinom with n = 1 and in your case: prob = 1-p hth Matthias >Hi > >I want to sample from a discrete random variable X, defined on >the non-negative integers, with > >prob(X=0) = (1-p) >prob(X=1) = (1-p)*p >... >prob(X=i)=(1-p)*p^i >... > > > >Before reinventing the wheel, >is there a ready-made R function to give me such a thing? > > >-- >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 > > -- StaMatS - Statistik + Mathematik Service Dr. rer. nat. Matthias Kohl Gottlieb-Keim-Straße 60 95448 Bayreuth Phone: +49 921 50736457 E-Mail: [EMAIL PROTECTED] Home: www.stamats.de __ R-help@stat.math.ethz.ch 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] New user: Custom probability distribution
Hi, you can use our package distr respectively distrEx. require(distrEx) D1 <- DiscreteDistribution(supp=0:3, prob = 12/(25*(1:4))) plot(D1) r(D1)(50) hth Matthias >Hello, > > Given a probability function: p(x) = 12 / (25(x+1)) , x=0, 1, 2, 3 we >generate the following values: > > C1 C2 > 0 0.48 > 1 0.24 > 2 0.16 > 3 0.12 > > Now, I'm supposed to create 50 random values using this table. In >MiniTab, I simply selected Calc -> Random Data -> Discrete, and selected >the columns, and it created 50 random values in a new column.[1] > >How do I do the same thing in R? > > [1] You may be wondering why I'm telling you this. Well, it's because >if I were in your shoes, I'd think "Oh, someone wants me to solve his >homework.". I have already solved it using MiniTab, but I want to be >able to use R instead (since I prefer NetBSD). > > > > > >__ >R-help@stat.math.ethz.ch mailing list >https://stat.ethz.ch/mailman/listinfo/r-help >PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html > -- StaMatS - Statistik + Mathematik Service Dr. rer. nat. Matthias Kohl Gottlieb-Keim-Straße 60 95448 Bayreuth Phone: +49 921 50736457 E-Mail: [EMAIL PROTECTED] Home: www.stamats.de __ R-help@stat.math.ethz.ch 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] Extending a data frame with S4
the help page on "setOldClass" might help you. In particular the section
"Register or Convert?".
Matthias
hadley wickham schrieb:
>I'm trying to create an extension to data.frame with more complex row
>and column names, and have run into some problems:
>
>
>
>>setClass("new-data.frame", representation("data.frame"))
>>
>>
>[1] "new-data.frame"
>Warning message:
>old-style ('S3') class "data.frame" supplied as a superclass of
>"new-data.frame", but no automatic conversion will be peformed for S3
>classes in: .validDataPartClass(clDef, name)
>
>Do I need to be worried about this?
>
>
>
>>new("new-data.frame", data.frame())
>>
>>
>Error in initialize(value, ...) : initialize method returned an object
>of class "data.frame" instead of the required class "new-data.frame"
>
>I guess this is related to the warning above. I presume I can fix
>this with an initialize function, but I'm not sure how to go about
>referring to the data frame that is the object.
>Is there a way to extend a data.frame, or do I need to create an
>object that contains the data frame in a slot?
>
>Thanks for your help,
>
>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
>
>
--
StaMatS - Statistik + Mathematik Service
Dipl.Math.(Univ.) Matthias Kohl
www.stamats.de
__
R-help@stat.math.ethz.ch 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] convolution of the double exponential distribution
Hi David,
you can even increase the precision of these computations by changing
the variables
"DefaultNrFFTGridPointsExponent" (default: 12 -> 2^12 grid points) and
"TruncQuantile"
(default: 1e-5) in our package "distr"; see
distroptions()
You can for instance use
distroptions("DefaultNrFFTGridPointsExponent", 14)
and
distroptions("TruncQuantile", 1e-8)
We checked our algorithm in case of Binomial, Poisson, Normal and
Exponential distributions and found a very high precision; confer
Section 5 of
http://www.stamats.de/comp.pdf
Moreover, for n-fold convolutions you can also use the function
"convpow" which can be found under
http://www.uni-bayreuth.de/departments/math/org/mathe7/DISTR/nFoldConvolution.R
hth,
Matthias
>Ravi, Duncan, and Matthias,
>
>Thank you very much for the helpful replies. For convolutions with 2 or
>3 copies, I found that the CDFs from the distr package closely match the
>analytic results from this paper:
>K. Singh, M. Xie, and W. E. Strawderman, "Combining Information From
>Independent Sources Through Confidence Distributions," Annals of
>Statistics 33, no. 1 (2005): 159-183.
>
>That gives me confidence that the package will also work well for higher
>copy numbers. At least for me, using the package is much more convenient
>than programming all the needed integrals into R.
>
>David
>___
>David R. Bickel http://davidbickel.com
>Research Scientist
>Pioneer Hi-Bred International (DuPont)
>Bioinformatics and Exploratory Research
>7200 NW 62nd Ave.; PO Box 184
>Johnston, IA 50131-0184
>515-334-4739 Tel
>515-334-4473 Fax
>[EMAIL PROTECTED], [EMAIL PROTECTED]
>
>
>| -Original Message-
>| From: Matthias Kohl [mailto:[EMAIL PROTECTED]
>| Sent: Friday, December 23, 2005 9:09 AM
>| To: Bickel, David
>| Cc: Duncan Murdoch; r-help@stat.math.ethz.ch
>| Subject: Re: [R] convolution of the double exponential distribution
>|
>| Duncan Murdoch schrieb:
>|
>| >On 12/22/2005 7:56 PM, Bickel, David wrote:
>| >
>| >
>| >>Is there any R function that computes the convolution of the double
>| >>exponential distribution?
>| >>
>| >>If not, is there a good way to integrate ((q+x)^n)*exp(-2x)
>| over x from
>| >>0 to Inf for any value of q and for any positive integer n?
>| I need to
>| >>perform the integration within a function with q and n as
>| arguments. The
>| >>function integrate() is giving me this message:
>| >>
>| >>"evaluation of function gave a result of wrong length"
>| >>
>| >>
>| >
>| >Under the substitution of y = q+x, that looks like a gamma integral.
>| >The x = 0 to Inf range translates into y = q to Inf, so
>| you'll need an
>| >incomplete gamma function, such as pgamma. Be careful to get the
>| >constant multiplier right.
>| >
>| >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
>| >
>| >
>|
>| Hi,
>|
>| you can use our package "distr".
>|
>| require(distr)
>| ## define double exponential distribution
>| loc <- 0 # location parameter
>| sca <- 1 # scale parameter
>|
>| rfun <- function(n){ loc + scale * ifelse(runif(n) > 0.5, 1,
>| -1) * rexp(n) }
>| body(rfun) <- substitute({ loc + scale * ifelse(runif(n) >
>| 0.5, 1, -1) *
>| rexp(n) },
>| list(loc = loc, scale = sca))
>|
>| dfun <- function(x){ exp(-abs(x-loc)/scale)/(2*scale) }
>| body(dfun) <- substitute({ exp(-abs(x-loc)/scale)/(2*scale)
>| }, list(loc
>| = loc, scale = sca))
>|
>| pfun <- function(x){ 0.5*(1 +
>| sign(x-loc)*(1-exp(-abs(x-loc)/scale))) }
>| body(pfun) <- substitute({ 0.5*(1 +
>| sign(x-loc)*(1-exp(-abs(x-loc)/scale))) },
>| list(loc = loc, scale = sca))
>|
>| qfun <- function(x){ loc - scale*sign(x-0.5)*log(1 - 2*abs(x-0.5)) }
>| body(qfun) <- substitute({ loc - scale*sign(x-0.5)*log(1 -
>| 2*abs(x-0.5)) },
>| list(loc = loc, scale = sca))
>|
>| D1 <- new("AbscontDistribution", r = rfun, d = dfun, p =
>| pfun, q = qfun)
>| plot(D1)
>|
>| D2 <- D1 + D1 # convolution based on FFT
>| plot(D2)
>|
>| hth,
>| Matthias
>|
>| --
>| StaMatS - Statistik + Mathematik Service
>| Dipl.Math.(Univ.) Matthias Kohl
>| www.stamats.de
>|
>
Re: [R] convolution of the double exponential distribution
Duncan Murdoch schrieb:
>On 12/22/2005 7:56 PM, Bickel, David wrote:
>
>
>>Is there any R function that computes the convolution of the double
>>exponential distribution?
>>
>>If not, is there a good way to integrate ((q+x)^n)*exp(-2x) over x from
>>0 to Inf for any value of q and for any positive integer n? I need to
>>perform the integration within a function with q and n as arguments. The
>>function integrate() is giving me this message:
>>
>>"evaluation of function gave a result of wrong length"
>>
>>
>
>Under the substitution of y = q+x, that looks like a gamma integral.
>The x = 0 to Inf range translates into y = q to Inf, so you'll need an
>incomplete gamma function, such as pgamma. Be careful to get the
>constant multiplier right.
>
>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
>
>
Hi,
you can use our package "distr".
require(distr)
## define double exponential distribution
loc <- 0 # location parameter
sca <- 1 # scale parameter
rfun <- function(n){ loc + scale * ifelse(runif(n) > 0.5, 1, -1) * rexp(n) }
body(rfun) <- substitute({ loc + scale * ifelse(runif(n) > 0.5, 1, -1) *
rexp(n) },
list(loc = loc, scale = sca))
dfun <- function(x){ exp(-abs(x-loc)/scale)/(2*scale) }
body(dfun) <- substitute({ exp(-abs(x-loc)/scale)/(2*scale) }, list(loc
= loc, scale = sca))
pfun <- function(x){ 0.5*(1 + sign(x-loc)*(1-exp(-abs(x-loc)/scale))) }
body(pfun) <- substitute({ 0.5*(1 +
sign(x-loc)*(1-exp(-abs(x-loc)/scale))) },
list(loc = loc, scale = sca))
qfun <- function(x){ loc - scale*sign(x-0.5)*log(1 - 2*abs(x-0.5)) }
body(qfun) <- substitute({ loc - scale*sign(x-0.5)*log(1 - 2*abs(x-0.5)) },
list(loc = loc, scale = sca))
D1 <- new("AbscontDistribution", r = rfun, d = dfun, p = pfun, q = qfun)
plot(D1)
D2 <- D1 + D1 # convolution based on FFT
plot(D2)
hth,
Matthias
--
StaMatS - Statistik + Mathematik Service
Dipl.Math.(Univ.) Matthias Kohl
www.stamats.de
__
R-help@stat.math.ethz.ch 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] Loading namespaces
BBK schrieb:
>Just noticed the mssing ) at the end of the setClass statement, it is there
>in the orginal
>
>Phineas
>
>-Original Message-
>From: [EMAIL PROTECTED]
>[mailto:[EMAIL PROTECTED] Behalf Of BBK
>Sent: Thursday, December 08, 2005 8:18 PM
>To: 'R-Help
>Subject: [R] Loading namespaces
>
>
>I'm creating a package for my own use that uses some S4 classes but no
>methods.
>
>I have a file called NAMESPACE it contains the line:
>
>exportClasses("foo")
>
>and at the top of the R file I have
>
>setClass("foo", representation(x="numeric")
>
>and the line:
>
>.onLoad<-function(libname,pkgname)
>
>
Do you mean
.onLoad <- function(lib, pkg) require(methods)
as given in Section 1.6.6 of "Writing R Extensions"
>When I run R CMD check I get Syntax error in the only R file. If I comment
>out the .onLoad function I get a package/namespace load failed error.
>
>Are libname and pkgname parameters for the function .onLoad that need to
>explicitly stated, or does R populate them when the package is loaded?
>
>Does .onLoad as defined above do enough to ensure that the namesapce is
>loaded?
>
>
see Section 1.6.6 (ibid.): "There needs to be an .onLoad action to
ensure that the methods package is loaded and attached"
hth
Matthias
>Phineas Campbell
>
>
>
>
>>version
>>
>>
> _
>platform sparc-sun-solaris2.9
>arch sparc
>os solaris2.9
>system sparc, solaris2.9
>status
>major2
>minor1.0
>year 2005
>month04
>day 18
>language R
> _
>
>__
>R-help@stat.math.ethz.ch mailing list
>https://stat.ethz.ch/mailman/listinfo/r-help
>PLEASE do read the posting guide!
>http://www.R-project.org/posting-guide.html
>
>______
>R-help@stat.math.ethz.ch mailing list
>https://stat.ethz.ch/mailman/listinfo/r-help
>PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
>
>
--
StaMatS - Statistik + Mathematik Service
Dipl.Math.(Univ.) Matthias Kohl
www.stamats.de
__
R-help@stat.math.ethz.ch 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] Vectors of S4 Classes
Phineas schrieb:
>I have a function from which I wish to return two vectors of equal length of
>a class
>
>eg
>
>
>>setClass("ClassOne", representation(x="numeric"))
>>
>>
>[1] "ClassOne"
>
>
>>first<-new("ClassOne", x=1)
>>second<-new("ClassOne",x=2)
>>
>>
Do you want
list(first = first, second = second)
or something like this
setClass("ClassOneList", contains = "list")
new("ClassOneList", list(first = first, second = second))
hth
Matthias
>
>
>>first<-rbind(first,second)
>>first
>>
>>
>
>first
>second
>
>Is it possible to create vector or list of an S4 class?
>
>Phineas
>
>
>
>>version
>>
>>
> _
>platform i386-pc-mingw32
>arch i386
>os mingw32
>system i386, mingw32
>status
>major2
>minor1.0
>year 2005
>month04
>day 18
>language 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
>
>
--
StaMatS - Statistik + Mathematik Service
Dipl.Math.(Univ.) Matthias Kohl
Gottlieb-Keim-Straße 60
95448 Bayreuth
Phone: +49 921 50736 457
E-Mail: [EMAIL PROTECTED]
www.stamats.de
__
R-help@stat.math.ethz.ch 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] solution of convolution equation
Anna Oganyan wrote:
>Hello,
>May be somebody can help me...
>I am trying to find a solution of a convolution equation using fft (and
>unfortunately I do not have a good background for this).
>So I am just trying to figure out how it can be implemented in R. I have
>two multidimensional independent variables X and Z
>and I know their densities fx and fz, which are multidimensional arrays.
>So I have to find the density of Y, such that X+Y=Z.
>So, first I tried on a simple example, where the variables are
>one-dimensional, say X is N(0,1) and Z is N(7,3).
>So I want to find the density of Y (which should be N(7, sqrt(8)).
>I did:
>x <- seq(-6, 20, length=500)
>fx <- dnorm(x)
>z <- seq(-6, 20, length=500)
>fz <- dnorm(z, mean=7, sd=3)
>ffty <- fft(fz)/fft(fx)
>fy <- fft(ffty, inverse=T)/length(ffty)
>plot(Re(fy))
>
>But the plot is far from being normal. May be the problem is with the
>lengths of fx and fz? should they be of different lengths and fx padded
>with zeros? May be somebody could point out that…? Thanks!
>Anna
>
>__
>R-help@stat.math.ethz.ch mailing list
>https://stat.ethz.ch/mailman/listinfo/r-help
>PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
>
>
Hello Anna,
in our R package "distr" (on CRAN) we have implemented a convolution
algorithm via fft.
See also:
http://www.uni-bayreuth.de/departments/math/org/mathe7/KOHL/pubs/comp.pdf
respectively
library(distr)
getMethod("+", signature=c("AbscontDistribution","AbscontDistribution"))
(or see the sources)
Unfortunatelly, we haven't implemented our algorithms for
multidimensional distributions so far.
hope that helps,
Matthias
__
R-help@stat.math.ethz.ch 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] two-tailed exact binomail test
Peter Dalgaard wrote:
>katrina smith <[EMAIL PROTECTED]> writes:
>
>
>
>>I am trying to find a definition for the two-tailed exact binomial
>>test but have been unsuccessful. Can you help?
>>
>>
>
>Just read binom.test. The relevant bit is this:
>(m is the mean == n*p)
>
>else if (x < m) {
>i <- seq(from = ceiling(m), to = n)
>y <- sum(dbinom(i, n, p) <= d * relErr)
>pbinom(x, n, p) + pbinom(n - y, n, p, lower = FALSE)
>}
>
>i.e. we take the lower tail, including the value observed + the part
>of the upper tail where the binomial density is less than or equal to
>that of x (with a little fuzz added in). Symmetrically for observations
>in the upper tail of course.
>
>If you were looking for an "official" definition of the two sided
>exact test, I don't think one exists. R's version is equivalent to the
>likelihood ratio test, but there are alternatives (tail-balancing,
>doubling the one-sided p, and maybe more).
>
>
>
there is a reference:
Section 2.4.2 ("Zweiseitige Tests in einparametrigen
Exponentialfamilien" - two sided tests in one-parameter exponential
families) in
H. Witting (1985): Mathematische Statistik I. Teubner. Stuttgart
confer Satz 2.70, Korollar 2.73 (in case of symmetry)
and Beispiel 2.74 (application of Korollar 2.73 to binomial model for p=
0.5).
Matthias
__
R-help@stat.math.ethz.ch 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] rlm/M/MM-estimator questions
Christian Hennig wrote: >Hi list, > >1) How can the MM-estimator method="MM" in function rlm be tuned to 85% >efficiency? It seems that there is a default tuning to 95%. I presume, but >am not sure, that the MM-estimator uses phi=phi.bisquare as default and >the tuning constant could be set by adding a parameter c=... >Is this true? Which value to use for 85%? >(In principle I should be able to figure that out theoretically, but it >would be much easier if somebody already knew the constant or a >straightforward way to compute it.) > > Hi Christian, I have not calculated the efficiency myself ... But the thesis of Matias Salibian-Barrera (SB 2000) might help you to find the answer (cf. Chapter 4). See: http://mathstat.math.carleton.ca:16080/~matias/thesis.pdf As far as I understand the choice k0=1.548 is to obtain a breakdown point 0.5 whereas k0=1.988 leads to a breakdown point of 0.4 - at least in the location case; confer p. 60 of SB 2000. In the article "Optimal robust $M$-estimates of location" by Fraiman, Yohai and Zamar (Ann. Stat. 29(1): 194 - 223) which is, of course, concerned with the location case, the authors recommend to use k0=1.988 instead of k0=1.548 (cf. p. 206). Hope that helps! Matthias >2) The M-estimator with bisquare uses "rescaled MAD of the residuals" as >scale estimator according to the rlm help page. Does this mean that a >scale estimator is used which is computed from least squares residuals? Are >M-estimator and residual scale estimator iterated until convergence of >them both? (Does this converge?) Or what else? What does "rescaled" mean? > >Thank you, >Christian > > >*** NEW ADDRESS! *** >Christian Hennig >University College London, Department of Statistical Science >Gower St., London WC1E 6BT, phone +44 207 679 1698 >[EMAIL PROTECTED], www.homepages.ucl.ac.uk/~ucakche > >__ >R-help@stat.math.ethz.ch mailing list >https://stat.ethz.ch/mailman/listinfo/r-help >PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html > > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] User-defined random variable
Achim Zeileis wrote:
On Sat, 30 Apr 2005, Peter Dalgaard wrote:
Paul Smith <[EMAIL PROTECTED]> writes:
Dear All
I would like to know whether it is possible with R to define a
discrete random variable different from the ones already defined
inside R and generate random numbers from that user-defined
distribution.
Yes. One generic way is to specify the quantile function (as in
qpois() etc.) and do qfun(runif(N)).
If the support discrete but also finite, you can also use sample(), e.g.
sample(myset, N, replace = TRUE, prob = myprob)
Z
--
O__ Peter Dalgaard Blegdamsvej 3
c/ /'_ --- Dept. of Biostatistics 2200 Cph. N
(*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918
~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Hi,
one can also use our R package "distr" to generate discrete random
variables. The subsequent code provides a function which generates an
object of class "DiscreteDistribution" based on a finite support "supp".
If "prob" is missing all elements in "supp" are equally weighted.
Matthias
# generating function
DiscreteDistribution <- function(supp, prob){
if(!is.numeric(supp))
stop("'supp' is no numeric vector")
if(any(!is.finite(supp))) # admit +/- Inf?
stop("inifinite or missing values in supp")
len <- length(supp)
if(missing(prob)){
prob <- rep(1/len, len)
}else{
if(len != length(prob))
stop("'supp' and 'prob' must have equal lengths")
if(any(!is.finite(prob)))
stop("inifinite or missing values in prob")
if(!identical(all.equal(sum(prob), 1,
tolerance = 2*distr::TruncQuantile), TRUE))
stop("sum of 'prob' has to be (approximately) 1")
if(!all(prob >= 0))
stop("'prob' contains values < 0")
}
if(length(usupp <- unique(supp)) < len){
warning("collapsing to unique support values")
prob <- as.vector(tapply(prob, supp, sum))
supp <- sort(usupp)
len <- length(supp)
}else{
o <- order(supp)
supp <- supp[o]
prob <- prob[o]
}
if(len > 1){
if(min(supp[2:len] - supp[1:(len - 1)]) < distr::DistrResolution)
stop("grid too narrow --> change DistrResolution")
}
rfun <- function(n){
sample(x = supp, size = n, replace = TRUE, prob = prob)
}
intervall <- distr::DistrResolution / 2
supp.grid <- as.numeric(matrix(
rbind(supp - intervall, supp + intervall), nrow = 1))
prob.grid <- c(as.numeric(matrix(rbind(0, prob), nrow = 1)), 0)
dfun <- function(x){ stepfun(x = supp.grid, y = prob.grid)(x) }
cumprob <- cumsum(prob)
pfun <- function(x){ stepfun(x = supp, y = c(0, cumprob))(x) }
qfun <- function(x){ supp[sum(cumprob
return(new("DiscreteDistribution", r = rfun, d = dfun, p = pfun,
q = qfun, support = supp))
}
# some examples
supp <- rpois(20, lambda=7) # some support vector
D1 <- DiscreteDistribution(supp = supp) # prob is missing
r(D1)(10) # 10 random numbers of Distribution D1
support(D1) # support
d(D1)(support(D1)) # pdf
p(D1)(5) # cdf
q(D1)(0.5) # quantile (here: median)
plot(D1) # plot of pdf, cdf and quantile
D2 <- lgamma(D1) # apply member of group generic "Math"
plot(D2)
D3 <- D1 + Binom(size=10) # convolution with object of class
"DiscreteDistribution"
plot(D3)
D4 <- D1 + Norm() # convolution with object of class "AbscontDistribution"
plot(D4)
__
R-help@stat.math.ethz.ch 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] Using R to illustrate the Central Limit Theorem
(Ted Harding) wrote:
On 21-Apr-05 [EMAIL PROTECTED] wrote:
Here's a bit of a refinement on Ted's first suggestion.
[ corrected from runif(M*k), N, k) to runif(N*k), N, k) ]
N <- 1
graphics.off()
par(mfrow = c(1,2), pty = "s")
for(k in 1:20) {
m <- (rowMeans(matrix(runif(N*k), N, k)) - 0.5)*sqrt(12*k)
hist(m, breaks = "FD", xlim = c(-4,4), main = k,
prob = TRUE, ylim = c(0,0.5), col = "lemonchiffon")
pu <- par("usr")[1:2]
x <- seq(pu[1], pu[2], len = 500)
lines(x, dnorm(x), col = "red")
qqnorm(m, ylim = c(-4,4), xlim = c(-4,4), pch = ".", col = "blue")
abline(0, 1, col = "red")
Sys.sleep(1)
}
Very nice! (I can better keep up with it mentally, though, with
Sys.sleep(2) or Sys.sleep(3), which moght be better for classroom
demo).
One thing occurred to me, watching it: people might say "Yes,
we can see how the distribution -> Normal, nice and smooth,
especially in the tails and side-arms; but the peaks always look
a bit rough."
Which could be the cue for introducing "SD(ni) = sqrt(E[ni])",
and the following hack of the above code seems to show this OK
in the "rootograms":
N <- 1
graphics.off()
par(mfrow = c(1,2), pty = "s")
for(k in 1:20) {
m <- (rowMeans(matrix(runif(N*k), N, k)) - 0.5)*sqrt(12*k)
hm <- hist(m, breaks = "FD", xlim = c(-4,4), main = k, plot=FALSE,
prob = TRUE, ylim = c(0,0.5), col = "lemonchiffon")
hm$counts<-sqrt(hm$counts) ;
plot(hm,xlim = c(-4,4),main = k,ylab="sqrt(Frequency)")
pu <- par("usr")[1:2]
x <- seq(pu[1], pu[2], len = 500)
lines(x, sqrt(N*dnorm(x)*(hm$breaks[2]-hm$breaks[1])), col = "red")
qqnorm(m, ylim = c(-4,4), xlim = c(-4,4), pch = ".", col = "blue")
abline(0, 1, col = "red")
Sys.sleep(2)
}
(and also shows clearly how the tails of the sample move outwards
into the tails of the Normal, as in fact you expect from the finite
range of mean(runif(k)), especially initially: very visible for
k up to about 5, and not really settled down for k<10).
Next stop: Hanging rootograms!
Best wishes,
Ted.
E-Mail: (Ted Harding) <[EMAIL PROTECTED]>
Fax-to-email: +44 (0)870 094 0861
Date: 22-Apr-05 Time: 13:10:19
-- XFMail --
__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Hello,
here is another idea to illustrate the Central Limit Theorem which is
based on our package "distr".
# Illustration of the Central Limit Theorem
# using package "distr"
require(distr)
CLT <- function(Distr, n, sleep = 1){
# Distr: object of class "AbscontDistribution"
# n: iterations
# sleep: time to sleep
graphics.off()
par(mfrow = c(1,2))
# expectation of Distr
fun1 <- function(x, Distr){x*d(Distr)(x)}
E <- try(integrate(fun1, lower = q(Distr)(0), upper = q(Distr)(1),
Distr = Distr)$value,
silent = TRUE)
if(!is.numeric(E))
E <- try(integrate(fun1, lower = q(Distr)(distr::TruncQuantile),
upper = q(Distr)(1-distr::TruncQuantile),
Distr = Distr)$value,
silent = TRUE)
# standard deviation of Distr
fun2 <- function(x, Distr){x^2*d(Distr)(x)}
E2 <- try(integrate(fun2, lower = q(Distr)(0), upper = q(Distr)(1),
Distr = Distr)$value,
silent = TRUE)
if(!is.numeric(E2))
E2 <- try(integrate(fun2, lower = q(Distr)(distr::TruncQuantile),
upper = q(Distr)(1-distr::TruncQuantile),
Distr = Distr)$value,
silent = TRUE)
std <- sqrt(E2 - E^2)
Sn <- 0
N <- Norm()
for(k in 1:n) {
Sn <- Sn + Distr
Tn <- (Sn - k*E)/(std*sqrt(k))
x <- seq(-5,5,0.01)
dTn <- d(Tn)(x)
ymax <- max(1/sqrt(2*pi), dTn)
plot(x, d(Tn)(x), ylim = c(0, ymax), type = "l", ylab =
"densities", main = k, lwd = 4)
lines(x, d(N)(x), col = "orange", lwd = 2)
plot(x, p(Tn)(x), ylim = c(0, 1), type = "l", ylab = "cdfs",
main = k, lwd = 4)
lines(x, p(N)(x), col = "orange", lwd = 2)
Sys.sleep(sleep)
}
}
# some examples
distroptions("DefaultNrFFTGridPointsExponent", 13)
CLT(Distr = Unif(), n = 20, sleep = 0)
CLT(Distr = Exp(), n = 20, sleep = 0)
CLT(Distr = Chisq(), n = 20, sleep = 0)
CLT(Distr = Td(df = 5), n = 20, sleep = 0)
CLT(Distr = Beta(), n = 20, sleep = 0)
distroptions("DefaultNrFFTGridPointsExponent", 14)
CLT(Distr = Lnorm(), n = 20, sleep = 0)
__
R-help@stat.math.ethz.ch 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] usage and behavior of 'setIs'
Thank you,
Matthias
> Hi Matthias,
>
> A similar problem to yours (with one level of inheritance less) was
> disccussed this month on the r-devel list.
> You find an answer from JChambers here:
>
> https://stat.ethz.ch/pipermail/r-devel/2004-October/030980.html
>
> And yes specifying _setAs_ to each _setIs_ with the coerce and replace
> is a _hack_ which is with this version of methods necessary when
> inherting from Old Classes.
>
> /E
>
>
> [EMAIL PROTECTED] wrote:
>
>>Hello,
>>
>>am I using 'setIs' in the correct way in the subsequent (artifical)
>> example?
>>
>>Do I have to specify explicit 'setAs' for 'list' and 'vector' or
>>should this work automatically, since "getClass("List1")" states
>>an explicit coerce also for these classes.
>>
>>I'm working with R 2.0.0 Patched (2004-10-06) on windows 2000.
>>
>>Thanks for your advice,
>>Matthias
>>
>>
>># example
>>setClass(Class = "List1", representation(List = "list"))
>>setClass(Class = "List2", contains = "list")
>>
>>setIs(class1 = "List1", class2 = "List2",
>>coerce = function(obj){ new("List2", [EMAIL PROTECTED]) },
>>replace = function(obj, value){
>>[EMAIL PROTECTED] <- value
>>})
>>
>>getClass("List1")
>># states explicit coerce for 'list' and 'vector'
>>getClass("List2")
>>L1 <- new("List1", List = list("a"))
>>
>># all TRUE
>>is(L1, "List2")
>>is(L1, "list")
>>is(L1, "vector")
>>
>>as(L1, "List2") # works
>>
>># both return 'list()'
>># why not a 'list' with entry "a"?
>># Is there an additional 'setAs' needed?
>>as(L1, "list")
>>as(L1, "vector")
>>
>>L2 <- as(L1, "List2")
>>as(L2, "list") # works
>>as(L2, "vector") # works
>>
>>__
>>[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
>>
>>
>>
>
>
> --
> Dipl. bio-chem. Witold Eryk Wolski
> MPI-Moleculare Genetic
> Ihnestrasse 63-73 14195 Berlin
> tel: 0049-30-83875219 __("<_
> http://www.molgen.mpg.de/~wolski \__/'v'
> http://r4proteomics.sourceforge.net||/ \
> mail: [EMAIL PROTECTED]^^ m m
> [EMAIL PROTECTED]
>
> __
> [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
__
[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] usage and behavior of 'setIs'
Hello,
am I using 'setIs' in the correct way in the subsequent (artifical) example?
Do I have to specify explicit 'setAs' for 'list' and 'vector' or
should this work automatically, since "getClass("List1")" states
an explicit coerce also for these classes.
I'm working with R 2.0.0 Patched (2004-10-06) on windows 2000.
Thanks for your advice,
Matthias
# example
setClass(Class = "List1", representation(List = "list"))
setClass(Class = "List2", contains = "list")
setIs(class1 = "List1", class2 = "List2",
coerce = function(obj){ new("List2", [EMAIL PROTECTED]) },
replace = function(obj, value){
[EMAIL PROTECTED] <- value
})
getClass("List1")
# states explicit coerce for 'list' and 'vector'
getClass("List2")
L1 <- new("List1", List = list("a"))
# all TRUE
is(L1, "List2")
is(L1, "list")
is(L1, "vector")
as(L1, "List2") # works
# both return 'list()'
# why not a 'list' with entry "a"?
# Is there an additional 'setAs' needed?
as(L1, "list")
as(L1, "vector")
L2 <- as(L1, "List2")
as(L2, "list") # works
as(L2, "vector") # works
__
[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] Gumbel distribution
please, take a look at package "evd" Matthias > Does R have built in Gumbel distribution (pdf, ecdf, hazard, parameters > estimation) for the minimum case? Thanks > > Anne > [[alternative HTML version deleted]] > > __ > [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 __ [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] setClassUnion
sorry, please replace exportClass by exportClasses
> Hello,
>
> I have a question concerning "setClassUnion".
> I'm working with R 2.0.0 Patched (2004-10-06) on windows 2000.
>
> I tried to use "setClassUnion" in a package I am currently working on.
> The situation is similar to the following example:
>
> The DESCRIPTION file has entries:
> Depends: R (>= 2.0.0), methods
> Imports: methods
> LazyLoad: yes
>
> The NAMESPACE file has entries:
> importClassesFrom("methods", "NULL", "numeric")
> exportClass("OptionalNumeric", "class1", "class2")
exportClasses("OptionalNumeric", "class1", "class2")
>
> The example R code is:
> .onLoad <- function(lib, pkg){
> require("methods", character = TRUE, quietly = TRUE)
> }
>
> setClassUnion("OptionalNumeric", c("numeric", "NULL"))
>
> setClass("class1",
> representation(test1 = "OptionalNumeric"),
> prototype(test1 = numeric(1)))
>
> # why does this not work?
> # The error I get is:
> # Error in makePrototypeFromClassDef(properties, ClassDef, immediate, #
> where) :In making the prototype for class "class1" elements of the #
> prototype failed to match the corresponding slot class: test1
> # (class "OptionalNumeric ")
> # Sourcing this into R gives no error for me
>
> # but instead using
> prototype(test1 = NULL)
> # works
>
> # Moreover, using the second version (with test1 = NULL)
> # the following works, too
> setClass("class2",
> representation(test2 = "class1"),
> prototype(test2 = new("class1", test1 = numeric(1
>
> What am I doing wrong?
> Can someone please explain this to me?
>
> Thanks for your help,
> Matthias
>
> __
> [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
__
[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] setClassUnion
Hello,
I have a question concerning "setClassUnion".
I'm working with R 2.0.0 Patched (2004-10-06) on windows 2000.
I tried to use "setClassUnion" in a package I am currently working on. The
situation is similar to the following example:
The DESCRIPTION file has entries:
Depends: R (>= 2.0.0), methods
Imports: methods
LazyLoad: yes
The NAMESPACE file has entries:
importClassesFrom("methods", "NULL", "numeric")
exportClass("OptionalNumeric", "class1", "class2")
The example R code is:
.onLoad <- function(lib, pkg){
require("methods", character = TRUE, quietly = TRUE)
}
setClassUnion("OptionalNumeric", c("numeric", "NULL"))
setClass("class1",
representation(test1 = "OptionalNumeric"),
prototype(test1 = numeric(1)))
# why does this not work?
# The error I get is:
# Error in makePrototypeFromClassDef(properties, ClassDef, immediate,
# where) :In making the prototype for class "class1" elements of the
# prototype failed to match the corresponding slot class: test1
# (class "OptionalNumeric ")
# Sourcing this into R gives no error for me
# but instead using
prototype(test1 = NULL)
# works
# Moreover, using the second version (with test1 = NULL)
# the following works, too
setClass("class2",
representation(test2 = "class1"),
prototype(test2 = new("class1", test1 = numeric(1
What am I doing wrong?
Can someone please explain this to me?
Thanks for your help,
Matthias
__
[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] How to pass strings to functions? [once once more, now With content I hope...]
> Dear expeRts,
>
> I fail to succesfully pass strings to functions. It comes down to the
> observation that
>
> > plot(someVariable,anotherVariable)
>
> works fine, but
>
> > x <- "someVariable"
> > y <- "anotherVariable"
> > plot(x,y)
>
> does not.
>
> Does this have something to do with the returned value of x being
> /"someVariable"/ and not /someVariable/, i.e. without the quotation
> marks? Is there any way to work around this?
>
> Ultimately I'd like to make multiple graphs by looping throught the
> values in vectors. Something like:
> > var<-c(var1,var2...n)
> > for (v in var)
> >{
> > plot(var, x))
> >}
>
what about plot(get(v), x)?
> I've looked around for help on this but am stuck.
>
> Hope you can help,
> Gijs Plomp
>
> __
> [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
[R] questions about setMethod("Arith", ...)
Hi,
we have some questions concerning the definition of new arithmetic methods.
In our package "distr" (on CRAN) we define some new arithmetic methods
for "+", "-", "*", "/".
After loading "distr" the corresponding arithmetic methods work. Now, if
we define a new
class and also a new method for one of the arithmetic methods "+", "-",
"*", "/", still everything works fine. But, if we then define new
methods for the whole group "Arith", the arithmetic methods of "distr"
no longer work. (for more details see example code below)
What are we missing?
Moreover, there is a certain precendence; i.e., if one defines a single
arithmetic method (e.g., "/") and alterwards defines a method for the
whole group "Arith", the "old" method "/" remains valid.
However, if we first define a method for the whole group "Arith" and
afterwards define a new single arithmetic method (e.g., "+") the new one
is valid. (for more details see example code below).
Is this intended?
Thanks for your help,
Matthias, Thomas
###
## Example code
###
require(distr)
getMethods("/") # shows the corresponding methods of "distr"
## now define a new class "track" (see Chambers (1998))
## and define "/"
setClass("track", representation(x = "numeric", y = "numeric"))
setMethod("/", signature("track", "numeric"),
function(e1, e2){ [EMAIL PROTECTED] = [EMAIL PROTECTED]/e2; e1 })
getMethods("/") # shows the corresponding methods
# of "distr" and class "track"
(N <- Norm()) # creates an object of standard normal distribution
(N1 <- N/3) # works
(tr <- new("track", x = 1:3, y = 4:6))
(tr1 <- tr/3) # works
## now define new methods for "Arith"
setMethod("Arith", signature("track", "numeric"),
function(e1, e2){
[EMAIL PROTECTED] = callGeneric([EMAIL PROTECTED], e2)
e1
})
getMethods("/") # "/" for "distr" is lost
N2 <- N/3 # fails
(tr2 <- tr/3) # works, "but" still the "old" method
tr + 2 # works
## now a new method "+"
setMethod("+", signature("track", "numeric"),
function(e1, e2){ [EMAIL PROTECTED] = [EMAIL PROTECTED]; e1 })
tr + 2 # works, "but" with the "new" method
__
[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] Fw: Evaluating strings as variables
> "Robin Gruna" <[EMAIL PROTECTED]> writes:
>
>> I have the following problem: I have a list as follows,
>>
>> > values <- list(red = 1, yellow = 2, blue = 3)
>>
>> > values
>> $red
>> [1] 1
>>
>> $yellow
>> [1] 2
>>
>> $blue
>> [1] 3
>>
>> There is also a vector containing the diffrent "colors" as character
>> strings:
>>
>> > colors <- c("red", "red", "blue", "yellow", "red", "blue")
>> > colors
>> [1] "red""red""blue" "yellow" "red""blue"
>>
>> Now i can attach the list values to R:
>>
>> > attach(values)
>> > red
>> [1] 1
>>
>> etc...
>>
>> Now to my problem: How can I make R in a simple way to evaluate the
>> strings "red", "blue" etc. as variables, returning their numeric
>> values ? As result I want to get a vector containing the values of the
>> colors like this one:
>>
>> > values.colors
>> [1] 1 1 3 2 1 3
>
> Forget about the attach() business, and do
>
> values <- unlist(values) # or values <- c(red = 1, yellow = 2, blue = 3)
> values.colors <- values[colors]
>
> If you insist on going via variables, try
>
> sapply(colors, get)
>
or?
unlist(mget(colors, envir = as.environment(-1), inherits = TRUE))
Matthias
> --
>O__ Peter Dalgaard Blegdamsvej 3
> c/ /'_ --- Dept. of Biostatistics 2200 Cph. N
> (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918
> ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907
>
> __
> [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
[R] setValidity() changes "Extends"?
Hi,
I'm using Version 1.9.0 (2004-04-12) on Windows NT/98/2000 and found
the following difference between using setClass(..., valdity = ),
respectively using setValidity() afterwards:
setValidity() changes the "Extends"-part of a derived class, is this
intended or a bug or am I missing something?
##
## Example code
##
setClass("Class1", representation("name" = "character"))
if(!isGeneric("name")) setGeneric("name", function(object)
standardGeneric("name"))
setMethod("name", "Class1", function(object) [EMAIL PROTECTED])
setClass("Class2", representation("Class1"))
setClass("Class3", representation("Class2"))
getClass("Class3") # as I expected
#Slots:
#
#Name: name
#Class: character
#
#Extends:
#Class "Class2", directly
#Class "Class1", by class "Class2"
validClass3 <- function(object){TRUE}
setValidity("Class3", validClass3)
#Slots:
#
#Name: name
#Class: character
#
#Extends: "Class2" # has been changed???
getClass("Class3") # why does setValidity change "Extends"?
# am I missing something?
# This doesn't happen if I use
# setClass(..., validity = )
# It, of course, also works if I
explicitly use
# setClass("Class3", contains =
c("Class2", "Class1")
C32 <- new("Class3")
name(C32) # generates an error?!
#Error in name(C32) : No direct or inherited method for function "name"
for this call
is(C32, "Class1") # however
#TRUE
__
[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] "privileged slots"
Martin Maechler schrieb:
"Matthias" == Matthias Kohl <[EMAIL PROTECTED]>
on Thu, 27 May 2004 14:01:51 +0100 writes:
Matthias> Hi all, in the help for RClassUtils I found the
Matthias> expression "privileged slots" in function
Matthias> "checkSlotAssignment" with the explanation:
Matthias> /privileged slots (those that can only be set by
Matthias> accesor functions defined along with the class
Matthias> itself)/
RClassUtils ???
> help.search("RClassUtils")
your right, sorry
but, at least a R Site search in "Functions"
gives me one match: "Utilities for Managing Class Definitions"
which hast the "title": RClassUtils{methods} R Documentation
No help files found with alias or concept or title matching
'RClassUtils' using fuzzy matching.
-
So I guess that's not something in a standard R document.
You should rather keep to the 'official documentation' ...
I thought this is a official documentation ...
Matthias> I thought all slots of a (not private) class can
Matthias> be a accessed by a user via the @ Operator.
I tend to agree with your thoughts...
Matthias> Is there a way to make a single slot of a class (not
Matthias> the whole class) private, so that you can access
Matthias> this slot only via an accessor function (not via @)?
I'd rather guess not.
Matthias> Thanks, for your help Matthias
Martin
__
[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
[R] "privileged slots"
Hi all, in the help for RClassUtils I found the expression "privileged slots" in function "checkSlotAssignment" with the explanation: /privileged slots (those that can only be set by accesor functions defined along with the class itself)/ I thought all slots of a (not private) class can be a accessed by a user via the @ Operator. Is there a way to make a single slot of a class (not the whole class) private, so that you can access this slot only via an accessor function (not via @)? Thanks, for your help Matthias __ [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] (offtopic) I need two sets of 5 different color scales
maybe, the RColorBrewer package does what you want? see also: ColorBrewer.org > Hi, > > I am plotting a policy function (result from a dynamic stochastic > optimization problem, discretized approximation). The policy function > maps from an 2 x 2 x 2 x 3 x B x F state space to a B x F state space (B > and F are usually between 4-6, and represent domestic and foreign > savings. The other variables are income (Y), inflation (Pi), domestic > and foreign interest rates (R and Z)). I actually wrote a plotting > function to represent all this, the result is attached -- please have a > look at it and help me... > > I need advice in the following: I need two sets of colors for B and F > which are easy to distinguish (when printed on a color laser printer), > represent cardinality (ie have an intuitive mapping to an interval) or > at least ordinality. > > I have experimented with the following: > > Bcolors <- hsv(.6, seq(0.2, 1, length=5), 1) > Fcolors <- hsv(seq(.1,0, length=5), seq(0.2, 1, length=5) > > this is what you see in the plot. What colors would you use? Do you > think that varying both brightness and hue helps to distinguish > colors? Should I change saturation, too? > > Thanks, > > Tamas > > PS.: The plot is simply gzipped. If you need a zipped version, or the > source code, contact me. > > -- > Tamás K. Papp > E-mail: [EMAIL PROTECTED] > Please try to send only (latin-2) plain text, not HTML or other garbage. __ [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] Exactness of ppois
Hello,
by checking the precision of a convolution algorithm, we found the
following "inexactness":
We work with R Version 1.8.1 (2003-11-21) on Windows systems (NT, 2000,
XP).
Try the code:
## Kolmogorov distance between two methods to
## determine P(Poisson(lambda)<=x)
Kolm.dist <- function(lam, eps){
x <- seq(0,qpois(1-eps, lambda=lam), by=1)
max(abs(ppois(x, lambda=lam)-cumsum(dpois(x, lambda=lam
}
erg<-optimize(Kolm.dist, lower=900, upper=1000, maximum=TRUE, eps=1e-15)
erg
Kolm1.dist <- function(lam, eps){
x <- seq(0,qpois(1-eps, lambda=lam), by=1)
which.max(abs(ppois(x, lambda=lam)-cumsum(dpois(x, lambda=lam
}
Kolm1.dist(lam=erg$max, eps=1e-15)
So for lambda=977.8 and x=1001 we get a distance of about 5.2e-06.
(This inexactness seems to hold for all lambda values greater than about
900.)
BUT, summing about 1000 terms of exactness around 1e-16,
we would expect an error of order 1e-13.
We suspect algorithm AS 239 to cause that flaw.
Do you think this could cause other problems apart from
that admittedly extreme example?
Thanks for your attention!
Matthias
__
[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] var on a vector
[EMAIL PROTECTED] schrieb: Hello I am an "R" newbie. I have a problem with computing a variance on a vector. data(cars) variance <- function (x) mean(x^2)-mean(x)^2; variance(cars[,1]) [1] 27.4 var(cars[,1]) [1] 27.95918 What did I assume/understand wrong ? TIA Hello, help(var) says: The denominator n - 1 is used which gives an unbiased estimator of the (co)variance for i.i.d. observations. Try: n <- length(cars[,1]) var(cars[,1])*(n-1)/n Matthias Kohl __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
