Re: [R] sequence with start and stop positions
Try:
unlist(mapply(seq, c(1,20,50), c(7,25,53)))
[1] 1 2 3 4 5 6 7 20 21 22 23 24 25 50 51 52 53
-Christos
-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Chris Oldmeadow
Sent: Tuesday, August 26, 2008 12:42 AM
To: r-help@r-project.org
Subject: [R] sequence with start and stop positions
Hi,
I have a vector of start positions, and another vector of
stop positions,
eg start-c(1,20,50)
stop-c(7,25,53)
Is there a quick way to create a sequence from these vectors?
new-c(1,2,3,4,5,6,7,20,21,22,23,24,25,50,51,52,53)
the way Im doing it at the moment is
pos-seq(start[1],stop[1])
for (i in 2:length(start)){
new-seq(start[i],stop[i])
pos-c(pos,new)
}
This works on small data, but its very inefficient on large
vectors, and is taking forever!
Anybody no a better way?
many thanks,
Chris
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] lattice plotting character woes
Hi Rolf, Hi Mark, Hi List,
I have not digested Rolf's response yet. It may well answer my problems.
In the meantime I have some reproducible code which actually shows my
problem:
patches -
structure(list(URBAN_AREA = structure(c(2L, 19L, 23L, 2L, 19L,
23L, 2L, 19L, 23L, 2L, 19L, 23L), .Label = c(CENTRAL AUCKLAND ZONE,
CHRISTCHURCH, DUNEDIN, HAMILTON ZONE, HASTINGS ZONE,
INVERCARGILL, LOWER HUTT ZONE, mean, NAPIER ZONE, NELSON,
NEW PLYMOUTH, NORTHERN AUCKLAND ZONE, PALMERSTON NORTH,
PORIRUA ZONE, ROTORUA, SD, SE, SOUTHERN AUCKLAND ZONE,
TAURANGA, WANGANUI, WELLINGTON ZONE, WESTERN AUCKLAND ZONE,
WHANGAREI), class = factor), NO_PATCHES = c(11L, 16L, 21L,
87L, 192L, 324L, 164L, 417L, 773L, 679L, 757L, 3083L), MEAN_AREA =
c(9.623631225,
15.29089619, 149.2063532, 14.1676, 247.5262, 28.611, 11.5698,
221.0022, 37.3725, 11.918, 133.5804, 25.6759), AREA.ZONE = c(13683,
3666, 1558, 64830, 41103, 22581, 123819, 90107, 57627, 264735,
223963, 174456), Buffer.zone = c(0L, 0L, 0L, 5L, 5L, 5L, 10L,
10L, 10L, 20L, 20L, 20L)), .Names = c(URBAN_AREA, NO_PATCHES,
MEAN_AREA, AREA.ZONE, Buffer.zone), class = data.frame,
row.names = c(2L,
15L, 19L, 22L, 36L, 40L, 42L, 56L, 60L, 62L, 76L, 80L))
library(lattice)
Region = factor(patches$URBAN_AREA)
lpatches = patches
for (i in 2:4) lpatches[,i] = log10(patches[,i])
# zone not transformed
lpatches.pca = princomp(lpatches[,-c(1,17)],cor=FALSE)
x = lpatches.pca$scores[,1]
y = lpatches.pca$scores[,2]
zz = as.character(patches$Buffer.zone/5)
table(zz)
plsy - trellis.par.get(plot.symbol)
# only 0 or 1 used as plotting symbol
plsy$pch = as.character(rep(1:6,2))
trellis.par.set(plot.symbol,plsy)
xyplot(y ~ x |Region)
# only 1,2,3,4 used as plotting symbol
I actually wish 0,1,2, or 4 to be used - to indicate the zone coded in zz.
Cheers, Murray
PS The xyplots produced on R 2.7.0 for Mac OS X, the first one also on
an older Windows version.
Rolf Turner wrote:
Murray:
I'm not at all sure that I understand what you're driving at --- but
does this do something like what you want?
require(lattice)
set.seed(260808)
n = 50
x = rnorm(n)
y = rnorm(n)
z = ceiling(runif(n,0,4))
g = runif(n,0,6)
G = factor(ceiling(g))
print(xyplot(y ~ x | G,pch=G,
panel=function(x,y,...,subscripts,pch) {
panel.xyplot(x,y,pch=pch[subscripts])
}
))
cheers,
Rolf
##
Attention: This e-mail message is privileged and confidential. If you
are not the intended recipient please delete the message and notify
the sender. Any views or opinions presented are solely those of the
author.
This e-mail has been scanned and cleared by MailMarshal
www.marshalsoftware.com
##
--
Dr Murray Jorgensen http://www.stats.waikato.ac.nz/Staff/maj.html
Department of Statistics, University of Waikato, Hamilton, New Zealand
Email: [EMAIL PROTECTED][EMAIL PROTECTED] Fax 7 838 4155
Phone +64 7 838 4773 wkHome +64 7 825 0441Mobile 021 139 5862
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Displaying Equations in Documentation
On Tue, 26 Aug 2008, Rolf Turner wrote:
On 26/08/2008, at 4:50 AM, Jarrett Byrnes wrote:
I'm currently working on writing up some documentation for some of my code,
but am having the darndest time coding in equations. For example, the
equation in the following:
\details{ Calculated the R Squared for observed endogenous variables in a
structural equation model, as well as several other useful summary
statistics about the error in thoe variables.
R Squared values are calculated as
\deqn{R^{2} = 1-\frac{estimated variance}{observed variance}}
Standardized error coefficients are then calculated as sqrt(1 - R^2).
}
While it shows normally using R CMD Rd2dvi, when I actually compile and
load the package, displays as follows:
R^{2} = 1-frac{estimated variance}{observed variance}
Displays, presumably, in a plain text or html environment.
As my late Uncle Stanley would have said, ``What the hell do you
expect?''
Plain text and html, unlike LaTeX, do not have advanced mathematical
display
capabilities.
That's not the whole story for HTML, as there is MathML. At the time the
decisions were made for Rd rendering, MathML was pretty much unsupported.
That has changed slowly, and nowadays some browsers (e.g. Firefox, Opera)
do support Presentation MathML -- unfortunately they often are short of
suitable fonts. AFAIK, Safari and IE (including Compiled HTML) do not yet
support it.
So someone interested could rewrite the Rdconv code to make use of MathML,
but for now the subset of R users who could benefit from it would be
small.
The R documentation facility --- RTFM!!! (section 2.6 Mathematics)
--- provides a way
to allow for both possibilities. You can do:
\deqn{R^2 = 1 - \frac{estimated variance}{observed
variance}}{R-squared = 1 - (estimated variance)/(observed variance)}
Then in the pdf manual for your package you'll get the sexy
mathematical display, but when you call
for help on line and get a plain text or html version you'll see
R-squared = 1 - (estimated variance)/(observed variance)
which is the best that can be done under the circumstances
I have also tried
\deqn{R^{2} = 1-\frac{{estimated variance}{observed variance}}}
and
\deqn{R^{2} = \frac{1-{estimated variance}{observed variance}}}
with the same result - Rd2dvi is happy, but the display is still wonky in
practice. I've also tried subbing in \eqn{R^{2}} in the rest of the text
in a few places, but, again, it shows as R^{2}. Is there something I'm
missing about inserting equations into R documentation?
Yes.
cheers,
Rolf Turner
##
Attention:\ This e-mail message is privileged and confid...{{dropped:9}}
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
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@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] problem using cat within Sweave code chunks when fig=TRUE
On Mon, 25 Aug 2008, Wolfgang Huber wrote: Hi Paul, You could label=codechunk1a,fig=TRUE,include=FALSE= plot(1:10) @ label=codechunk1b= cat(figure1.eps caption goes here \n,file=readme.txt,append=T) cat(figure1.pdf caption goes here \n,file=readme.txt,append=T) @ You can suppress one of the three executions by setting eps=FALSE in the code chunk options or in the document-wide \SweaveOpts settings. I would also be interested what the third execution is good for. Sweave runs the code once (with the default graphics device), then once more for each 'fig' output. Also, if you set both eps=FALSE and pdf=FALSE, the code chunk is not shown in the (.tex) output, not graphic output file is produced, but the chunk still seems to be executed once. I think you will find graphics output *is* produced, e.g. on Rplots.pdf in non-interactive use. Look for 'eval' in makeRweaveLatexCodeRunner() to see why. best wishes Wolfgang [EMAIL PROTECTED] wrote: Hello R list I was intending to use a cat statement within Sweave code chunks that generate greaphs to generate a readme.txt file listing all the figures generated with a brief caption along the lines of: desired format of readme.txt _ figure1.eps caption for figure1 figure1.pdf caption for figure1 figure2.eps caption for figure2 figure2.pdf caption for figure2 _ As an example, this block of sweave code replicates what I would like to do in principle: label=codechunk1,fig=TRUE,include=FALSE= plot(1:10) cat(figure1.eps caption goes here \n,file=readme.txt,append=T) cat(figure1.pdf caption goes here \n,file=readme.txt,append=T) @ label=codechunk2,fig=TRUE,include=FALSE= plot(11:20) cat(figure2.eps caption goes here \n,file=readme.txt,append=T) cat(figure2.pdf caption goes here \n,file=readme.txt,append=T) @ which I originally though would produce the desired result, however, the cat statement appears to get executed three times producing: readme.txt--- figure1.eps caption goes here figure1.pdf caption goes here figure1.eps caption goes here figure1.pdf caption goes here figure1.eps caption goes here figure1.pdf caption goes here figure2.eps caption goes here figure2.pdf caption goes here figure2.eps caption goes here figure2.pdf caption goes here figure2.eps caption goes here figure2.pdf caption goes here I have figured out that any sweave code with fig=TRUE appears to be executed multiple times (up to three), presumably to write to both eps and pdf graphics devices (not sure what the first/last execution does though...). Does anyone have any ideas about how to only execute the cat statements the first time around so that the output looks like what I specified at the top of this message? Paul Paul Rustomji Rivers and Estuaries CSIRO Land and Water GPO Box 1666 Canberra ACT 2601 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- 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@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] parse and eval character vector
Dear R-help, I have a character vector, some elements will be numeric, some not, and some even empty. E.g.: temp1 - c(abcd, 2 ,) I'm only interested in the numeric elements, the rest I can just throw away. It is easy enough to loop through the vector: temp - try(eval(parse(text=temp1[1])), silent=TRUE); class(temp) # try-error temp - try(eval(parse(text=temp1[2])), silent=TRUE); class(temp) # numeric temp - try(eval(parse(text=temp1[3])), silent=TRUE); class(temp) # NULL and then throw away the non-numeric/NULL stuff. But, as this vector will be long, I would really like to speed things up by not using a loop, and I thought that lapply might do the trick. However: temp.fn - function(x) try(eval(parse(text=x)), silent=TRUE) temp2 - lapply(temp1, FUN=temp.fn) class(temp2[2]) # list, for elements 1, 2, and 3 and I don't know how to extract the numeric elements from here. So, can I either use lapply as above and somehow get the information I need out of temp2 (I've tried using unlist but had no success), or is there some other function that I can apply to my character vector to avoid looping? Rob. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] parse and eval character vector
Just use as.numeric. Non numeric will be NA. So the solution of your problem is na.omit(as.numeric(temp1)) HTH, Thierry ir. Thierry Onkelinx Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest Cel biometrie, methodologie en kwaliteitszorg / Section biometrics, methodology and quality assurance Gaverstraat 4 9500 Geraardsbergen Belgium tel. + 32 54/436 185 [EMAIL PROTECTED] www.inbo.be To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher The plural of anecdote is not data. ~ Roger Brinner The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. ~ John Tukey -Oorspronkelijk bericht- Van: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Namens Rob Foxall Verzonden: dinsdag 26 augustus 2008 10:36 Aan: r-help@r-project.org Onderwerp: [R] parse and eval character vector Dear R-help, I have a character vector, some elements will be numeric, some not, and some even empty. E.g.: temp1 - c(abcd, 2 ,) I'm only interested in the numeric elements, the rest I can just throw away. It is easy enough to loop through the vector: temp - try(eval(parse(text=temp1[1])), silent=TRUE); class(temp) # try-error temp - try(eval(parse(text=temp1[2])), silent=TRUE); class(temp) # numeric temp - try(eval(parse(text=temp1[3])), silent=TRUE); class(temp) # NULL and then throw away the non-numeric/NULL stuff. But, as this vector will be long, I would really like to speed things up by not using a loop, and I thought that lapply might do the trick. However: temp.fn - function(x) try(eval(parse(text=x)), silent=TRUE) temp2 - lapply(temp1, FUN=temp.fn) class(temp2[2]) # list, for elements 1, 2, and 3 and I don't know how to extract the numeric elements from here. So, can I either use lapply as above and somehow get the information I need out of temp2 (I've tried using unlist but had no success), or is there some other function that I can apply to my character vector to avoid looping? Rob. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is door een geldig ondertekend document. The views expressed in this message and any annex are purely those of the writer and may not be regarded as stating an official position of INBO, as long as the message is not confirmed by a duly signed document __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] parse and eval character vector
That's great, thanks. I can live with the warnings! Cheers, Rob On Tue, Aug 26, 2008 at 4:49 PM, ONKELINX, Thierry [EMAIL PROTECTED] wrote: Just use as.numeric. Non numeric will be NA. So the solution of your problem is na.omit(as.numeric(temp1)) HTH, Thierry ir. Thierry Onkelinx Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest Cel biometrie, methodologie en kwaliteitszorg / Section biometrics, methodology and quality assurance Gaverstraat 4 9500 Geraardsbergen Belgium tel. + 32 54/436 185 [EMAIL PROTECTED] www.inbo.be To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher The plural of anecdote is not data. ~ Roger Brinner The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. ~ John Tukey -Oorspronkelijk bericht- Van: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Namens Rob Foxall Verzonden: dinsdag 26 augustus 2008 10:36 Aan: r-help@r-project.org Onderwerp: [R] parse and eval character vector Dear R-help, I have a character vector, some elements will be numeric, some not, and some even empty. E.g.: temp1 - c(abcd, 2 ,) I'm only interested in the numeric elements, the rest I can just throw away. It is easy enough to loop through the vector: temp - try(eval(parse(text=temp1[1])), silent=TRUE); class(temp) # try-error temp - try(eval(parse(text=temp1[2])), silent=TRUE); class(temp) # numeric temp - try(eval(parse(text=temp1[3])), silent=TRUE); class(temp) # NULL and then throw away the non-numeric/NULL stuff. But, as this vector will be long, I would really like to speed things up by not using a loop, and I thought that lapply might do the trick. However: temp.fn - function(x) try(eval(parse(text=x)), silent=TRUE) temp2 - lapply(temp1, FUN=temp.fn) class(temp2[2]) # list, for elements 1, 2, and 3 and I don't know how to extract the numeric elements from here. So, can I either use lapply as above and somehow get the information I need out of temp2 (I've tried using unlist but had no success), or is there some other function that I can apply to my character vector to avoid looping? Rob. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is door een geldig ondertekend document. The views expressed in this message and any annex are purely those of the writer and may not be regarded as stating an official position of INBO, as long as the message is not confirmed by a duly signed document __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] environments
Hi there,
I try to understand the usage of environments but I'm not sure if I get it. I
wrote a test script like this:
testenv - new.env(environment())
myfun - function(x) {
print(testvar)
testenv$testvar_2 - 20
}
environment(myfun) - testenv
testenv$testvar - 10
myfun(hello)
ls(envir = testenv)
Now, I was wondering if there is any way to create new variables in my
environment without this testenv$ I know that I can access it that way if
I do an attach(testenv) before, but that does not help when creating new ones...
Do I completely misunderstand the concept?
I'm just looking for an elegant way to access objects of a graphical
userinterface from each handler-function and so on. And I thought it might be
good to pack them into an environment...
Antje
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] savePlot() does not save plot with format set
R-help, Whenever I try to save a plot with savePlot the file is not stored in my hard disk with the selected format. Several formats are set and none of them works. I just get the file name with missing extension and it can't be open with programs like Paint and Microsoft Photo Editor Th only one able to open it is Windows Picture and Fax Viewer plot(rnorm(10)) savePlot(test, type=png) savePlot(test, type=bmp) My platform is Windows XP SP3 version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major 2 minor 7.0 year 2008 month 04 day22 svn rev45424 language R version.string R version 2.7.0 (2008-04-22) Thanks in advanced __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Variance-covariance matrix
Dear R help forum, I am using the function 'coxph' to obtain hazard ratios for the comparison of a standard treatment to new treatments. This is easily obtained by fitting the relevant model and then calling exp(coef(fit1)) say. I now want to obtain the hazard ratio for the comparison of two non-standard treatments. From a statistical point of view, this can be achieved by dividing the exponentiated coefficients of 2 comparisions. E.g. to compared new treatment 1 (nt1) to new treatment 2 (nt2) we can fit 2 models: fit1 = standard treatment vs nt1 fit2 = standard treatment vs nt2. The required hazard ratio is therefore exp(coef(fit1))/exp(coef(fit2)) In order to obtain an associated confidence interval for this I require the covariance of this comparison. I know that R gives the variance-covariance matrix by the command 'fit$var'. However, this only gives the covariance matrix for non standard drugs and not the full covariance matrix. Can anyone tell me how to obtain the full covariance matrix? Thank you, Laura [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Output to Latex using Memisc almost works
Thank you Martin,
It turned out that my problem was that I used underscores _ in my variable
names. No problem for R, but Latex didn't like it.
Kind regards,
Eelke
Martin Elff wrote:
On Monday 25 August 2008 (15:36:50), Eelke wrote:
Hello,
I'm using memisc to output regression results to tables and latex. My
problem is that the output that Latex needs must be in between $ $ so
that
it is read as formula but memisc does not output the result between $ $.
For example, latex needs: $0.05^{***}$ and memisc outputs 0.05^{***} in
an
entry.
I am new to Latex and I imagine it could also be a latex 'problem' and
not
necessarily memisc.
I would suggest you put
\usepackage{dcolumn} into the preamble of your LaTeX file. This will
automagically insert the dollars into the columns.
memisc's mtable function works fine with that (at least in my
experience, I use it on a regular basis as one might suspect :-)
Cheers,
Martin
--
-
Dr. Martin Elff
Department of Social Sciences
University of Mannheim
A5, Room 328
68131 Mannheim
Germany
Phone: ++49-621-181-2093
Fax: ++49-621-181-2099
E-Mail: [EMAIL PROTECTED]
Homepage:
http://webrum.uni-mannheim.de/sowi/elff
http://www.sowi.uni-mannheim.de/lspwivs/
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
View this message in context:
http://www.nabble.com/Output-to-Latex-using-Memisc-almost-works-tp19144010p19158769.html
Sent from the R help mailing list archive at Nabble.com.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Variance-covariance matrix
Laura Bonnett wrote: Dear R help forum, I am using the function 'coxph' to obtain hazard ratios for the comparison of a standard treatment to new treatments. This is easily obtained by fitting the relevant model and then calling exp(coef(fit1)) say. I now want to obtain the hazard ratio for the comparison of two non-standard treatments. From a statistical point of view, this can be achieved by dividing the exponentiated coefficients of 2 comparisions. E.g. to compared new treatment 1 (nt1) to new treatment 2 (nt2) we can fit 2 models: fit1 = standard treatment vs nt1 fit2 = standard treatment vs nt2. The required hazard ratio is therefore exp(coef(fit1))/exp(coef(fit2)) In order to obtain an associated confidence interval for this I require the covariance of this comparison. I know that R gives the variance-covariance matrix by the command 'fit$var'. However, this only gives the covariance matrix for non standard drugs and not the full covariance matrix. Can anyone tell me how to obtain the full covariance matrix? What kind of data do you have? Is the standard treatment group the same in both comparisons? If so, why not just have a three-level treatment factor and compare nt1 to nt2 directly. If the control groups are completely separate, then the covariance between fits made on independent data is of course zero. Thank you, Laura [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] options(contrasts)
Code: options(contrasts) $contrasts factor ordered contr.treatment contr.poly I want to change the first entry ONLY, without retyping contr.poly. How do I do it? I have tried various possibilities and cannot get anything to work. I found out that the response to options(contrasts) has class list, but that doesn't help me, although I think it ought to help. Second question (metaquestion). How should I go about finding out the answer to a question like How does one change a single item in a list? My answer to the meta-meta-question is to post to this list. I hope that at least that part is correct. Thanks for any help. David Epstein -- View this message in context: http://www.nabble.com/options%28%22contrasts%22%29-tp19158786p19158786.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Variance-covariance matrix
The standard treatment is the same in both comparison. How do you do a three-level treatment factor? I thought you had to have a censoring indicator which took values 0 or 1 not 1, 2 or 3? Thanks, Laura On Tue, Aug 26, 2008 at 11:05 AM, Peter Dalgaard [EMAIL PROTECTED]wrote: Laura Bonnett wrote: Dear R help forum, I am using the function 'coxph' to obtain hazard ratios for the comparison of a standard treatment to new treatments. This is easily obtained by fitting the relevant model and then calling exp(coef(fit1)) say. I now want to obtain the hazard ratio for the comparison of two non-standard treatments. From a statistical point of view, this can be achieved by dividing the exponentiated coefficients of 2 comparisions. E.g. to compared new treatment 1 (nt1) to new treatment 2 (nt2) we can fit 2 models: fit1 = standard treatment vs nt1 fit2 = standard treatment vs nt2. The required hazard ratio is therefore exp(coef(fit1))/exp(coef(fit2)) In order to obtain an associated confidence interval for this I require the covariance of this comparison. I know that R gives the variance-covariance matrix by the command 'fit$var'. However, this only gives the covariance matrix for non standard drugs and not the full covariance matrix. Can anyone tell me how to obtain the full covariance matrix? What kind of data do you have? Is the standard treatment group the same in both comparisons? If so, why not just have a three-level treatment factor and compare nt1 to nt2 directly. If the control groups are completely separate, then the covariance between fits made on independent data is of course zero. Thank you, Laura [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] loop with splitted data
Hi to all,
seems to be simple, but I do not find the solution:
What must I write for the splitted to get
splitted$3$x and splitted$3$x
y = c(rep(2,5),rep(3,5))
ma - data.frame(x = 1:10, y=y )
splitted - split(ma, ma$y)
for (counter in (min(ma$y):max(ma$y)))
{
splitted$x
}
Regards Knut
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] paste: multiple collapse arguments?
Try this also: paste(apply(matrix(x, 2), 2, paste, collapse = ' '), collapse = \n) On Mon, Aug 25, 2008 at 8:36 PM, remko duursma [EMAIL PROTECTED] wrote: Dear R-helpers, I have a numeric vector, like: x - c(1,2,3,4,5,6) I make this into a string for output to a text file, separated by \n: paste(x, collapse=\n) Is there a way to alternate the collapse argument? So between the first two elements of x, I want to separate by , then by \n, and so forth. The result should then look like: 1 2\n3 4\n5 6 (This way I get 2 elements of x on each line using writeLines, instead of one or all). I could do this in some ugly loop, but surely there is a better way? thanks, Remko _ Get thousands of games on your PC, your mobile phone, and the web with Windows(R). [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Henrique Dallazuanna Curitiba-Paraná-Brasil 25° 25' 40 S 49° 16' 22 O __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] environments
I think you need assign, see ?assign for more details.
On Tue, Aug 26, 2008 at 6:02 AM, Antje [EMAIL PROTECTED] wrote:
Hi there,
I try to understand the usage of environments but I'm not sure if I get it.
I wrote a test script like this:
testenv - new.env(environment())
myfun - function(x) {
print(testvar)
testenv$testvar_2 - 20
}
environment(myfun) - testenv
testenv$testvar - 10
myfun(hello)
ls(envir = testenv)
Now, I was wondering if there is any way to create new variables in my
environment without this testenv$ I know that I can access it that way
if I do an attach(testenv) before, but that does not help when creating new
ones...
Do I completely misunderstand the concept?
I'm just looking for an elegant way to access objects of a graphical
userinterface from each handler-function and so on. And I thought it might
be good to pack them into an environment...
Antje
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Henrique Dallazuanna
Curitiba-Paraná-Brasil
25° 25' 40 S 49° 16' 22 O
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] savePlot() does not save plot with format set
You need type the extension of the file: plot(rnorm(10)) savePlot(test.png, type=png) savePlot(test.bmp, type=bmp) On Tue, Aug 26, 2008 at 6:29 AM, Luis Ridao Cruz [EMAIL PROTECTED] wrote: R-help, Whenever I try to save a plot with savePlot the file is not stored in my hard disk with the selected format. Several formats are set and none of them works. I just get the file name with missing extension and it can't be open with programs like Paint and Microsoft Photo Editor Th only one able to open it is Windows Picture and Fax Viewer plot(rnorm(10)) savePlot(test, type=png) savePlot(test, type=bmp) My platform is Windows XP SP3 version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major 2 minor 7.0 year 2008 month 04 day22 svn rev45424 language R version.string R version 2.7.0 (2008-04-22) Thanks in advanced __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Henrique Dallazuanna Curitiba-Paraná-Brasil 25° 25' 40 S 49° 16' 22 O __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] options(contrasts)
On 26/08/2008 6:30 AM, David Epstein wrote: Code: options(contrasts) $contrasts factor ordered contr.treatment contr.poly I want to change the first entry ONLY, without retyping contr.poly. How do I do it? I have tried various possibilities and cannot get anything to work. This doesn't really save typing, but you could do this: conts - options(contrasts)$contrasts conts[1] - contr.poly options(contrasts = conts) (One thing I wonder about: which version of R are you using? Mine shows the names of the contrasts as unordered and ordered.) I found out that the response to options(contrasts) has class list, but that doesn't help me, although I think it ought to help. Second question (metaquestion). How should I go about finding out the answer to a question like How does one change a single item in a list? Read the Introduction to R manual, in particular the bits on replacement functions. The general idea is that you can sometimes index the target of an assignment, as I did above. I could also have written it as conts - options(contrasts) conts$contrasts[1] - contr.poly options(conts) One might guess that options(contrasts)[1] - contr.poly would work, but one would be wrong. There is no options- replacement function. (There might be one in some contributed package; I was just looking in the base packages.) Duncan Murdoch My answer to the meta-meta-question is to post to this list. I hope that at least that part is correct. Thanks for any help. David Epstein __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to do a meta-analysis plot
At 19:44 25/08/2008, Jorge Ivan Velez wrote:
Dear R-list,
I'd like to do a meta-analysis plot similar to
Since these plots are known as forest plots
?forestplot
might help you.
As it says you have to do a bit more work but you do get much more flexibility
install.packages('rmeta')
require(rmeta)
data(catheter)
a - meta.MH(n.trt, n.ctrl, col.trt, col.ctrl, data=catheter,
names=Name, subset=c(13,6,5,3,7,12,4,11,1,8,10,2))
summary(a)
plot(a)
(see attached file) by using my own OR (Odds Ratio) and 95% Confidence
Interval data set, which looks like
mydata=data.frame(OR=c(2.04545454545, 1.10434782609, 1.22588104401,
1.14102564103,
1.20579245527, 1.375, 1.16535433071),
L95=c(1.22839621997, 0.858106819302, 1.0964802088, 0.841934120955,
0.969786886818, 1.01498537023, 0.919391492382),
U95=c(3.40546755139, 1.42122051928, 1.37055308613, 1.54632513827,
1.49917372998, 1.86258857302, 1.47707220868)
)
rownames(mydata)=c(paste(Study,1:6,sep=),'Summary')
mydata
My problem is that I don't have the raw data as rmeta _requires_ and, even
when I have my data set in the _same_ (?) format that summary(a), when I
tried plot(mydata) it doesn't work. Another approach I used was to change
the class of my object but it didn't work either. I'm running XP SP2 on a
2.4 GHz Intel-Core 2 Duo processor and my R-session info is the following:
R version 2.7.2 RC (2008-08-18 r46388)
i386-pc-mingw32
locale:
LC_COLLATE=English_United States.1252;LC_CTYPE=English_United
States.1252;LC_MONETARY=English_United
States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252
attached base packages:
[1] grid stats graphics grDevices utils datasets methods
base
other attached packages:
[1] rmeta_2.14 RODBC_1.2-3
loaded via a namespace (and not attached):
[1] tools_2.7.2
I would greatly appreciate any ideas about how should I proceed.
Thanks in advance,
Jorge
Content-Type: application/pdf; name=example MH.pdf
X-Attachment-Id: f_fkbfhzhb0
Content-Disposition: attachment; filename=example MH.pdf
Michael Dewey
http://www.aghmed.fsnet.co.uk
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
Hi Mario, If I understand your problem statement, the random choice of one of the two envelopes ensures that the probability of choosing one or the other is 0.5. If you find 10 (units, I can't find the pound symbol on this keyboard), there is an equal likelihood that the other envelope contains 5 or 20 units. Your expected value for the swap is thus 7.5*0.5 + 15*0.5 = 11.25 units. The trick is in the twice as much, half as much statement. If it was x units more or less, a symmetric distribution, swapping wouldn't make any difference. Jim __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] String search: Return closest match
Hi, I have to match names where names can be recorded with errors or additions. Now I am searching for a string search function which returns always the closest match. E.g. searching for Washington it should return only Washington but not Washington, D.C. But it also could be that the list contains only Hamburg but the record I am searching for is Hamburg OTM and then we still want to find Hamburg. Or maybe the list contains Hamburg and Hamberg but we are searching for Hamburg and thus only this should this one should be returned. agrep() returns all close matches but unfortunately does not return the degree of closeness otherwise selection would be easy. Is there such a function already implemented? Thanks a million for your help, Werner __ verfügt über einen herausragenden Schutz gegen Massenmails. http://mail.yahoo.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
Hi again, Oops, I meant the expected value of the swap is: 5*0.5 + 20*0.5 = 12.5 Too late, must get to bed. Jim __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to do a meta-analysis plot
Another way to do it (based on summary data): http://finzi.psych.upenn.edu/R/Rhelp02a/archive/131342.html -- Jonathan Baron, Professor of Psychology, University of Pennsylvania Home page: http://www.sas.upenn.edu/~baron __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
On 26/08/2008 7:54 AM, Jim Lemon wrote: Hi again, Oops, I meant the expected value of the swap is: 5*0.5 + 20*0.5 = 12.5 Too late, must get to bed. But that is still wrong. You want a conditional expectation, conditional on the observed value (10 in this case). The answer depends on the distribution of the amount X, where the envelopes contain X and 2X. For example, if you knew that X was at most 5, you would know you had just observed 2X, and switching would be a bad idea. The paradox arises because people want to put a nonsensical Unif(0, infinity) distribution on X. The Wikipedia article points out that it can also arise in cases where the distribution on X has infinite mean: a mathematically valid but still nonsensical possibility. Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] savePlot() does not save plot with format set
Strictly, you need to type the extension _if your filename or path
include a period (.)_ but not otherwise.
The filename test alone will be saved as paste(test,type). So will,
for example, /plots/test. But filenames such as
test.it or ../plots/test will not include the extension
automatically.
Steve Ellison
Henrique Dallazuanna [EMAIL PROTECTED] 26/08/2008 12:12
You need type the extension of the file:
plot(rnorm(10))
savePlot(test.png, type=png)
savePlot(test.bmp, type=bmp)
On Tue, Aug 26, 2008 at 6:29 AM, Luis Ridao Cruz [EMAIL PROTECTED] wrote:
R-help,
Whenever I try to save a plot with savePlot
the file is not stored in my hard disk with the selected
format.
***
This email and any attachments are confidential. Any use...{{dropped:8}}
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Using interactive plots to get information about data points
jcarmichael jcarmichael314 at gmail.com writes: Thank you for the GGobi reference, it is a very handy tool! My main goal, however, is to be able to identify univariate outliers (with boxplots for example), and I'm having a hard time getting the rggobi package to do that. Hmmm... I guess because GGobi is designed for the visualization of high-dimensional data, boxplots of individual vectors isn't a priority. Sorry if that was a bad lead. (Still cool software, though!) You can, however, use identify. An RSiteSearch for identify boxplot will give you a number of leads. I often recommend John Fox's R Commander GUI as a tool to help with learning syntax. In this case, the Rcmdr dialogue for boxplots includes a checkbox for identify outliers with mouse. You can use this dialogue with this option checked, and then examine the syntax to see how it was done: data(Angell, package=car) boxplot(Angell$hetero, ylab=hetero) identify(rep(1, length(Angell$hetero)), Angell$hetero, rownames(Angell)) The results of the RSiteSearch will explain how/why this works. (Hint: the rep(1,length(z)),z pattern essentially defines x,y coordinates of all the points that make up the boxplot - for a univariate boxplot, all have an 'x' coordinate of 1). In a similar vein, the latticist GUI (package playWith - need GTK libraries or runtime (windows) installed) has similar functionality with parallel boxplots. Hope this helps. Michael Bibo Queensland Health __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] SQL Primer for R
Sorry, chaps. I need one more: dbDisconnect(con.in) Error in sqliteCloseConnection(conn, ...) : RS-DBI driver: (close the pending result sets before closing this connection) I am pretty sure I have fetched everything there is to be fetched. I am not sure what I need to do to say goodbye (or to find out what is still pending). ?dbDisconnect doesn't tell me. PS: the documentation for dbConnect should probably add dbDisconnect to its 'See also' section. regards, /iaw Really irrelevant PS: the by function could keep the number of observations that go into each category. I know it can be computed separately, which is what I am doing now. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] environments
On Tue, Aug 26, 2008 at 6:07 AM, Henrique Dallazuanna [EMAIL PROTECTED] wrote:
I think you need assign, see ?assign for more details.
On Tue, Aug 26, 2008 at 6:02 AM, Antje [EMAIL PROTECTED] wrote:
Hi there,
I try to understand the usage of environments but I'm not sure if I get it.
I wrote a test script like this:
testenv - new.env(environment())
myfun - function(x) {
print(testvar)
testenv$testvar_2 - 20
}
environment(myfun) - testenv
testenv$testvar - 10
As Henrique said, the canonical way of assigning a value within an
environment is the assign. A more obscure, but also more effective,
approach is evalq which quotes an expression then evaluates it in the
given environment. For example
env - new.env()
evalq({aa - 1:3; bb - LETTERS[1:9]; cc - list(A = aa, B = bb)}, env)
objects(env)
[1] aa bb cc
env$aa
[1] 1 2 3
myfun(hello)
ls(envir = testenv)
Now, I was wondering if there is any way to create new variables in my
environment without this testenv$ I know that I can access it that way
if I do an attach(testenv) before, but that does not help when creating new
ones...
Do I completely misunderstand the concept?
I'm just looking for an elegant way to access objects of a graphical
userinterface from each handler-function and so on. And I thought it might
be good to pack them into an environment...
Antje
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Henrique Dallazuanna
Curitiba-Paraná-Brasil
25° 25' 40 S 49° 16' 22 O
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] String search: Return closest match
I have to match names where names can be recorded with errors or
additions.
Now I am searching for a string search function which returns always
the closest match. E.g. searching for Washington it should
return only Washington but not Washington, D.C. But it also could be
that the list contains only Hamburg but the record I am searching
for is Hamburg OTM and then we still want to find Hamburg. Or
maybe the list contains Hamburg and Hamberg but we are searching
for Hamburg and thus only this should this one should be returned.
agrep() returns all close matches but unfortunately does not
return the degree of closeness otherwise selection would be easy.
Is there such a function already implemented?
The Levenshtein distance is a common metric for determining how close two
string are (in fact, agrep uses this). There's a function to calculate it
on the R wiki.
http://wiki.r-project.org/rwiki/doku.php?id=tips:data-strings:levenshtein
You can use this to find the closest string. (If your set of cities is
large, it may be quickest to use agrep to narrow the selection first,
since the pure R implementation of levenshtein is likely to be slow.)
Regards,
Richie.
Mathematical Sciences Unit
HSL
ATTENTION:
This message contains privileged and confidential inform...{{dropped:20}}
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] environments
Douglas Bates wrote:
As Henrique said, the canonical way of assigning a value within an
environment is the assign. A more obscure, but also more effective,
approach is evalq which quotes an expression then evaluates it in the
given environment. For example
env - new.env()
evalq({aa - 1:3; bb - LETTERS[1:9]; cc - list(A = aa, B = bb)}, env)
objects(env)
[1] aa bb cc
env$aa
[1] 1 2 3
Yes, and the with() construct works similarly. You do have to be careful
to note that also the right hand side of the assignment is evaluated in
env. This can have unexpected consequences.
Notice also the possibility of lexical scoping and superassignment
-. See for instance
file.show(system.file(demo/scoping.R, package=base))
--
O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B
c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
(*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918
~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] environments
Okay, I see, there is no really easy way (I was wondering whether I can set
an environment as default for new created variables). Is there any difference
if I call
myenv$myvar - 10
or
assign(myvar,10, env=myenv)
?
Antje
Douglas Bates schrieb:
On Tue, Aug 26, 2008 at 6:07 AM, Henrique Dallazuanna [EMAIL PROTECTED] wrote:
I think you need assign, see ?assign for more details.
On Tue, Aug 26, 2008 at 6:02 AM, Antje [EMAIL PROTECTED] wrote:
Hi there,
I try to understand the usage of environments but I'm not sure if I get it.
I wrote a test script like this:
testenv - new.env(environment())
myfun - function(x) {
print(testvar)
testenv$testvar_2 - 20
}
environment(myfun) - testenv
testenv$testvar - 10
As Henrique said, the canonical way of assigning a value within an
environment is the assign. A more obscure, but also more effective,
approach is evalq which quotes an expression then evaluates it in the
given environment. For example
env - new.env()
evalq({aa - 1:3; bb - LETTERS[1:9]; cc - list(A = aa, B = bb)}, env)
objects(env)
[1] aa bb cc
env$aa
[1] 1 2 3
myfun(hello)
ls(envir = testenv)
Now, I was wondering if there is any way to create new variables in my
environment without this testenv$ I know that I can access it that way
if I do an attach(testenv) before, but that does not help when creating new
ones...
Do I completely misunderstand the concept?
I'm just looking for an elegant way to access objects of a graphical
userinterface from each handler-function and so on. And I thought it might
be good to pack them into an environment...
Antje
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Henrique Dallazuanna
Curitiba-Paraná-Brasil
25° 25' 40 S 49° 16' 22 O
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] svymeans question
I have the following code which produces the output below it clus1 - svydesign(ids = ~schid, data = lower_dat) items - as.formula(paste( ~ , paste(lset, collapse= +))) rr1 - svymean(items, clus1, deff='replace', na.rm=TRUE) rr1 mean SE DEff W525209 0.719748 0.015606 2.4932 W525223 0.508228 0.027570 6.2802 W525035 0.827202 0.014060 2.8561 W525131 0.805421 0.015425 3.1350 W525033 0.242982 0.020074 4.5239 W525163 0.904647 0.013905 4.6289 W525165 0.439981 0.020029 3.3620 W525167 0.148112 0.013047 2.7860 W525177 0.865924 0.014977 3.9898 W525179 0.409003 0.020956 3.7515 W525181 0.634076 0.022076 4.3372 W525183 0.242498 0.019073 4.0894 W525401 0.262343 0.021830 3.4354 W525059 0.854792 0.016551 4.5576 W525251 0.691191 0.025010 6.0512 W525083 0.433204 0.017310 2.5200 W525289 0.634560 0.012762 1.4504 W524763 0.791868 0.014478 2.6265 W524765 0.223621 0.019627 4.5818 W524951 0.242982 0.016796 3.1669 W524769 0.820910 0.016786 3.9579 W524771 0.872701 0.015853 4.6712 W524839 0.518877 0.026433 5.7794 W525374 1.209584 0.043065 5.1572 W524885 0.585673 0.027780 6.5674 W525377 1.100678 0.050093 5.8851 W524787 0.839303 0.012994 2.5852 W524789 0.339787 0.019230 3.4041 W524791 0.847047 0.012885 2.6461 W524825 0.500968 0.021988 3.9935 W524795 0.868345 0.014951 4.0377 W524895 0.864472 0.013872 3.3917 W524897 0.804937 0.020070 5.2977 W524967 0.475799 0.032137 8.5511 W525009 0.681994 0.018670 3.3188 However, when I do the following: svymean(~W524787, clus1, deff='replace', na.rm=TRUE) mean SE DEff W524787 0.855547 0.011365 4.1158 Compare this to the value in the row 9 up from the bottom to see it is different. Computing the mean of the item by itself with svymeans agrees with the sample mean mean(lower_dat$W524787, na.rm=T) [1] 0.8555471 Now, I know that there is a covariance between the variables, but I was under the impression that the sample mean was still of pragmatic utility, but to account for sample design only the standard error is affected. In the work I am doing, it is important for the means of the items from svymeans to be the same as the sample mean when it is computed by itself. It's a bit of a story as to why, and I can provide that info if relevant. I don't see an argument in svydesign or in svymean that would allow for me to treat the variables as being independent. But, maybe I am missing something else and would welcome any reactions. Harold __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] environments
Antje wrote: Okay, I see, there is no really easy way (I was wondering whether I can set an environment as default for new created variables). Is there any difference if I call myenv$myvar - 10 or assign(myvar,10, env=myenv) ? No. And with(myenv, myvar - 10) is also the same. However, you do have to be careful with with(myenv, myvar - x) vs. myenv$myvar - x because it might be different x. -- O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] options(contrasts)
On 26-Aug-08 10:30:30, David Epstein wrote: Code: options(contrasts) $contrasts factor ordered contr.treatment contr.poly I want to change the first entry ONLY, without retyping contr.poly. How do I do it? I have tried various possibilities and cannot get anything to work. I found out that the response to options(contrasts) has class list, but that doesn't help me, although I think it ought to help. Second question (metaquestion). How should I go about finding out the answer to a question like How does one change a single item in a list? In view of your meta-meta-strategy, here is a response to the meta-question: If you sijmply want to replace a given component (say $C) of a list L, then use code like: L$C - your.replacement If you want to change the contents of a component of a list, then what you need to do will depend on the nature of that component (number, vector, array, anova table, list ... ). Simple example: L-list(A=A,B=B,C=Z,D=D) L # $A # [1] A # $B # [1] B # $C # [1] Z # $D # [1] D C-L$C ## extract $C from L C # [1] Z C-C ## change it L$C-C ## put it back L # $A # [1] A # $B # [1] B # $C # [1] C # $D # [1] D My answer to the meta-meta-question is to post to this list. I hope that at least that part is correct. It has been known to work ... Ted. E-Mail: (Ted Harding) [EMAIL PROTECTED] Fax-to-email: +44 (0)870 094 0861 Date: 26-Aug-08 Time: 14:42:29 -- XFMail -- __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] processing subset lists and then plot(density())
d - structure(list(Site = structure(c(8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L), .Label = c(119, 148, 179, 185, 190, 198, 202, 215, 61, BC, HC, SC), class = factor), EPT.Taxa = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1L, NA, 3L, 1L, 5L, 7L, 3L, 11L, 3L, 14L, 12L, 12L, 0L, 0L, 4L, 0L, 5L, 3L, 2L, 6L, 1L, 8L, 6L, 9L, 1L, 0L, 2L, 0L, 5L, 2L, 1L, 0L, 2L, 4L, 6L, 8L, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 2L, 4L, 3L, 2L, 7L, 6L, 4L, 8L, 7L, 11L, 8L, 11L, 1L, 3L, 0L, 2L, 7L, 8L, 2L, 7L, 6L, 11L, 6L, 12L, 1L, 1L, 0L, 0L, 0L, 1L, 2L, 9L, 6L, 16L, 6L, 10L, 2L, 3L, 3L, 2L, 5L, 2L, 0L, 3L, 6L, 10L, NA, 10L, 1L, 0L, 3L, 1L, 4L, 2L, 3L, 2L, 3L, 11L, 10L, 6L)), .Names = c(Site, EPT.Taxa), class = data.frame, row.names = c(NA, -144L)) subset(d, Site==215) #I would like to plot(density()) of each of the Sites #I tried list-levels(d$Site) lapply(d, FUN=subset, Site==list) #I would like to be able to make a list of the subsets based on Site and then plot them all on one graphics window # Am I working in the right direction- I have just discovered lists (I know I know, I am a little slow) -- Stephen Sefick Research Scientist Southeastern Natural Sciences Academy Let's not spend our time and resources thinking about things that are so little or so large that all they really do for us is puff us up and make us feel like gods. We are mammals, and have not exhausted the annoying little problems of being mammals. -K. Mullis __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
Duncan: I think I see what you're saying but the strange thing is that if you use the utility function log(x) rather than x, then the expected values are equal. Somehow, if you are correct and I think you are, then taking the log , fixes the distribution of x which is kind of odd to me. I'm sorry to belabor this non R related discussion and I won't say anything more about it but I worked/talked on this with someone for about a month a few years ago and we gave up so it's interesting for me to see this again. Mark -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Duncan Murdoch Sent: Tuesday, August 26, 2008 8:15 AM To: Jim Lemon Cc: r-help@r-project.org; Mario Subject: Re: [R] Two envelopes problem On 26/08/2008 7:54 AM, Jim Lemon wrote: Hi again, Oops, I meant the expected value of the swap is: 5*0.5 + 20*0.5 = 12.5 Too late, must get to bed. But that is still wrong. You want a conditional expectation, conditional on the observed value (10 in this case). The answer depends on the distribution of the amount X, where the envelopes contain X and 2X. For example, if you knew that X was at most 5, you would know you had just observed 2X, and switching would be a bad idea. The paradox arises because people want to put a nonsensical Unif(0, infinity) distribution on X. The Wikipedia article points out that it can also arise in cases where the distribution on X has infinite mean: a mathematically valid but still nonsensical possibility. Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] String search: Return closest match
That works perfectly, great. Thanks a lot for that Richard! Werner - Ursprüngliche Mail Von: [EMAIL PROTECTED] [EMAIL PROTECTED] An: Werner Wernersen [EMAIL PROTECTED] CC: [EMAIL PROTECTED]; [EMAIL PROTECTED] Gesendet: Dienstag, den 26. August 2008, 14:10:11 Uhr Betreff: Re: [R] String search: Return closest match I have to match names where names can be recorded with errors or additions. Now I am searching for a string search function which returns always the closest match. E.g. searching for Washington it should return only Washington but not Washington, D.C. But it also could be that the list contains only Hamburg but the record I am searching for is Hamburg OTM and then we still want to find Hamburg. Or maybe the list contains Hamburg and Hamberg but we are searching for Hamburg and thus only this should this one should be returned. agrep() returns all close matches but unfortunately does not return the degree of closeness otherwise selection would be easy. Is there such a function already implemented? The Levenshtein distance is a common metric for determining how close two string are (in fact, agrep uses this). There's a function to calculate it on the R wiki. http://wiki.r-project.org/rwiki/doku.php?id=tips:data-strings:levenshtein You can use this to find the closest string. (If your set of cities is large, it may be quickest to use agrep to narrow the selection first, since the pure R implementation of levenshtein is likely to be slow.) Regards, Richie. Mathematical Sciences Unit HSL ATTENTION: This message contains privileged and confidential information intended for the addressee(s) only. If this message was sent to you in error, you must not disseminate, copy or take any action in reliance on it and we request that you notify the sender immediately by return email. Opinions expressed in this message and any attachments are not necessarily those held by the Health and Safety Laboratory or any person connected with the organisation, save those by whom the opinions were expressed. Please note that any messages sent or received by the Health and Safety Laboratory email system may be monitored and stored in an information retrieval system. Scanned by MailMarshal - Marshal's comprehensive email content security solution. Download a free evaluation of MailMarshal at www.marshal.com __ gt über einen herausragenden Schutz gegen Massenmails. http://mail.yahoo.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] processing subset lists and then plot(density())
Here's a solution with ggplot2 library(ggplot2) ggplot(na.omit(d), aes(x = EPT.Taxa, colour = Site)) + geom_density() ggplot(na.omit(d), aes(x = EPT.Taxa)) + geom_density() + facet_grid(Site ~ .) HTH, Thierry ir. Thierry Onkelinx Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest Cel biometrie, methodologie en kwaliteitszorg / Section biometrics, methodology and quality assurance Gaverstraat 4 9500 Geraardsbergen Belgium tel. + 32 54/436 185 [EMAIL PROTECTED] www.inbo.be To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher The plural of anecdote is not data. ~ Roger Brinner The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. ~ John Tukey -Oorspronkelijk bericht- Van: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Namens stephen sefick Verzonden: dinsdag 26 augustus 2008 15:51 Aan: R Help Onderwerp: [R] processing subset lists and then plot(density()) d - structure(list(Site = structure(c(8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L), .Label = c(119, 148, 179, 185, 190, 198, 202, 215, 61, BC, HC, SC), class = factor), EPT.Taxa = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1L, NA, 3L, 1L, 5L, 7L, 3L, 11L, 3L, 14L, 12L, 12L, 0L, 0L, 4L, 0L, 5L, 3L, 2L, 6L, 1L, 8L, 6L, 9L, 1L, 0L, 2L, 0L, 5L, 2L, 1L, 0L, 2L, 4L, 6L, 8L, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 2L, 4L, 3L, 2L, 7L, 6L, 4L, 8L, 7L, 11L, 8L, 11L, 1L, 3L, 0L, 2L, 7L, 8L, 2L, 7L, 6L, 11L, 6L, 12L, 1L, 1L, 0L, 0L, 0L, 1L, 2L, 9L, 6L, 16L, 6L, 10L, 2L, 3L, 3L, 2L, 5L, 2L, 0L, 3L, 6L, 10L, NA, 10L, 1L, 0L, 3L, 1L, 4L, 2L, 3L, 2L, 3L, 11L, 10L, 6L)), .Names = c(Site, EPT.Taxa), class = data.frame, row.names = c(NA, -144L)) subset(d, Site==215) #I would like to plot(density()) of each of the Sites #I tried list-levels(d$Site) lapply(d, FUN=subset, Site==list) #I would like to be able to make a list of the subsets based on Site and then plot them all on one graphics window # Am I working in the right direction- I have just discovered lists (I know I know, I am a little slow) -- Stephen Sefick Research Scientist Southeastern Natural Sciences Academy Let's not spend our time and resources thinking about things that are so little or so large that all they really do for us is puff us up and make us feel like gods. We are mammals, and have not exhausted the annoying little problems of being mammals. -K. Mullis __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is door een geldig ondertekend document. The views expressed in this message and any annex are purely those of the writer and may not be regarded as stating an official position of INBO, as long as the message is not confirmed by a duly signed document __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
... To pick up on what Mark has said: it strikes me that this is related to the simplex, where the bounded nature of the vector space means that normal arithmetical operations (i.e. Euclidean) don't work---that is, they can be used, but the results are wrong. Covariances and correlations for instance, are wrong, something that Pearson noted more than a century ago. Taking logs of ratios of numbers (say a number divdided by geometric mean of the other numbers, as in Aitchison's set of transformations) in this space transfers everything to Euclidean space, so squaring the problem. This is why taking logs fixes things ?? Mark. statquant wrote: Duncan: I think I see what you're saying but the strange thing is that if you use the utility function log(x) rather than x, then the expected values are equal. Somehow, if you are correct and I think you are, then taking the log , fixes the distribution of x which is kind of odd to me. I'm sorry to belabor this non R related discussion and I won't say anything more about it but I worked/talked on this with someone for about a month a few years ago and we gave up so it's interesting for me to see this again. Mark -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Duncan Murdoch Sent: Tuesday, August 26, 2008 8:15 AM To: Jim Lemon Cc: r-help@r-project.org; Mario Subject: Re: [R] Two envelopes problem On 26/08/2008 7:54 AM, Jim Lemon wrote: Hi again, Oops, I meant the expected value of the swap is: 5*0.5 + 20*0.5 = 12.5 Too late, must get to bed. But that is still wrong. You want a conditional expectation, conditional on the observed value (10 in this case). The answer depends on the distribution of the amount X, where the envelopes contain X and 2X. For example, if you knew that X was at most 5, you would know you had just observed 2X, and switching would be a bad idea. The paradox arises because people want to put a nonsensical Unif(0, infinity) distribution on X. The Wikipedia article points out that it can also arise in cases where the distribution on X has infinite mean: a mathematically valid but still nonsensical possibility. Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- View this message in context: http://www.nabble.com/Two-envelopes-problem-tp19150357p19163195.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] savePlot() does not save plot with format set
Did you try: savePlot(test.bmp, type = bmp) Monica Message: 118 Date: Tue, 26 Aug 2008 10:29:13 +0100 From: Luis Ridao Cruz Subject: [R] savePlot() does not save plot with format set To: Message-ID: Content-Type: text/plain; charset=US-ASCII R-help, Whenever I try to save a plot with savePlot the file is not stored in my hard disk with the selected format. Several formats are set and none of them works. I just get the file name with missing extension Editor plot(rnorm(10)) savePlot(test, type=png) savePlot(test, type=bmp) version _ platform i386-pc-mingw32 arch i386 os mingw32 system i386, mingw32 status major 2 minor 7.0 year 2008 month 04 day 22 svn rev 45424 language R version.string R version 2.7.0 (2008-04-22) Thanks in advanced __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] processing subset lists and then plot(density())
definitly need to look into ggplot- very cool. thanks On Tue, Aug 26, 2008 at 10:20 AM, ONKELINX, Thierry [EMAIL PROTECTED] wrote: Here's a solution with ggplot2 library(ggplot2) ggplot(na.omit(d), aes(x = EPT.Taxa, colour = Site)) + geom_density() ggplot(na.omit(d), aes(x = EPT.Taxa)) + geom_density() + facet_grid(Site ~ .) HTH, Thierry ir. Thierry Onkelinx Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest Cel biometrie, methodologie en kwaliteitszorg / Section biometrics, methodology and quality assurance Gaverstraat 4 9500 Geraardsbergen Belgium tel. + 32 54/436 185 [EMAIL PROTECTED] www.inbo.be To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher The plural of anecdote is not data. ~ Roger Brinner The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. ~ John Tukey -Oorspronkelijk bericht- Van: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Namens stephen sefick Verzonden: dinsdag 26 augustus 2008 15:51 Aan: R Help Onderwerp: [R] processing subset lists and then plot(density()) d - structure(list(Site = structure(c(8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L, 8L, 12L, 7L, 6L, 11L, 5L, 10L, 4L, 3L, 2L, 1L, 9L), .Label = c(119, 148, 179, 185, 190, 198, 202, 215, 61, BC, HC, SC), class = factor), EPT.Taxa = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1L, NA, 3L, 1L, 5L, 7L, 3L, 11L, 3L, 14L, 12L, 12L, 0L, 0L, 4L, 0L, 5L, 3L, 2L, 6L, 1L, 8L, 6L, 9L, 1L, 0L, 2L, 0L, 5L, 2L, 1L, 0L, 2L, 4L, 6L, 8L, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 2L, 4L, 3L, 2L, 7L, 6L, 4L, 8L, 7L, 11L, 8L, 11L, 1L, 3L, 0L, 2L, 7L, 8L, 2L, 7L, 6L, 11L, 6L, 12L, 1L, 1L, 0L, 0L, 0L, 1L, 2L, 9L, 6L, 16L, 6L, 10L, 2L, 3L, 3L, 2L, 5L, 2L, 0L, 3L, 6L, 10L, NA, 10L, 1L, 0L, 3L, 1L, 4L, 2L, 3L, 2L, 3L, 11L, 10L, 6L)), .Names = c(Site, EPT.Taxa), class = data.frame, row.names = c(NA, -144L)) subset(d, Site==215) #I would like to plot(density()) of each of the Sites #I tried list-levels(d$Site) lapply(d, FUN=subset, Site==list) #I would like to be able to make a list of the subsets based on Site and then plot them all on one graphics window # Am I working in the right direction- I have just discovered lists (I know I know, I am a little slow) -- Stephen Sefick Research Scientist Southeastern Natural Sciences Academy Let's not spend our time and resources thinking about things that are so little or so large that all they really do for us is puff us up and make us feel like gods. We are mammals, and have not exhausted the annoying little problems of being mammals. -K. Mullis __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is door een geldig ondertekend document. The views expressed in this message and any annex are purely those of the writer and may not be regarded as stating an official position of INBO, as long as the message is not confirmed by a duly signed document -- Stephen Sefick Research Scientist Southeastern Natural Sciences Academy Let's not spend our time and resources thinking about things that are so little or so large that all they really do for us is puff us up and make us feel like gods. We are mammals, and have not exhausted the annoying little problems of being mammals. -K. Mullis __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How to do a meta-analysis plot
On Mon, 25 Aug 2008, Jorge Ivan Velez wrote: My problem is that I don't have the raw data as rmeta _requires_ and, even when I have my data set in the _same_ (?) format that summary(a), when I tried plot(mydata) it doesn't work. No, it doesn't _require_ that. If you just want a forest plot then use forestplot() or metaplot(). If you want to do a meta-analysis from summary data use meta.summaries(). -thomas __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] SQL Primer for R
On Tue, 26 Aug 2008, ivo welch wrote: Sorry, chaps. I need one more: dbDisconnect(con.in) Error in sqliteCloseConnection(conn, ...) : RS-DBI driver: (close the pending result sets before closing this connection) I am pretty sure I have fetched everything there is to be fetched. dbClearResult -thomas Thomas Lumley Assoc. Professor, Biostatistics [EMAIL PROTECTED] University of Washington, Seattle __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
Mark My experience was similarly frustrating. Maybe formulating the problem a bit differently could help to clarify it. State it like this: Someone chooses an amount of money x. He puts 2x/3 of it in one envelope and x/3 in an other. There is no assumption about the distribution of x. If you choose one envelope your expectation is x/2 and changing may lead to a gain or a loss of x/6. In my view there is no basis for a frequentist conditional expectation, conditional on the amount in the first envelope. Of course, after opening the first envelope and finding a, you know for sure that x can only be 3a or 3a/2, but to me there seems to be no basis to assign probabilities to these two alternatives. I am aware of the long lasting discussion and of course this will not end it. Heinz At 14:51 26.08.2008, Mark Leeds wrote: Duncan: I think I see what you're saying but the strange thing is that if you use the utility function log(x) rather than x, then the expected values are equal. Somehow, if you are correct and I think you are, then taking the log , fixes the distribution of x which is kind of odd to me. I'm sorry to belabor this non R related discussion and I won't say anything more about it but I worked/talked on this with someone for about a month a few years ago and we gave up so it's interesting for me to see this again. Mark -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Duncan Murdoch Sent: Tuesday, August 26, 2008 8:15 AM To: Jim Lemon Cc: r-help@r-project.org; Mario Subject: Re: [R] Two envelopes problem On 26/08/2008 7:54 AM, Jim Lemon wrote: Hi again, Oops, I meant the expected value of the swap is: 5*0.5 + 20*0.5 = 12.5 Too late, must get to bed. But that is still wrong. You want a conditional expectation, conditional on the observed value (10 in this case). The answer depends on the distribution of the amount X, where the envelopes contain X and 2X. For example, if you knew that X was at most 5, you would know you had just observed 2X, and switching would be a bad idea. The paradox arises because people want to put a nonsensical Unif(0, infinity) distribution on X. The Wikipedia article points out that it can also arise in cases where the distribution on X has infinite mean: a mathematically valid but still nonsensical possibility. Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] more dot plots on 1 page
Hi,
I would like to have six dot plots on one page. I make the plots in a
'for loop' because it is six times the same graph but for different
subjects (species).
I tried it with par(mfrow=c(3,2), oma=c(3,3,2,0), mar=c(2,2,2,2)); but
this does not work for dot plots apparently.
Then I tried with print(). But then I had to give the dot plots names.
This also does not work.
Any one who has some more ideas? What do I do wrong?
species.vector-c(Lp,Pp,Pl,Ra,Lc,Ml)
for(i in 1:6){
- inputdot[inputdot$sp_==species.vector[i],]
str()
XX7-subset(,chamber==7)
XX8-subset(,chamber==8)
XX9-subset(,chamber==9)
XX10-subset(,chamber==10)
species.vector[i] -
dotplot(XX7$BG_dry~XX7$stress,ylim=c(0,20),ylab=biomass,scales=list(ti
ck.number=10),
panel = function (x, y) {
panel.abline(v=c(1:8),lty=2,col=gray)
panel.xyplot(x, y, pch = 1, col = blue, cex = .75)
panel.xyplot(XX8$stress,XX8$BG_dry, pch = 16, col = red, cex = .6)
panel.xyplot(XX9$stress,XX9$BG_dry, pch = 2, col = blue, cex = .6)
panel.xyplot(XX10$stress,XX10$BG_dry, pch = 17, col = red, cex = .6)
},
key = list(text = list(c(ch7, ch8,ch9,ch10), cex = .75),
points = list(pch = c(1, 16,2,17), col = c(blue,red,blue,red),
cex = .75),
space = bottom, border = T),main=species.vector[i]) )
}
print(Lp, split=c(1,1,2,3), more=TRUE)
print(Pp, split=c(2,1,2,3), more=TRUE)
print(Pl, split=c(1,2,2,3), more=TRUE)
print(Ra, split=c(2,2,2,3), more=TRUE)
print(Lc, split=c(1,3,2,3), more=TRUE)
print(Ml, split=c(2,3,2,3))
Thanks a lot!
Joke
Joke Van den Berge, PhDs
University of Antwerp, Campus Drie Eiken
Department of Biology
Research Group of Plant and Vegetation Ecology
Universiteitsplein 1
B-2610 Wilrijk, Belgium
email: [EMAIL PROTECTED] mailto:[EMAIL PROTECTED]
tel: +32 3 820 22 72 fax: +32 3 820 22 71
[[alternative HTML version deleted]]
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
On 8/26/2008 9:51 AM, Mark Leeds wrote: Duncan: I think I see what you're saying but the strange thing is that if you use the utility function log(x) rather than x, then the expected values are equal. I think that's more or less a coincidence. If I tell you that the two envelopes contain X and 2X, and I also tell you that X is 1,2,3,4, or 5, and you open one and observe 10, then you know that X=5 is the content of the other envelope. The expected utility of switching is negative using any increasing utility function. On the other hand, if we know X is one of 6,7,8,9,10, and you observe a 10, then you know that you got X, so the other envelope contains 2X = 20, and the expected utility is positive. As Heinz says, the problem does not give enough information to come to a decision. The decision *must* depend on the assumed distribution of X, and the problem statement gives no basis for choosing one. There are probably some subjective Bayesians who would assume a particular default prior in a situation like that, but I wouldn't. Duncan Murdoch Somehow, if you are correct and I think you are, then taking the log , fixes the distribution of x which is kind of odd to me. I'm sorry to belabor this non R related discussion and I won't say anything more about it but I worked/talked on this with someone for about a month a few years ago and we gave up so it's interesting for me to see this again. Mark -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Duncan Murdoch Sent: Tuesday, August 26, 2008 8:15 AM To: Jim Lemon Cc: r-help@r-project.org; Mario Subject: Re: [R] Two envelopes problem On 26/08/2008 7:54 AM, Jim Lemon wrote: Hi again, Oops, I meant the expected value of the swap is: 5*0.5 + 20*0.5 = 12.5 Too late, must get to bed. But that is still wrong. You want a conditional expectation, conditional on the observed value (10 in this case). The answer depends on the distribution of the amount X, where the envelopes contain X and 2X. For example, if you knew that X was at most 5, you would know you had just observed 2X, and switching would be a bad idea. The paradox arises because people want to put a nonsensical Unif(0, infinity) distribution on X. The Wikipedia article points out that it can also arise in cases where the distribution on X has infinite mean: a mathematically valid but still nonsensical possibility. Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] embedded examples
I am working on embedding R into some visualization research programs. Can any point me to a collection of embedded and standalone R/C/C++ examples? The documentation is to terse for me to figure out how to develop this and I am looking for some simple examples to study. Thanks and best regards, EBo -- __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] environments
On 8/26/2008 9:44 AM, Peter Dalgaard wrote: Antje wrote: Okay, I see, there is no really easy way (I was wondering whether I can set an environment as default for new created variables). Is there any difference if I call myenv$myvar - 10 or assign(myvar,10, env=myenv) ? No. And with(myenv, myvar - 10) is also the same. However, you do have to be careful with with(myenv, myvar - x) vs. myenv$myvar - x because it might be different x. It's a little unfortunate that with(myenv, myvar - x) works when myenv is an environment but not a list, while mylist - within(mylist, myvar - x) works when mylist is a list but not an environment. Is this something that should be fixed? Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
Hi Duncan: I think I get you. Once one takes expectations, there is an underlying assumption about the distribution of X and , in this problem, we don't have one so taking expectations has no meaning. If the log utility fixing the problem is purely just a coincidence, then it's surely an odd one because log(utility) is often used in economics for expressing how investors view the notion of accumulating capital versus the risk of losing it. I'm not a economist but it's common for them to use log utility to prove theorems about optimal consumption etc. Thanks because I think I see it now by your example below. Mark -Original Message- From: Duncan Murdoch [mailto:[EMAIL PROTECTED] Sent: Tuesday, August 26, 2008 11:26 AM To: Mark Leeds Cc: r-help@r-project.org Subject: Re: [R] Two envelopes problem On 8/26/2008 9:51 AM, Mark Leeds wrote: Duncan: I think I see what you're saying but the strange thing is that if you use the utility function log(x) rather than x, then the expected values are equal. I think that's more or less a coincidence. If I tell you that the two envelopes contain X and 2X, and I also tell you that X is 1,2,3,4, or 5, and you open one and observe 10, then you know that X=5 is the content of the other envelope. The expected utility of switching is negative using any increasing utility function. On the other hand, if we know X is one of 6,7,8,9,10, and you observe a 10, then you know that you got X, so the other envelope contains 2X = 20, and the expected utility is positive. As Heinz says, the problem does not give enough information to come to a decision. The decision *must* depend on the assumed distribution of X, and the problem statement gives no basis for choosing one. There are probably some subjective Bayesians who would assume a particular default prior in a situation like that, but I wouldn't. Duncan Murdoch Somehow, if you are correct and I think you are, then taking the log , fixes the distribution of x which is kind of odd to me. I'm sorry to belabor this non R related discussion and I won't say anything more about it but I worked/talked on this with someone for about a month a few years ago and we gave up so it's interesting for me to see this again. Mark -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Duncan Murdoch Sent: Tuesday, August 26, 2008 8:15 AM To: Jim Lemon Cc: r-help@r-project.org; Mario Subject: Re: [R] Two envelopes problem On 26/08/2008 7:54 AM, Jim Lemon wrote: Hi again, Oops, I meant the expected value of the swap is: 5*0.5 + 20*0.5 = 12.5 Too late, must get to bed. But that is still wrong. You want a conditional expectation, conditional on the observed value (10 in this case). The answer depends on the distribution of the amount X, where the envelopes contain X and 2X. For example, if you knew that X was at most 5, you would know you had just observed 2X, and switching would be a bad idea. The paradox arises because people want to put a nonsensical Unif(0, infinity) distribution on X. The Wikipedia article points out that it can also arise in cases where the distribution on X has infinite mean: a mathematically valid but still nonsensical possibility. Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] embedded examples
You might also want to look at existing visualisation applications that connect with R: * http://ggobi.org * http://rosuda.org/mondrian * http://rosuda.org/software/Gauguin/gauguin.html to name a few. Hadley On Tue, Aug 26, 2008 at 10:31 AM, EBo [EMAIL PROTECTED] wrote: I am working on embedding R into some visualization research programs. Can any point me to a collection of embedded and standalone R/C/C++ examples? The documentation is to terse for me to figure out how to develop this and I am looking for some simple examples to study. Thanks and best regards, EBo -- __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- http://had.co.nz/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] permutation cycles
I'm replying in public to this,
since it relates to an R-help thread of April 15-16, 2008
(and am also adding back 'References:' to one of those postings
of mine).
Mark Segal mark at biostat ucsf
on Mon, 25 Aug 2008 11:25:05 -0700 writes:
Hi Martin,
You posted the following code for obtaining cycles of a permutation:
## MM: Peter's function is so compact and elegant,
## now try to optimize for speed() :
cycLengths2 - function(p) {
n - length(p)
x - integer(n)
ii - seq_len(n)
for (i in ii) {
z - ii[!x][1] # index of first unmarked x[] entry
if (is.na(z)) break
repeat { ## mark x[] - i for those in i-th cycle
x[z] - i
z - p[z]
if (x[z]) break
}
}
## Now, table(x) gives the cycle lengths,
## where split(seq_len(n), x) would give the cycles list
tabulate(x, i - 1L) ## is quite a bit faster than the equivalent
## table(x)
}
But, while it is the case that
split(seq_len(n), x)
gives the members of the respective cycle lists, their cycle ordering
is destroyed: it is replaced with increasing order.
Do you have an efficient way of obtaining (preserved) cycle order?
Yes, indeed, the above code (a version of an original function
by Peter Dalgaard) does not very easily lend itself to get the
cycles in order.
Well, as a matter of fact, I've spent too much time achieving
that but the resulting solution was slower than
the function that Barry Rowlingson had sent me back in April
which exactly solves the problem of finding the
cycle decomposition of permutation vector p[]:
cycles - function(p) {
## From: Barry Rowlingson at Lancaster ac UK
## To: Martin Maechler at ETH Zurich
## Subject: Re: [R] sign(permutation) in R ?
## Date: Tue, 15 Apr 2008 18:15:19 +0100
## with minor tweaks by MM
cycles - list()
ii - seq_along(p)
while(!all(is.na(p))) {
start - ii[!is.na(p)][1]# == which(!is.na(p))[1]
cycle - start
i - start
repeat {
nextV - p[i]
p[i] - NA
if(nextV == start) # cycle is closed
break
## else
cycle - c(cycle,nextV)
i - nextV
}
cycles - c(cycles,list(cycle))
}
return(cycles)
}
You can check that it is working
(and also find that it's reasonably fast),
e.g., with
p - sample(12); rbind(seq_along(p), p)
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12]
123456789101112
p 11453629 121 810 7
str(cycles(p))
List of 2
$ : int [1:7] 1 11 10 8 12 7 9
$ : int [1:5] 2 4 3 5 6
---
Martin Maechler, ETH Zurich
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problem with Integrate for NEF-HS distribution
Hi,
Simply re-write your integrand as follows:
integrand - function(X) {0.5 * cos(theta) * 2 / (exp(-pi*X/2 - X*theta) +
exp(pi*X/2 - X*theta)) }
Now the problem goes away.
theta - -1
integrate(integrand, -Inf, 1)
0.9842315 with absolute error 1.2e-05
This would work for all theta such that abs(theta) -pi/2.
Ravi.
---
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: [EMAIL PROTECTED]
Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On
Behalf Of xia wang
Sent: Tuesday, August 26, 2008 1:01 AM
To: r-help@r-project.org
Subject: [R] Problem with Integrate for NEF-HS distribution
I need to calcuate the cumulative probability for the Natural Exponential
Family - Hyperbolic secant distribution with a parameter theta between -pi/2
and pi/2. The integration should be between 0 and 1 as it is a probability.
The function integrate works fine when the absolute value of theta is not
too large. That is, the NEF-HS distribution is not too skewed. However,
once the theta gets large in absolute value, such as -1 as shown in the
example below, integrate keeps give me error message for non-finite
function when I put the lower bound as -Inf. I suspect that it is caused
by the very heavy tail of the distribution.
Is there any way that I can get around of this and let integrate work for
the skewed distribution? Or is there any other function for integrating in
R-package? Thanks a lot for your advice in advance!
theta--1
sech-function(X) 2/(exp(-X)+exp(X))
integrand - function(X) {0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
integrate(integrand, -3,1)
0.8134389 with absolute error 7e-09
integrate(integrand, -10,1)
0.9810894 with absolute error 5.9e-06
integrate(integrand, -15,1)
0.9840505 with absolute error 7e-09
integrate(integrand, -50,1)
0.9842315 with absolute error 4.4e-05
integrate(integrand, -100,1)
0.9842315 with absolute error 3.2e-05
integrate(integrand, -Inf,1)
Error in integrate(integrand, -Inf, 1) : non-finite function value
Xia
_
Be the filmmaker you always wanted to be-learn how to burn a DVD with
WindowsR.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] R binaries for linux on G4 Powerpc?
Hi Everyone I just installed Yellow Dog Linux (v 5.02) on an old G4 Powerpc replacing the Mac OSX. The system is running beautifully. I would like to install R on this system. I am not sure how to go about doing it. I have installed R on my Ubuntu system as well as on Mac OSX but there are binaries available for those systems from the R Project site. The Yellow Dog Linux is rpm based and has a package management tool called yum. The distribution is derived from Red Hat Linux but is specially designed for the G3/G4/G5 Powerpc architecture. As far as I know, there are no repositories for this particular distro to install R. I was wondering if there is a way to get access to other repositories to install R on my system. I would appreciate any suggestions and guidance. Regards, Tariq [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] package install failure
Hello, the arm package is failing to install, it says I need gnu make but it's not clear how I specify which make to use in the R environment. For example I already have gmake installed. How to I instruct R to use it? * Installing *source* package 'Matrix' ... usage: make [-BeikNnqrstWX] [-D variable] [-d flags] [-f makefile] [-I directory] [-J private] [-j max_jobs] [-m directory] [-T file] [-V variable] [variable=value] [target ...] *** You need GNU make for building this package from the sources ERROR: configuration failed for package 'Matrix' ** Removing '/usr/local/nb40i386/R-2.7.1/lib/R/library/Matrix' best regards, George __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Dramatic slowdown of R 2.7.2?
Dear R users/developers,
simple comparison of code execution time of R 2.7.1 and R 2.7.2 shows a
dramatic slowdown of the newer version. Rprof() identifies .Call
function as a main cause (see the code below). What happened with R
2.7.2?
Kind regards
Marek Wielgosz
Bayes Consulting
# Probably useful info ###
### CPU: Core2Duo T 7300, 2 GB RAM
### WIN XP
### both standard Rblass.dll files
### both pre-compiled binary versions of R
# R 2.7.1
#R version 2.7.1 (2008-06-23)
#i386-pc-mingw32
#
#locale:
#LC_COLLATE=Polish_Poland.1250;LC_CTYPE=Polish_Poland.1250;LC_MONETARY=Polish_Poland.1250;LC_NUMERIC=C;LC_TIME=Polish_Poland.1250
#
#attached base packages:
#[1] stats graphics grDevices utils datasets methods base
# R 2.7.2
#R version 2.7.2 (2008-08-25)
#i386-pc-mingw32
#
#locale:
#LC_COLLATE=Polish_Poland.1250;LC_CTYPE=Polish_Poland.1250;LC_MONETARY=Polish_Poland.1250;LC_NUMERIC=C;LC_TIME=Polish_Poland.1250
#
#attached base packages:
#[1] stats graphics grDevices utils datasets methods base
The following code has been executed (with both R 2.7.1 and R 2.7.2):
Rprof(tmp - tempfile())
sample_size=10
f1=function(k) {solve(matrix(rnorm(4*k^2),2*k,2*k))}
out=vector(,sample_size)
for (i in 1:sample_size) {out[i]=system.time(f1(10^3))[[3]]}
summaryRprof(tmp)
out_summary=matrix(,1,6,dimnames=list(value,c(dim,min,mean,med,sd,max)))
out_summary[1,1]=sample_size
out_summary[1,2]=min(out)
out_summary[1,3]=mean(out)
out_summary[1,4]=median(out)
out_summary[1,5]=sd(out)
out_summary[1,6]=max(out)
out_summary
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] no output when run densityplot...
Hi,
I have downloaded a R script from
http://www.wessa.net/rwasp_edauni.wasp#output.
This script produces a densityplot graphic, amongst others, when is
executed from the web page.
However, when I run it in my machine the *densityplot* function produces
any output, I mean a blank graphic.
But, it's interesting if I run the following lines in the R interactive
console:
y - as.ts(x$pending)
bitmap(file=tmp)
densityplot(~y,col=black,main=x)
dev.off()
Basically, they are the same lines inside the script, I got the expected
output.
I appreciate you advice,
The script code:
totalgraphics - function(server,status) {
par2 = '0'
par1 = '0'
filename - paste(server,status,sep='.')
filename - paste('./pre-processed/',filename,sep='')
v - read.csv(file=filename,sep=',')
x - v[[status]]
par1 - as.numeric(par1)
par2 - as.numeric(par2)
x - as.ts(x)
library(lattice)
filename - paste(server,'_',status,'-','density','.png',sep='')
bitmap(file=filename)
title - paste('Density Plot bw = ',par1,' for
',server,'[',status,']',sep=' ')
*if (par1 == 0)
{
densityplot(~x,col='black',main=title)
} else {
if (par1 0) {
densityplot(~x,col='black',main=title,bw=par1)
}
}*
dev.off()
filename - paste(server,'_',status,'-','sequence','.png',sep='')
bitmap(file=filename)
title - paste('Run Sequence Plot for',server,'[',status,']',sep=' ')
plot(x,type='l',main=title,xlab='time or index',ylab='value')
grid()
dev.off()
filename - paste(server,'_',status,'-','hist','.png',sep='')
bitmap(file=filename)
hist(x)
grid()
dev.off()
filename - paste(server,'_',status,'-','quantile','.png',sep='')
bitmap(file=filename)
qqnorm(x)
qqline(x)
grid()
dev.off()
if (par2 0)
{
filename - paste(server,'_',status,'-','lagplot1','.png',sep='')
bitmap(file=filename)
dum - cbind(lag(x,k=1),x)
dum
dum1 - dum[2:length(x),]
dum1
z - as.data.frame(dum1)
z
title - paste('Lag plot (k=1), lowess, and regression line for
',server,'[',status,']',sep=' ')
plot(z,main=title)
lines(lowess(z))
abline(lm(z))
dev.off()
if (par2 1) {
filename - paste(server,'_',status,'-','lagplot2','.png',sep='')
bitmap(file=filename)
dum - cbind(lag(x,k=par2),x)
dum
dum1 - dum[(par2+1):length(x),]
dum1
z - as.data.frame(dum1)
z
mylagtitle - 'Lag plot (k='
mylagtitle - paste(mylagtitle,par2,sep='')
mylagtitle - paste(mylagtitle,'), and lowess',sep='')
title - paste(mylagtitle,' for ',server,'[',status,']',sep=' ')
plot(z,main=title)
lines(lowess(z))
dev.off()
}
filename - paste(server,'_',status,'-','autocorre','.png',sep='')
bitmap(file=filename)
title - paste('Autocorrelation Function',' for
',server,'[',status,']',sep=' ')
acf(x,lag.max=par2,main=title)
grid()
dev.off()
}
summary(x);
}
servers -
c(komolongma.ece.uprm.edu,sakura.hpcc.jp,fsvc001.asc.hpcc.jp,rocks-52.sdsc.edu,rocks-153.sdsc.edu)
status - c(unsubmitted,active,pending)
for (i in servers) {
for (j in status) {
totalgraphics(i,j)
}
}
[[alternative HTML version deleted]]
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problem with Integrate for NEF-HS distribution
Thanks so much. It works well in my MCMC sampling. May I know why re-writing
the integrand can solve the problem? I have been thinking it was the skewness
of the distribution brought the error. Thanks!
Xia
From: [EMAIL PROTECTED]
To: [EMAIL PROTECTED]; r-help@r-project.org
Subject: RE: [R] Problem with Integrate for NEF-HS distribution
Date: Tue, 26 Aug 2008 11:59:44 -0400
Hi,
Simply re-write your integrand as follows:
integrand - function(X) {0.5 * cos(theta) * 2 / (exp(-pi*X/2 - X*theta) +
exp(pi*X/2 - X*theta)) }
Now the problem goes away.
theta - -1
integrate(integrand, -Inf, 1)
0.9842315 with absolute error 1.2e-05
This would work for all theta such that abs(theta) -pi/2.
Ravi.
---
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: [EMAIL PROTECTED]
Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On
Behalf Of xia wang
Sent: Tuesday, August 26, 2008 1:01 AM
To: r-help@r-project.org
Subject: [R] Problem with Integrate for NEF-HS distribution
I need to calcuate the cumulative probability for the Natural Exponential
Family - Hyperbolic secant distribution with a parameter theta between -pi/2
and pi/2. The integration should be between 0 and 1 as it is a probability.
The function integrate works fine when the absolute value of theta is not
too large. That is, the NEF-HS distribution is not too skewed. However,
once the theta gets large in absolute value, such as -1 as shown in the
example below, integrate keeps give me error message for non-finite
function when I put the lower bound as -Inf. I suspect that it is caused
by the very heavy tail of the distribution.
Is there any way that I can get around of this and let integrate work for
the skewed distribution? Or is there any other function for integrating in
R-package? Thanks a lot for your advice in advance!
theta--1
sech-function(X) 2/(exp(-X)+exp(X))
integrand - function(X) {0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
integrate(integrand, -3,1)
0.8134389 with absolute error 7e-09
integrate(integrand, -10,1)
0.9810894 with absolute error 5.9e-06
integrate(integrand, -15,1)
0.9840505 with absolute error 7e-09
integrate(integrand, -50,1)
0.9842315 with absolute error 4.4e-05
integrate(integrand, -100,1)
0.9842315 with absolute error 3.2e-05
integrate(integrand, -Inf,1)
Error in integrate(integrand, -Inf, 1) : non-finite function value
Xia
_
Be the filmmaker you always wanted to be-learn how to burn a DVD with
WindowsR.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
_
[[alternative HTML version deleted]]
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] plot3d origin
Hi all, I am trying to do a 3d plot where the x,y,z axes intersects with the origin (0,0,0) using the plot3d() funtion in the rgl package without success. I looked back at the past archives on this subject and someone suggested using djmrgl package. I searched and found it, installed it but when I try to load it I get the error ... Error in inDL(x, as.logical(local), as.logical(now), ...) : unable to load shared library 'C:/PROGRA~1/R/R-27~1.0/library/djmrgl/libs/djmrgl.dll': LoadLibrary failure: The specified procedure could not be found. (Using Windows XP) The file is in C:\Program Files\R\R-2.7.0\library\djmrgl\libs\ and I have tried putting this in the PATH variable but it still doesn't make a difference. I would appreciate it if someone could help me on this. Kind Regards Chibisi [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] sequence with start and stop positions
On Tue, 26 Aug 2008, Chris Oldmeadow wrote:
Hi,
I have a vector of start positions, and another vector of stop positions,
eg start-c(1,20,50)
stop-c(7,25,53)
Is there a quick way to create a sequence from these vectors?
new-c(1,2,3,4,5,6,7,20,21,22,23,24,25,50,51,52,53)
Vectorize the process.
start2 - rep( start, stop-start+1 )
lens - stop-start+1
offset - rep(cumsum(c( 0, lens )),c( lens ,0 ))
seq(along=offset)-offset+start2-1
[1] 1 2 3 4 5 6 7 20 21 22 23 24 25 50 51 52 53
HTH,
Chuck
the way Im doing it at the moment is
pos-seq(start[1],stop[1])
for (i in 2:length(start)){
new-seq(start[i],stop[i])
pos-c(pos,new)
}
This works on small data, but its very inefficient on large vectors, and is
taking forever!
Anybody no a better way?
many thanks,
Chris
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Charles C. Berry(858) 534-2098
Dept of Family/Preventive Medicine
E mailto:[EMAIL PROTECTED] UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
Duncan: Just one more thing which Heinz alerted me to. Suppose that we changed the game so that instead of being double or half of X, we said that one envelope will contain X + 5 and the other contains X-5. So someone opens it and sees 10 dollars. Now, their remaining choices are 5 and 15 so the expectation of switching is the same = 10. So, in this case, we don't know the distribution of X and yet the game is fair. This is why , although I like your example, I still think the issue has something to do with percentages versus additions. In finance, the cumulative return is ( in non continuous time ) over some horizon is a productive of the individual returns over whatever intervals you want to break the horizon into. In order to make things nicer statistically ( and for other reasons too like making the assumption of normaility somkewhat more plausible ) , finance people take the log of this product in order to to transform the cumulative return into an additive measure. So, I think there's still something going on with units as far as adding versus multiplying ? but I'm not sure what and I do still see what you're saying in your example. Thanks. Mark On Tue, Aug 26, 2008 at 11:44 AM, Mark Leeds wrote: Hi Duncan: I think I get you. Once one takes expectations, there is an underlying assumption about the distribution of X and , in this problem, we don't have one so taking expectations has no meaning. If the log utility fixing the problem is purely just a coincidence, then it's surely an odd one because log(utility) is often used in economics for expressing how investors view the notion of accumulating capital versus the risk of losing it. I'm not a economist but it's common for them to use log utility to prove theorems about optimal consumption etc. Thanks because I think I see it now by your example below. Mark -Original Message- From: Duncan Murdoch [mailto:[EMAIL PROTECTED] Sent: Tuesday, August 26, 2008 11:26 AM To: Mark Leeds Cc: r-help@r-project.org Subject: Re: [R] Two envelopes problem On 8/26/2008 9:51 AM, Mark Leeds wrote: Duncan: I think I see what you're saying but the strange thing is that if you use the utility function log(x) rather than x, then the expected values are equal. I think that's more or less a coincidence. If I tell you that the two envelopes contain X and 2X, and I also tell you that X is 1,2,3,4, or 5, and you open one and observe 10, then you know that X=5 is the content of the other envelope. The expected utility of switching is negative using any increasing utility function. On the other hand, if we know X is one of 6,7,8,9,10, and you observe a 10, then you know that you got X, so the other envelope contains 2X = 20, and the expected utility is positive. As Heinz says, the problem does not give enough information to come to a decision. The decision *must* depend on the assumed distribution of X, and the problem statement gives no basis for choosing one. There are probably some subjective Bayesians who would assume a particular default prior in a situation like that, but I wouldn't. Duncan Murdoch Somehow, if you are correct and I think you are, then taking the log , fixes the distribution of x which is kind of odd to me. I'm sorry to belabor this non R related discussion and I won't say anything more about it but I worked/talked on this with someone for about a month a few years ago and we gave up so it's interesting for me to see this again. Mark -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Duncan Murdoch Sent: Tuesday, August 26, 2008 8:15 AM To: Jim Lemon Cc: r-help@r-project.org; Mario Subject: Re: [R] Two envelopes problem On 26/08/2008 7:54 AM, Jim Lemon wrote: Hi again, Oops, I meant the expected value of the swap is: 5*0.5 + 20*0.5 = 12.5 Too late, must get to bed. But that is still wrong. You want a conditional expectation, conditional on the observed value (10 in this case). The answer depends on the distribution of the amount X, where the envelopes contain X and 2X. For example, if you knew that X was at most 5, you would know you had just observed 2X, and switching would be a bad idea. The paradox arises because people want to put a nonsensical Unif(0, infinity) distribution on X. The Wikipedia article points out that it can also arise in cases where the distribution on X has infinite mean: a mathematically valid but still nonsensical possibility. Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list
Re: [R] Two envelopes problem
So if I put $10 and $20 into the envelopes, then told you that the values were multiples of $10, it would be wrong for you to assess probabilities on $100, $1,000,000 and so on? :-) But what if you reasoned that there were far more multiples of 10 above 20 than below 20? What if I was really evil and the multiples of 10 that I always put in the envelopes were 0x$10 and 0x$10. When you opened an empty envelope, what are the odds the other has $1,000,000? What are the odds it has $10? It seems the information given can be misleading. Robert Farley Metro www.Metro.net -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Duncan Murdoch Sent: Tuesday, August 26, 2008 05:15 To: Jim Lemon Cc: r-help@r-project.org; Mario Subject: Re: [R] Two envelopes problem On 26/08/2008 7:54 AM, Jim Lemon wrote: Hi again, Oops, I meant the expected value of the swap is: 5*0.5 + 20*0.5 = 12.5 Too late, must get to bed. But that is still wrong. You want a conditional expectation, conditional on the observed value (10 in this case). The answer depends on the distribution of the amount X, where the envelopes contain X and 2X. For example, if you knew that X was at most 5, you would know you had just observed 2X, and switching would be a bad idea. The paradox arises because people want to put a nonsensical Unif(0, infinity) distribution on X. The Wikipedia article points out that it can also arise in cases where the distribution on X has infinite mean: a mathematically valid but still nonsensical possibility. Duncan Murdoch __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
Duncan Murdoch [EMAIL PROTECTED] 26/08/2008 16:17:34
If this is indeed the case, switch; the expected gain is
positive because _you already have the information that you hold the
median value of the three possibilities_. The tendency when
presented
with the problem is to reason as if this is the case.
No, you don't know that A is the median. That doesn't even make
sense,
based on the context of the question: there is no 3-valued random
variable here.
This is my point; apologies if I put it badly. The fact is that you
_don't_ hold the median value and that this is indeed a silly way of
looking at it. My assertion was that this is the way many folk DO look
at it, and that this results in an apparent paradox.
In fact, you inadvertently gave an example when you said
The unopened envelope can hold only two values, given
that yours contains A.
True - for a rather restricted meaning of 'true'.
As written, it implicitly allows three values; A for the envelope you
hold, and two more (2A and 0.5A) for the alternatives you permit. The
usual (and incorrect) expected gain calculation uses all three; 2A-A for
one choice and 0.5-A for the other. To do that at all, we must be
playing with three possible values for the contents of an envelope.
This clearly cannot be the case if there are only two possible values,
as we are told in posing the problem. The situation is that you hold
_either_ A _or_ 2A and the other envelope holds (respectively) 2A or A.
We just don't know what A is until we open both envelopes. So for a
given pair of envelopes, it is the choice of coefficient (1 or 2) that
is random.
If I were to describe this in terms of a random variable I would have
to assume an unknown but - for this run - fixed value A multiplied by a
two-valued random variable with possible values 1 and 2, would I not? We
surely can't have both 0.5 and 2 in our distribution at the same time,
because the original proposition said there were only two possibilities
and they differ by a factor of 2, not 4.
You have no basis for putting a probability distribution on those
values, because you don't know the distribution of X.
But we do know; or at least we know the distribution of the
coefficient. We have two envelopes of which one is selected at random,
and the ratio of values is 2. On that basis, assigning a 50:50
probability of ending up with A or 2A on first selection seems
uncontroversial.
But I'm more than willing to have my intuition corrected - possibly
off-line, of course, since this stopped being R a while back!
Steve E
***
This email and any attachments are confidential. Any use...{{dropped:8}}
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Maintaining repeated ID numbers when transposing with reshape
Thank you for your suggestion, I will play around with it. I guess my concern is that I need each test result to occupy its own cell rather than have one or more in the same row. Adaikalavan Ramasamy-2 wrote: There might be a more elegant way of doing this but here is a way of doing it without reshape(). df - data.frame( ID=c(1,1,1,1,2,2), TEST=c(A,A,B,C,B,B), RESULT=c(17,12,15,12,8,9) ) df.s - split( df, df$ID ) out - sapply( df.s, function(m) tapply( m$RESULT, m$TEST, paste, collapse=, ) ) t(out) A B C 1 17,12 15 12 2 NA 8,9 NA Not the same output as you wanted. This makes more sense unless you have a reason to priotize 17 instead of 12 in the first row. Regards, Adai jcarmichael wrote: I have a dataset in long format that looks something like this: ID TESTRESULT 1 A 17 1 A 12 1 B 15 1 C 12 2 B 8 2 B 9 Now what I would like to do is transpose it like so: IDTEST ATEST BTEST C 1 17 15 12 1 12.. 2 . 8. 2 . 9. When I try: reshape(mydata, v.names=result, idvar=id,timevar=test, direction=wide) It gives me only the first occurrence of each test for each subject. How can I transpose my dataset in this way without losing information about repeated tests? Any help or guidance would be appreciated! Thanks! __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- View this message in context: http://www.nabble.com/Maintaining-repeated-ID-numbers-when-transposing-with-reshape-tp19151853p19166910.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Latin squares in R
Dear R Gurus: What would be the best way to evaluate a Latin square problem, please? Does it work in Rcmdr, please? thanks, Edna Bell __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
On 8/26/2008 1:12 PM, [EMAIL PROTECTED] wrote: Duncan: Just one more thing which Heinz alerted me to. Suppose that we changed the game so that instead of being double or half of X, we said that one envelope will contain X + 5 and the other contains X-5. So someone opens it and sees 10 dollars. Now, their remaining choices are 5 and 15 so the expectation of switching is the same = 10. So, in this case, we don't know the distribution of X and yet the game is fair. No, you can't do that calculation without knowing the distribution. Take my first example again: if X is known in advance to be 1,2,3,4, or 5, then an observation of 10 must be X+5, so you'd expect 0 in the other envelope. The conditional expectation depends on the full distribution, and if you aren't willing to state that, you can't do the calculation properly. Duncan Murdoch This is why , although I like your example, I still think the issue has something to do with percentages versus additions. In finance, the cumulative return is ( in non continuous time ) over some horizon is a productive of the individual returns over whatever intervals you want to break the horizon into. In order to make things nicer statistically ( and for other reasons too like making the assumption of normaility somkewhat more plausible ) , finance people take the log of this product in order to to transform the cumulative return into an additive measure. So, I think there's still something going on with units as far as adding versus multiplying ? but I'm not sure what and I do still see what you're saying in your example. Thanks. Mark On Tue, Aug 26, 2008 at 11:44 AM, Mark Leeds wrote: Hi Duncan: I think I get you. Once one takes expectations, there is an underlying assumption about the distribution of X and , in this problem, we don't have one so taking expectations has no meaning. If the log utility fixing the problem is purely just a coincidence, then it's surely an odd one because log(utility) is often used in economics for expressing how investors view the notion of accumulating capital versus the risk of losing it. I'm not a economist but it's common for them to use log utility to prove theorems about optimal consumption etc. Thanks because I think I see it now by your example below. Mark -Original Message- From: Duncan Murdoch [mailto:[EMAIL PROTECTED] Sent: Tuesday, August 26, 2008 11:26 AM To: Mark Leeds Cc: r-help@r-project.org Subject: Re: [R] Two envelopes problem On 8/26/2008 9:51 AM, Mark Leeds wrote: Duncan: I think I see what you're saying but the strange thing is that if you use the utility function log(x) rather than x, then the expected values are equal. I think that's more or less a coincidence. If I tell you that the two envelopes contain X and 2X, and I also tell you that X is 1,2,3,4, or 5, and you open one and observe 10, then you know that X=5 is the content of the other envelope. The expected utility of switching is negative using any increasing utility function. On the other hand, if we know X is one of 6,7,8,9,10, and you observe a 10, then you know that you got X, so the other envelope contains 2X = 20, and the expected utility is positive. As Heinz says, the problem does not give enough information to come to a decision. The decision *must* depend on the assumed distribution of X, and the problem statement gives no basis for choosing one. There are probably some subjective Bayesians who would assume a particular default prior in a situation like that, but I wouldn't. Duncan Murdoch Somehow, if you are correct and I think you are, then taking the log , fixes the distribution of x which is kind of odd to me. I'm sorry to belabor this non R related discussion and I won't say anything more about it but I worked/talked on this with someone for about a month a few years ago and we gave up so it's interesting for me to see this again. Mark -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Duncan Murdoch Sent: Tuesday, August 26, 2008 8:15 AM To: Jim Lemon Cc: r-help@r-project.org; Mario Subject: Re: [R] Two envelopes problem On 26/08/2008 7:54 AM, Jim Lemon wrote: Hi again, Oops, I meant the expected value of the swap is: 5*0.5 + 20*0.5 = 12.5 Too late, must get to bed. But that is still wrong. You want a conditional expectation, conditional on the observed value (10 in this case). The answer depends on the distribution of the amount X, where the envelopes contain X and 2X. For example, if you knew that X was at most 5, you would know you had just observed 2X, and switching would be a bad idea. The paradox arises because people want to put a nonsensical Unif(0, infinity) distribution on X. The Wikipedia article points out that it can also arise in cases where the distribution on X has infinite
Re: [R] plot3d origin
On 8/26/2008 1:05 PM, Chibisi Chima-Okereke wrote:
Hi all,
I am trying to do a 3d plot where the x,y,z axes intersects with the origin
(0,0,0) using the plot3d() funtion in the rgl package without success. I
looked back at the past archives on this subject and someone suggested using
djmrgl package. I searched and found it, installed it but when I try to load
it I get the error ...
Don't use djmrgl. I don't support it any more.
Use rgl. If the automatic axes don't work, then use axis3d (or even
segments3d) to draw them exactly where you want.
For example,
points3d(rnorm(100),rnorm(100),rnorm(100))
axis3d('x', pos=c(0,0,0))
axis3d('y', pos=c(0,0,0))
axis3d('z', pos=c(0,0,0))
Duncan Murdoch
Error in inDL(x, as.logical(local), as.logical(now), ...) :
unable to load shared library
'C:/PROGRA~1/R/R-27~1.0/library/djmrgl/libs/djmrgl.dll':
LoadLibrary failure: The specified procedure could not be found.
(Using Windows XP)
The file is in C:\Program Files\R\R-2.7.0\library\djmrgl\libs\ and I have
tried putting this in the PATH variable but it still doesn't make a
difference.
I would appreciate it if someone could help me on this.
Kind Regards
Chibisi
[[alternative HTML version deleted]]
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Dramatic slowdown of R 2.7.2?
Marek Wielgosz wrote:
Dear R users/developers,
simple comparison of code execution time of R 2.7.1 and R 2.7.2 shows
a dramatic slowdown of the newer version. Rprof() identifies .Call
function as a main cause (see the code below). What happened with R
2.7.2?
I don't see much of a difference, running on Linux or under Wine on
Linux though, and scaled down to 3*10^2. With a large matrix inversion
as the workhorse, it is hardly surprising that .Call eats most of the
time, but one could easily get the idea that you aren't picking up the
same BLAS in both cases(?)
Kind regards
Marek Wielgosz
Bayes Consulting
# Probably useful info ###
### CPU: Core2Duo T 7300, 2 GB RAM
### WIN XP
### both standard Rblass.dll files
### both pre-compiled binary versions of R
# R 2.7.1
#R version 2.7.1 (2008-06-23)
#i386-pc-mingw32
#
#locale:
#LC_COLLATE=Polish_Poland.1250;LC_CTYPE=Polish_Poland.1250;LC_MONETARY=Polish_Poland.1250;LC_NUMERIC=C;LC_TIME=Polish_Poland.1250
#
#attached base packages:
#[1] stats graphics grDevices utils datasets methods base
# R 2.7.2
#R version 2.7.2 (2008-08-25)
#i386-pc-mingw32
#
#locale:
#LC_COLLATE=Polish_Poland.1250;LC_CTYPE=Polish_Poland.1250;LC_MONETARY=Polish_Poland.1250;LC_NUMERIC=C;LC_TIME=Polish_Poland.1250
#
#attached base packages:
#[1] stats graphics grDevices utils datasets methods base
The following code has been executed (with both R 2.7.1 and R 2.7.2):
Rprof(tmp - tempfile())
sample_size=10
f1=function(k) {solve(matrix(rnorm(4*k^2),2*k,2*k))}
out=vector(,sample_size)
for (i in 1:sample_size) {out[i]=system.time(f1(10^3))[[3]]}
summaryRprof(tmp)
out_summary=matrix(,1,6,dimnames=list(value,c(dim,min,mean,med,sd,max)))
out_summary[1,1]=sample_size
out_summary[1,2]=min(out)
out_summary[1,3]=mean(out)
out_summary[1,4]=median(out)
out_summary[1,5]=sd(out)
out_summary[1,6]=max(out)
out_summary
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B
c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
(*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918
~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
beautiful. now, i think i got it. the fact that the calculation works in the additive case doesn't make the calculation correct. the expected value calculation is kind of assuming that the person putting the things in the envelopes chooses what's in the second envelope AFTER knowing what's already in the other one. But, if he/she does know and still chooses half or double with 50/50 chance, then one can't use the original wealth as the reference point/ initial expected value because the envelope chooser isn't using it. I think it's best for me to think about it like that but your example is good to fall back on when I get confused again 2 years from now. Thank you very much for your patience and explanation. On Tue, Aug 26, 2008 at 2:06 PM, Duncan Murdoch wrote: On 8/26/2008 1:12 PM, [EMAIL PROTECTED] wrote: Duncan: Just one more thing which Heinz alerted me to. Suppose that we changed the game so that instead of being double or half of X, we said that one envelope will contain X + 5 and the other contains X-5. So someone opens it and sees 10 dollars. Now, their remaining choices are 5 and 15 so the expectation of switching is the same = 10. So, in this case, we don't know the distribution of X and yet the game is fair. No, you can't do that calculation without knowing the distribution. Take my first example again: if X is known in advance to be 1,2,3,4, or 5, then an observation of 10 must be X+5, so you'd expect 0 in the other envelope. The conditional expectation depends on the full distribution, and if you aren't willing to state that, you can't do the calculation properly. Duncan Murdoch This is why , although I like your example, I still think the issue has something to do with percentages versus additions. In finance, the cumulative return is ( in non continuous time ) over some horizon is a productive of the individual returns over whatever intervals you want to break the horizon into. In order to make things nicer statistically ( and for other reasons too like making the assumption of normaility somkewhat more plausible ) , finance people take the log of this product in order to to transform the cumulative return into an additive measure. So, I think there's still something going on with units as far as adding versus multiplying ? but I'm not sure what and I do still see what you're saying in your example. Thanks. Mark On Tue, Aug 26, 2008 at 11:44 AM, Mark Leeds wrote: Hi Duncan: I think I get you. Once one takes expectations, there is an underlying assumption about the distribution of X and , in this problem, we don't have one so taking expectations has no meaning. If the log utility fixing the problem is purely just a coincidence, then it's surely an odd one because log(utility) is often used in economics for expressing how investors view the notion of accumulating capital versus the risk of losing it. I'm not a economist but it's common for them to use log utility to prove theorems about optimal consumption etc. Thanks because I think I see it now by your example below. Mark -Original Message- From: Duncan Murdoch [mailto:[EMAIL PROTECTED] Sent: Tuesday, August 26, 2008 11:26 AM To: Mark Leeds Cc: r-help@r-project.org Subject: Re: [R] Two envelopes problem On 8/26/2008 9:51 AM, Mark Leeds wrote: Duncan: I think I see what you're saying but the strange thing is that if you use the utility function log(x) rather than x, then the expected values are equal. I think that's more or less a coincidence. If I tell you that the two envelopes contain X and 2X, and I also tell you that X is 1,2,3,4, or 5, and you open one and observe 10, then you know that X=5 is the content of the other envelope. The expected utility of switching is negative using any increasing utility function. On the other hand, if we know X is one of 6,7,8,9,10, and you observe a 10, then you know that you got X, so the other envelope contains 2X = 20, and the expected utility is positive. As Heinz says, the problem does not give enough information to come to a decision. The decision *must* depend on the assumed distribution of X, and the problem statement gives no basis for choosing one. There are probably some subjective Bayesians who would assume a particular default prior in a situation like that, but I wouldn't. Duncan Murdoch Somehow, if you are correct and I think you are, then taking the log , fixes the distribution of x which is kind of odd to me. I'm sorry to belabor this non R related discussion and I won't say anything more about it but I worked/talked on this with someone for about a month a few years ago and we gave up so it's interesting for me to see this again. Mark -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of Duncan Murdoch Sent: Tuesday, August
[R] stats tests on large tables
I have a large table of data which i can read in using read.table. I then want to conduct various tests on the data. Is there a simple way to conduct these tests by specifying the column headers which relate to different conditions of the experiment? e.g. data1 - read.table(test.table, header=TRUE) t.test(test~control, data = data1) Which doesn't work and results in the error Error in t.test.formula(test ~ control, data = data1) : grouping factor must have exactly 2 levels Many Thanks Richard [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problem with Integrate for NEF-HS distribution
Thanks. I revised the code but it doesn't make difference. The cause of the
error is just as stated in the ?integrate If the function is approximately
constant (in particular, zero) over nearly all
its range it is possible that the result and error estimate may be seriously
wrong. I have searched R-archive. It may not be feasible to solve the
problem by rectangle approximation as some postings suggested because the
integration is in fact within MCMC samplings as follows. The lower bound is
always -Inf. The upper bound and the value of theta changes for each sampling
in MCMC.
I tried to multiple the integrand by a large number ( like 10^16) and changes
the rel.tol. It does help for some range of theta but not all.
Xia
like.fun - function(beta,theta) {
integrand - function(X,theta) {0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
upper-x%*%beta
prob-matrix(NA, nrow(covariate),1)
for (i in 1:nrow(covariate))
{prob[i]-integrate(integrand,lower=-Inf, upper=upper[i],theta=theta)$value
}
likelihood - sum(log(dbinom(y,n,prob)))
return(likelihood)
}
Date: Tue, 26 Aug 2008 00:49:02 -0700
From: [EMAIL PROTECTED]
Subject: Re: [R] Problem with Integrate for NEF-HS distribution
To: [EMAIL PROTECTED]
Shouldn't you define
integrand - function(X,theta) ..
and not
integrand - function(X) .
--- On Tue, 26/8/08, xia wang [EMAIL PROTECTED] wrote:
From: xia wang [EMAIL PROTECTED]
Subject: [R] Problem with Integrate for NEF-HS distribution
To: r-help@r-project.org
Received: Tuesday, 26 August, 2008, 3:00 PM
I need to calcuate the cumulative probability for the
Natural Exponential Family - Hyperbolic secant distribution
with a parameter theta between -pi/2 and pi/2. The
integration should be between 0 and 1 as it is a
probability.
The function integrate works fine when the
absolute value of theta is not too large. That is, the
NEF-HS distribution is not too skewed. However, once the
theta gets large in absolute value, such as -1 as shown in
the example below, integrate keeps give me error
message for non-finite function when I put the
lower bound as -Inf. I suspect that it is caused by the
very heavy tail of the distribution.
Is there any way that I can get around of this and let
integrate work for the skewed distribution? Or
is there any other function for integrating in R-package?
Thanks a lot for your advice in advance!
theta--1
sech-function(X) 2/(exp(-X)+exp(X))
integrand - function(X)
{0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
integrate(integrand, -3,1)
0.8134389 with absolute error 7e-09
integrate(integrand, -10,1)
0.9810894 with absolute error 5.9e-06
integrate(integrand, -15,1)
0.9840505 with absolute error 7e-09
integrate(integrand, -50,1)
0.9842315 with absolute error 4.4e-05
integrate(integrand, -100,1)
0.9842315 with absolute error 3.2e-05
integrate(integrand, -Inf,1)
Error in integrate(integrand, -Inf, 1) : non-finite
function value
Xia
_
Be the filmmaker you always wanted to belearn how to
burn a DVD with Windows®.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained,
reproducible code.
_
Be the filmmaker you always wanted to belearn how to burn a DVD with Windows®.
[[alternative HTML version deleted]]
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] sequence with start and stop positions
And somewhat faster still (YMMV):
unlist(mapply(:, start, stop))
HTH,
Ray Brownrigg
Victoria University of Wellington
Christos Hatzis wrote:
Try:
unlist(mapply(seq, c(1,20,50), c(7,25,53)))
[1] 1 2 3 4 5 6 7 20 21 22 23 24 25 50 51 52 53
-Christos
-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Chris Oldmeadow
Sent: Tuesday, August 26, 2008 12:42 AM
To: r-help@r-project.org
Subject: [R] sequence with start and stop positions
Hi,
I have a vector of start positions, and another vector of
stop positions,
eg start-c(1,20,50)
stop-c(7,25,53)
Is there a quick way to create a sequence from these vectors?
new-c(1,2,3,4,5,6,7,20,21,22,23,24,25,50,51,52,53)
the way Im doing it at the moment is
pos-seq(start[1],stop[1])
for (i in 2:length(start)){
new-seq(start[i],stop[i])
pos-c(pos,new)
}
This works on small data, but its very inefficient on large
vectors, and is taking forever!
Anybody no a better way?
many thanks,
Chris
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problem with Integrate for NEF-HS distribution
Try the following:
sech-function(X) 2/(exp(-X)+exp(X))
your.integrand - function(X) {0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
my.integrand - function(X) {0.5 * cos(theta) * 2 / (exp(-pi*X/2 - X*theta)
+ exp(pi*X/2 - X*theta)) }
theta - -1
my.integrand(-800)
your.integrand(-800)
Do you see the problem?
Ravi.
---
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: [EMAIL PROTECTED]
Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
_
From: xia wang [mailto:[EMAIL PROTECTED]
Sent: Tuesday, August 26, 2008 12:59 PM
To: Ravi Varadhan
Cc: r-help@r-project.org
Subject: RE: [R] Problem with Integrate for NEF-HS distribution
Thanks so much. It works well in my MCMC sampling. May I know why
re-writing the integrand can solve the problem? I have been thinking it was
the skewness of the distribution brought the error. Thanks!
Xia
From: [EMAIL PROTECTED]
To: [EMAIL PROTECTED]; r-help@r-project.org
Subject: RE: [R] Problem with Integrate for NEF-HS distribution
Date: Tue, 26 Aug 2008 11:59:44 -0400
Hi,
Simply re-write your integrand as follows:
integrand - function(X) {0.5 * cos(theta) * 2 / (exp(-pi*X/2 - X*theta) +
exp(pi*X/2 - X*theta)) }
Now the problem goes away.
theta - -1
integrate(integrand, -Inf, 1)
0.9842315 with absolute error 1.2e-05
This would work for all theta such that abs(theta) -pi/2.
Ravi.
---
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: [EMAIL PROTECTED]
Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED]
On
Behalf Of xia wang
Sent: Tuesday, August 26, 2008 1:01 AM
To: r-help@r-project.org
Subject: [R] Problem with Integrate for NEF-HS distribution
I need to calcuate the cumulative probability for the Natural Exponential
Family - Hyperbolic secant distribution with a parameter theta between
-pi/2
and pi/2. The integration should be between 0 and 1 as it is a
probability.
The function integrate works fine when the absolute value of theta is
not
too large. That is, the NEF-HS distribution is not too skewed. However,
once the theta gets large in absolute value, such as -1 as shown in the
example below, integrate keeps give me error message for non-finite
function when I put the lower bound as -Inf. I suspect that it is caused
by the very heavy tail of the distribution.
Is there any way that I can get around of this and let integrate work
for
the skewed distribution? Or is there any other function for integrating in
R-package? Thanks a lot for your advice in advance!
theta--1
sech-function(X) 2/(exp(-X)+exp(X))
integrand - function(X) {0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
integrate(integrand, -3,1)
0.8134389 with absolute error 7e-09
integrate(integrand, -10,1)
0.9810894 with absolute error 5.9e-06
integrate(integrand, -15,1)
0.9840505 with absolute error 7e-09
integrate(integrand, -50,1)
0.9842315 with absolute error 4.4e-05
integrate(integrand, -100,1)
0.9842315 with absolute error 3.2e-05
integrate(integrand, -Inf,1)
Error in integrate(integrand, -Inf, 1) : non-finite function value
Xia
_
Be the filmmaker you always wanted to be-learn how to burn a DVD with
WindowsR.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
_
See what people are saying about Windows Live. Check out featured posts.
Check It Out!
http://www.windowslive.com/connect?ocid=TXT_TAGLM_WL_connect2_082008
[[alternative HTML version deleted]]
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] Two envelopes problem
A very important point is missing here. If there is x in one envelope
and 2x in the other the expected gain is 3x/2. If the idea is to
switch after observing the second envelope the expected gain is again
3x/2. In the case being put x will be either 5 or 10. But x is a
parameter and in this case does not have a probability distribution.
Then one can not take an expectation with respect to x.
John Frain
2008/8/26 S Ellison [EMAIL PROTECTED]:
Duncan Murdoch [EMAIL PROTECTED] 26/08/2008 16:17:34
If this is indeed the case, switch; the expected gain is
positive because _you already have the information that you hold the
median value of the three possibilities_. The tendency when
presented
with the problem is to reason as if this is the case.
No, you don't know that A is the median. That doesn't even make
sense,
based on the context of the question: there is no 3-valued random
variable here.
This is my point; apologies if I put it badly. The fact is that you
_don't_ hold the median value and that this is indeed a silly way of
looking at it. My assertion was that this is the way many folk DO look
at it, and that this results in an apparent paradox.
In fact, you inadvertently gave an example when you said
The unopened envelope can hold only two values, given
that yours contains A.
True - for a rather restricted meaning of 'true'.
As written, it implicitly allows three values; A for the envelope you
hold, and two more (2A and 0.5A) for the alternatives you permit. The
usual (and incorrect) expected gain calculation uses all three; 2A-A for
one choice and 0.5-A for the other. To do that at all, we must be
playing with three possible values for the contents of an envelope.
This clearly cannot be the case if there are only two possible values,
as we are told in posing the problem. The situation is that you hold
_either_ A _or_ 2A and the other envelope holds (respectively) 2A or A.
We just don't know what A is until we open both envelopes. So for a
given pair of envelopes, it is the choice of coefficient (1 or 2) that
is random.
If I were to describe this in terms of a random variable I would have
to assume an unknown but - for this run - fixed value A multiplied by a
two-valued random variable with possible values 1 and 2, would I not? We
surely can't have both 0.5 and 2 in our distribution at the same time,
because the original proposition said there were only two possibilities
and they differ by a factor of 2, not 4.
You have no basis for putting a probability distribution on those
values, because you don't know the distribution of X.
But we do know; or at least we know the distribution of the
coefficient. We have two envelopes of which one is selected at random,
and the ratio of values is 2. On that basis, assigning a 50:50
probability of ending up with A or 2A on first selection seems
uncontroversial.
But I'm more than willing to have my intuition corrected - possibly
off-line, of course, since this stopped being R a while back!
Steve E
***
This email and any attachments are confidential. Any u...{{dropped:20}}
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] error message in cor.dist
Hello, I am trying to calculate pairwise Pearson correlation distances, and using biodist package, function cor.dist. I start with a table of 4 rows and about 10 columns. (each of 4 samples, or rows have values in each of the 10 categories, no zeros or NAs). I am getting an error message: cor.dist(dmatrixD) Error in cor(t(x)) : missing observations in cov/cor In addition: Warning message: In cor(t(x)) : NAs introduced by coercion Does anyone know what that means? Also, are there other functions for Pearson distances? Thanks! -- Tanya. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] stats tests on large tables
Try t.test(data1[,1],data1[,2]) for the first two columns of data1. --- On Tue, 8/26/08, Richard Emes [EMAIL PROTECTED] wrote: From: Richard Emes [EMAIL PROTECTED] Subject: [R] stats tests on large tables To: r-help@r-project.org Received: Tuesday, August 26, 2008, 11:52 AM I have a large table of data which i can read in using read.table. I then want to conduct various tests on the data. Is there a simple way to conduct these tests by specifying the column headers which relate to different conditions of the experiment? e.g. data1 - read.table(test.table, header=TRUE) t.test(test~control, data = data1) Which doesn't work and results in the error Error in t.test.formula(test ~ control, data = data1) : grouping factor must have exactly 2 levels Many Thanks Richard [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ [[elided Yahoo spam]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problem with Integrate for NEF-HS distribution
Got it. Thanks so much!
Xia
From: [EMAIL PROTECTED]
To: [EMAIL PROTECTED]
CC: r-help@r-project.org
Subject: RE: [R] Problem with Integrate for NEF-HS distribution
Date: Tue, 26 Aug 2008 14:50:18 -0400
Try the following:
sech-function(X)
2/(exp(-X)+exp(X))
your.integrand - function(X)
{0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
my.integrand - function(X) {0.5 *
cos(theta) * 2 / (exp(-pi*X/2 - X*theta) + exp(pi*X/2 - X*theta))
}
theta -
-1
my.integrand(-800)
your.integrand(-800)
Do you see the problem?
Ravi.
---
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: [EMAIL PROTECTED]
Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
From: xia wang [mailto:[EMAIL PROTECTED]
Sent: Tuesday, August 26, 2008 12:59 PM
To: Ravi
Varadhan
Cc: r-help@r-project.org
Subject: RE: [R] Problem
with Integrate for NEF-HS distribution
Thanks so much. It works well in my MCMC sampling. May I know
why re-writing the integrand can solve the problem? I have been thinking it was
the skewness of the distribution brought the error.
Thanks!
Xia
From: [EMAIL PROTECTED]
To:
[EMAIL PROTECTED]; r-help@r-project.org
Subject: RE: [R] Problem with
Integrate for NEF-HS distribution
Date: Tue, 26 Aug 2008 11:59:44
-0400
Hi,
Simply re-write your integrand as
follows:
integrand - function(X) {0.5 * cos(theta)
* 2 / (exp(-pi*X/2 - X*theta) +
exp(pi*X/2 - X*theta)) }
Now the problem goes away.
theta - -1
integrate(integrand, -Inf, 1)
0.9842315 with absolute error
1.2e-05
This would work for all theta such that abs(theta)
-pi/2.
Ravi.
---
Ravi Varadhan, Ph.D.
Assistant
Professor, The Center on Aging and Health
Division of Geriatric
Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: [EMAIL PROTECTED]
Webpage:
-Original Message-
From:
[EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On
Behalf Of xia wang
Sent: Tuesday, August 26, 2008 1:01 AM
To:
r-help@r-project.org
Subject: [R] Problem with Integrate for NEF-HS
distribution
I need to calcuate the
cumulative probability for the Natural Exponential
Family - Hyperbolic
secant distribution with a parameter theta between -pi/2
and pi/2. The
integration should be between 0 and 1 as it is a probability.
The function integrate works fine when the absolute value of theta is
not
too large. That is, the NEF-HS distribution is not too skewed.
However,
once the theta gets large in absolute value, such as -1 as
shown in the
example below, integrate keeps give me error message for
non-finite
function when I put the lower bound as -Inf. I suspect that
it is caused
by the very heavy tail of the distribution.
Is there any way that I can get around of this and let integrate work
for
the skewed distribution? Or is there any other function for
integrating in
R-package? Thanks a lot for your advice in
advance!
theta--1
sech-function(X) 2/(exp(-X)+exp(X))
integrand - function(X)
{0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
integrate(integrand, -3,1)
0.8134389 with absolute error
7e-09
integrate(integrand, -10,1)
0.9810894 with absolute
error 5.9e-06
integrate(integrand, -15,1)
0.9840505
with absolute error 7e-09
integrate(integrand, -50,1)
0.9842315 with absolute error 4.4e-05
integrate(integrand,
-100,1)
0.9842315 with absolute error 3.2e-05
integrate(integrand, -Inf,1)
Error in integrate(integrand, -Inf, 1) :
non-finite function value
Xia
_
Be the
filmmaker you always wanted to be-learn how to burn a DVD with
WindowsR.
__
R-help@r-project.org
mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE
do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
eck It Out!
_
[[alternative HTML version deleted]]
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] How many parameters does my model (gls) have?
Hello, Is there a way to output the number of parameters in my model (gls)? I can count the number of estimates, but I'd like to use the number of parameters as an R object. Thanks, Mark __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How many parameters does my model (gls) have?
?extractAIC 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 On Aug 26, 2008, at 2:27 PM, Mark Na wrote: Hello, Is there a way to output the number of parameters in my model (gls)? I can count the number of estimates, but I'd like to use the number of parameters as an R object. Thanks, Mark __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] R for Windows GUI closes when I try to read.spss
Can't thank you enough! This worked! N'DOYE Souleymane [EMAIL PROTECTED] 8/26/2008 1:20 AM Hi Roderick, Let me suggest you to save your spss file in txt, and use the read.table or read.csv fonction. That is how I proceed. I hope it will help you, Best regards On Mon, Aug 25, 2008 at 7:39 PM, Roderick Harrison [EMAIL PROTECTED] wrote: ** High Priority ** I have been trying to read an SPSS file into R using read.spss (C:/Documents and Settings/Roderick Harrison/My Documents/RWORK/ihisdat.sav, use.value.labels = TRUE, to.data.frame = FALSE, max.value.labels = 500, trim.factor.names = FALSE, trim_values = TRUE, reencode = NA, use.missings = to.data.frame) Each time (at least 5 or 6 by now) I get the following Microsoft error message and the R-Console crashes. R for Windows GUI front-end has encountered a problem and needs to close. We are sorry for the inconvenience. I would greatly appreciate help with this as we need to use R to estimate confidence intervals for NHIS data, and our deliverable is due this week. Thanks in advance. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Souleymane N'DOYE, M.Sc. Statistician Engineer Decision Support System Consultant www.labstatconseil.com [EMAIL PROTECTED] P.O. Box 1601 * 00606 Sarit Center, Nairobi, Kenya Mobile : +254 736 842 478. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Potential Error/Bug?: addition of NULL during arithmetic
I encountered an error that does not make sense to me given my reading of the documentation and does not appear to be referenced in the list archives or online. The error occurred when a function received a NULL value rather than a numeric value as a parameter and then attempted to use the NULL value in the calculateion. The error: Error during wrapup: nothing to replace with can be easily recreated with the following code: 1 + NULL # note that the opperation succeeds and returns a numeric vector with dim(NULL) numeric(0) bar - 1 bar [1] 1 foo - 1 + NULL foo numeric(0) bar - bar + foo # note that here the assignment operation succeeds and 1 + (1 + NULL) - numeric(0) bar numeric(0) bar - c(1, 1) # however if the assignment is into a vector bar[1] - bar[1] + foo # note that the mathematical operation is identical, but the assignment fails Error during wrapup: nothing to replace with If this is the intended behavior, a more informative error message (e.g. 'attempt to assign NULL to vector element') would be useful. If it is not the intended behavior, should I log this as a bug? -eric sessionInfo() R version 2.7.1 (2008-06-23) powerpc-apple-darwin8.10.1 locale: en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8 attached base packages: [1] stats graphics grDevices utils datasets methods base __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Latin squares in R
Is this of any help http://cran.r-project.org/web/views/ExperimentalDesign.html --- On Tue, 8/26/08, Edna Bell [EMAIL PROTECTED] wrote: From: Edna Bell [EMAIL PROTECTED] Subject: [R] Latin squares in R To: [EMAIL PROTECTED] Received: Tuesday, August 26, 2008, 1:55 PM Dear R Gurus: What would be the best way to evaluate a Latin square problem, please? Does it work in Rcmdr, please? thanks, Edna Bell __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] stats tests on large tables
On 27/08/2008, at 3:52 AM, Richard Emes wrote:
I have a large table of data which i can read in using read.table.
I then want to conduct various tests on the data.
Is there a simple way to conduct these tests by specifying the
column headers which relate to different conditions of the experiment?
e.g. data1 - read.table(test.table, header=TRUE)
t.test(test~control, data = data1)
Which doesn't work and results in the error
Error in t.test.formula(test ~ control, data = data1) :
grouping factor must have exactly 2 levels
with(data1,t.test(test ~ control)) # You do not need to specify the
method explicitly.
cheers,
Rolf Turner
##
Attention:\ This e-mail message is privileged and confid...{{dropped:9}}
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] lattice: plotting an arbitrary number of panels, defining arbitrary groups
R Friends,
I'm running R2.7.1 on Windows XP.
I'm trying to get some lattice functionality which I have not seen
previously documented--I'd like to plot the exact same data in multiple
panels but changing the grouping variable each time so that each panel
highlights a different feature of the data set. The following code does
exactly that with a simple and fabricated air quality data set.
dataSet - data.frame(Pollutant=c(rep(Black Carbon,5),rep(PM10,5)),
Detector=c(1:5,1:5), Value=c(seq(50,10,-10),seq(100,60,-10)),
Class=Mass)
xyplot(
Value ~ Detector | Pollutant,
data=dataSet,
aspect = 1.0,
subscripts=TRUE,
panel = function(x,y,subscripts,...) {
if(panel.number() == 1)
panel.superpose(x=dataSet$Detector,y=dataSet$Value,1:nrow(dataSet),groups=dataSet$Pollutant);
if(panel.number() == 2)
panel.superpose(x=dataSet$Detector,y=dataSet$Value,1:nrow(dataSet),groups=normToEdge_dataSet$Class);
}
)
Although the panel labels indicate that only one type of pollutant is
displayed in each, I've instead forced all of the data to be plotted in
both. The first panel shows two colors, grouped by pollutant, the second
shows one color, grouped by class.
Here's where the problem comes, if I add an additional pollutant, instead
defining the data set as follows:
dataSet - data.frame(Pollutant=c(rep(Black
Carbon,5),rep(PM10,5),Ultrafines),
Detector=c(1:5,1:5,10),Value=c(seq(50,10,-10),seq(100,60,-10),75),Class=c(rep(Mass,10),Count))
and rerun the same plotting script, I obtain three panels. The one labeled
Black Carbon correctly displays all three pollutants in different colors.
PM10 however, displays all classes in one color when there should now be
two. Additionally, I now obtain a panel entitled Ultrafines which I'd like
to suppress.
The actual data set has a number of different pollutants, so what I'd
ideally like to do is arbitrarily define two panels with different grouping
variables. I've tried to set up dummy groups and to condition on those, but
with no luck. I think what I need to do is possible with viewports, but is
there no way to entice lattice to function in this way?
Any help would be appreciated.
cheers,
Alex Karner
Department of Civil and Environmental Engineering
University of California, Davis
[[alternative HTML version deleted]]
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] svymeans question
On 26/8/08 9:40 AM, Doran, Harold wrote: Computing the mean of the item by itself with svymeans agrees with the sample mean mean(lower_dat$W524787, na.rm=T) [1] 0.8555471 Compare this to the value in the row 9 up from the bottom to see it is different. You might be omitting more cases due to missing values than you expect. Does the following calculation give you the same results as in rr1? mean( lower_dat$W524787[ apply( lower_dat[lset], 1, function(x) !any(is.na(x)) ) ] ) James -- James Reilly Department of Statistics, University of Auckland Private Bag 92019, Auckland, New Zealand __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] lattice plotting character woes
[Rolf, this crosses with your reply. I will look at your email next.] I pasted the wrong code last time. The following code is supposed to illustrate my problem with lattice plotting character changes. patches - structure(list(URBAN_AREA = structure(c(2L, 19L, 23L, 2L, 19L, 23L, 2L, 19L, 23L, 2L, 19L, 23L), .Label = c(CENTRAL AUCKLAND ZONE, CHRISTCHURCH, DUNEDIN, HAMILTON ZONE, HASTINGS ZONE, INVERCARGILL, LOWER HUTT ZONE, mean, NAPIER ZONE, NELSON, NEW PLYMOUTH, NORTHERN AUCKLAND ZONE, PALMERSTON NORTH, PORIRUA ZONE, ROTORUA, SD, SE, SOUTHERN AUCKLAND ZONE, TAURANGA, WANGANUI, WELLINGTON ZONE, WESTERN AUCKLAND ZONE, WHANGAREI), class = factor), NO_PATCHES = c(11L, 16L, 21L, 87L, 192L, 324L, 164L, 417L, 773L, 679L, 757L, 3083L), MEAN_AREA = c(9.623631225, 15.29089619, 149.2063532, 14.1676, 247.5262, 28.611, 11.5698, 221.0022, 37.3725, 11.918, 133.5804, 25.6759), AREA.ZONE = c(13683, 3666, 1558, 64830, 41103, 22581, 123819, 90107, 57627, 264735, 223963, 174456), Buffer.zone = c(0L, 0L, 0L, 5L, 5L, 5L, 10L, 10L, 10L, 20L, 20L, 20L)), .Names = c(URBAN_AREA, NO_PATCHES, MEAN_AREA, AREA.ZONE, Buffer.zone), class = data.frame, row.names = c(2L, 15L, 19L, 22L, 36L, 40L, 42L, 56L, 60L, 62L, 76L, 80L)) library(lattice) Region = factor(patches$URBAN_AREA) lpatches = patches for (i in 2:4) lpatches[,i] = log10(patches[,i]) # zone not transformed lpatches.pca = princomp(lpatches[,-c(1,17)],cor=FALSE) x = lpatches.pca$scores[,1] y = lpatches.pca$scores[,2] zz = as.character(patches$Buffer.zone/5) table(zz) plsy - trellis.par.get(plot.symbol) # only 0 or 1 used as plotting symbol # I expected 0,1,2,4 plsy$pch = as.character(rep(1:6,2)) trellis.par.set(plot.symbol,plsy) xyplot(y ~ x |Region) # only 1,2,3,4 used as plotting symbol # I expected 1:6 Mark Leeds has pointed out that whatever numbers appear as plotting characters on the screen, they are replaced by circles when you save a pdf via pdf() xyplot() dev.off() I have reproduced this on my Mac. Cheers, Murray -- Dr Murray Jorgensen http://www.stats.waikato.ac.nz/Staff/maj.html Department of Statistics, University of Waikato, Hamilton, New Zealand Email: [EMAIL PROTECTED][EMAIL PROTECTED] Fax 7 838 4155 Phone +64 7 838 4773 wkHome +64 7 825 0441Mobile 021 139 5862 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] loop with splitted data
Is this what you want:
y = c(rep(2,5),rep(3,5))
ma - data.frame(x = 1:10, y=y )
splitted - split(ma, ma$y)
for (counter in (min(ma$y):max(ma$y)))
+ {
+ cat(counter, :, splitted[[as.character(counter)]]$x, '\n')
+ }
2 : 1 2 3 4 5
3 : 6 7 8 9 10
On Tue, Aug 26, 2008 at 6:37 AM, Knut Krueger [EMAIL PROTECTED] wrote:
Hi to all,
seems to be simple, but I do not find the solution:
What must I write for the splitted to get
splitted$3$x and splitted$3$x
y = c(rep(2,5),rep(3,5))
ma - data.frame(x = 1:10, y=y )
splitted - split(ma, ma$y)
for (counter in (min(ma$y):max(ma$y)))
{
splitted$x
}
Regards Knut
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
--
Jim Holtman
Cincinnati, OH
+1 513 646 9390
What is the problem that you are trying to solve?
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] lattice: plotting an arbitrary number of panels, defining arbitrary groups
On Tue, Aug 26, 2008 at 2:26 PM, Alex Karner [EMAIL PROTECTED] wrote:
R Friends,
I'm running R2.7.1 on Windows XP.
I'm trying to get some lattice functionality which I have not seen
previously documented--I'd like to plot the exact same data in multiple
panels but changing the grouping variable each time so that each panel
highlights a different feature of the data set. The following code does
exactly that with a simple and fabricated air quality data set.
dataSet - data.frame(Pollutant=c(rep(Black Carbon,5),rep(PM10,5)),
Detector=c(1:5,1:5), Value=c(seq(50,10,-10),seq(100,60,-10)),
Class=Mass)
xyplot(
Value ~ Detector | Pollutant,
data=dataSet,
aspect = 1.0,
subscripts=TRUE,
panel = function(x,y,subscripts,...) {
if(panel.number() == 1)
panel.superpose(x=dataSet$Detector,y=dataSet$Value,1:nrow(dataSet),groups=dataSet$Pollutant);
if(panel.number() == 2)
panel.superpose(x=dataSet$Detector,y=dataSet$Value,1:nrow(dataSet),groups=normToEdge_dataSet$Class);
}
)
Although the panel labels indicate that only one type of pollutant is
displayed in each, I've instead forced all of the data to be plotted in
both. The first panel shows two colors, grouped by pollutant, the second
shows one color, grouped by class.
Here's where the problem comes, if I add an additional pollutant, instead
defining the data set as follows:
dataSet - data.frame(Pollutant=c(rep(Black
Carbon,5),rep(PM10,5),Ultrafines),
Detector=c(1:5,1:5,10),Value=c(seq(50,10,-10),seq(100,60,-10),75),Class=c(rep(Mass,10),Count))
and rerun the same plotting script, I obtain three panels. The one labeled
Black Carbon correctly displays all three pollutants in different colors.
PM10 however, displays all classes in one color when there should now be
two. Additionally, I now obtain a panel entitled Ultrafines which I'd like
to suppress.
The actual data set has a number of different pollutants, so what I'd
ideally like to do is arbitrarily define two panels with different grouping
variables. I've tried to set up dummy groups and to condition on those, but
with no luck. I think what I need to do is possible with viewports, but is
there no way to entice lattice to function in this way?
Any help would be appreciated.
Panels can be repeated, using the standard R indexing interface; so,
for example,
## trellis object with one panel, but different groups depending on
panel.number()
p -
with(dataSet,
xyplot(Value ~ Detector,
group.list = list(Pollutant, Class),
aspect = 1.0,
subscripts = TRUE,
panel = function(..., group.list) {
panel.xyplot(...,
groups = group.list[[panel.number()]])
}))
## plot first panel twice
p[c(1, 1)]
## add a strip function
update(p[c(1, 1)],
strip = function(...) {
lab - c(Pollutant, Class)[panel.number()]
strip.default(1, 1, var.name = , factor.levels = lab)
})
-Deepayan
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] svymeans question
Indeed, missing data are the problem. The function removes any row with missing data and the means are based on the remaining rows. So, I wrote a function that just loops over each variable individually and organizes the data as I need it. -Original Message- From: James Reilly on behalf of James Reilly Sent: Tue 8/26/2008 5:39 PM To: Doran, Harold Cc: r-help@r-project.org Subject: Re: [R] svymeans question On 26/8/08 9:40 AM, Doran, Harold wrote: Computing the mean of the item by itself with svymeans agrees with the sample mean mean(lower_dat$W524787, na.rm=T) [1] 0.8555471 Compare this to the value in the row 9 up from the bottom to see it is different. You might be omitting more cases due to missing values than you expect. Does the following calculation give you the same results as in rr1? mean( lower_dat$W524787[ apply( lower_dat[lset], 1, function(x) !any(is.na(x)) ) ] ) James -- James Reilly Department of Statistics, University of Auckland Private Bag 92019, Auckland, New Zealand [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] loop with splitted data
jim holtman wrote:
Is this what you want:
y = c(rep(2,5),rep(3,5))
ma - data.frame(x = 1:10, y=y )
splitted - split(ma, ma$y)
for (counter in (min(ma$y):max(ma$y)))
+ {
+ cat(counter, :, splitted[[as.character(counter)]]$x, '\n')
+ }
2 : 1 2 3 4 5
3 : 6 7 8 9 10
But maybe this is what he really wanted:
lapply(splitted,[[, x)
$`2`
[1] 1 2 3 4 5
$`3`
[1] 6 7 8 9 10
or even
split(ma$x, ma$y)
$`2`
[1] 1 2 3 4 5
$`3`
[1] 6 7 8 9 10
--
O__ Peter Dalgaard Øster Farimagsgade 5, Entr.B
c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
(*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918
~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] awkward behavior with densityplot function
Hi,
I have the following script:
t.R ---
grafica - function() {
v - read.csv('preprocessed/komolongma.ece.uprm.edu.active',sep=',')
x - as.ts(v$active)
bitmap(file=output.png)
densityplot(~x,col='blue',main='Density Plot')
dev.off()
}
grafica()
t.R ---
When I sourced it from R prompt, it quietly runs. However the
output.png generated contains any visible data. I said that, because
the file is created and it has 2062 bytes.
--- execution ---
library(lattice)
source(file=t.R)
--- execution ---
Now, if I sequentially run the lines above in the R prompt, the
output.png file contains a graphical representation of the data.
--- execution 2
v - read.csv('preprocessed/komolongma.ece.uprm.edu.active',sep=',')
x - as.ts(v$active)
bitmap(file=output.png)
densityplot(~x,col='blue',main='Density Plot')
dev.off()
execution 2 ---
What is wrong with the densityplot function that produces any output
when is invoked from a script?
Someone has a script invoking the densityplot function?
Thanks a lot for your help.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] awkward behavior with densityplot function
Lattice functions do not create plots. In this respect they are quite
different to functions like plot(...), hist(...) c, which do create plots.
If you want to see the plot from a lattice function, you need to print the
object it creates. It's the action of printing that creates the graphical
image.
This happens automatically at the command line, but if you want it to happen
inside the function you need to use print(...) explicitly.
So your function should be as in the following:
t.R ---
grafica - function() {
v - read.csv('preprocessed/komolongma.ece.uprm.edu.active',sep=',')
x - as.ts(v$active)
bitmap(file=output.png)
print(densityplot(~x,col='blue',main='Density Plot')) ###---
dev.off()
invisible(v) ### further suggestion.
}
grafica()
t.R ---
Bill Venables
http://www.cmis.csiro.au/bill.venables/
-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On Behalf Of John Sanabria
Sent: Wednesday, 27 August 2008 9:08 AM
To: r-help@r-project.org
Subject: [R] awkward behavior with densityplot function
Hi,
I have the following script:
t.R ---
grafica - function() {
v - read.csv('preprocessed/komolongma.ece.uprm.edu.active',sep=',')
x - as.ts(v$active)
bitmap(file=output.png)
densityplot(~x,col='blue',main='Density Plot')
dev.off()
}
grafica()
t.R ---
When I sourced it from R prompt, it quietly runs. However the
output.png generated contains any visible data. I said that, because
the file is created and it has 2062 bytes.
--- execution ---
library(lattice)
source(file=t.R)
--- execution ---
Now, if I sequentially run the lines above in the R prompt, the
output.png file contains a graphical representation of the data.
--- execution 2
v - read.csv('preprocessed/komolongma.ece.uprm.edu.active',sep=',')
x - as.ts(v$active)
bitmap(file=output.png)
densityplot(~x,col='blue',main='Density Plot')
dev.off()
execution 2 ---
What is wrong with the densityplot function that produces any output
when is invoked from a script?
Someone has a script invoking the densityplot function?
Thanks a lot for your help.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] awkward behavior with densityplot function
See FAQ 7.22
On 27/08/2008, at 11:07 AM, John Sanabria wrote:
Hi,
I have the following script:
t.R ---
grafica - function() {
v - read.csv('preprocessed/komolongma.ece.uprm.edu.active',sep=',')
x - as.ts(v$active)
bitmap(file=output.png)
densityplot(~x,col='blue',main='Density Plot')
dev.off()
}
grafica()
t.R ---
When I sourced it from R prompt, it quietly runs. However the
output.png generated contains any visible data. I said that,
because the file is created and it has 2062 bytes.
--- execution ---
library(lattice)
source(file=t.R)
--- execution ---
Now, if I sequentially run the lines above in the R prompt, the
output.png file contains a graphical representation of the data.
--- execution 2
v - read.csv('preprocessed/komolongma.ece.uprm.edu.active',sep=',')
x - as.ts(v$active)
bitmap(file=output.png)
densityplot(~x,col='blue',main='Density Plot')
dev.off()
execution 2 ---
What is wrong with the densityplot function that produces any
output when is invoked from a script?
Someone has a script invoking the densityplot function?
##
Attention:\ This e-mail message is privileged and confid...{{dropped:9}}
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
[R] A Tip: lm, glm, and retained cases
Hi Folks, This tip is probably lurking somewhere already, but I've just discovered it the hard way, so it is probably worth passing on for the benefit of those who might otherwise hack their way along the same path. Say (for example) you want to do a logistic regression of a binary response Y on variables X1, X2, X3, X4: GLM - glm(Y ~ X1 + X2 + X3 + X4) Say there are 1000 cases in the data. Because of missing values (NAs) in the variables, the number of complete cases retained for the regression is, say, 600. glm() does this automatically. QUESTION: Which cases are they? You can of course find out by hand on the lines of ix - which( (!is.na(Y))(!is.na(X1))...(!is.na(X4)) ) but one feels that GLM already knows -- so how to get it to talk? ANSWER: (e.g.) ix - as.integer(names(GLM$fit)) Reason: When glm(Y~X1+...) picks up the data passed to it, it assigns[*] to each element of Y a name which is its integer position in the variable, expressed as a character string (1, 2, 3, ... ). [*] Assuming (as is usually the case) that the elements didn't have names in the first place. Otherwise these names are used; modify the above approach accordingly. These names are retained during the computation, and when incomplete cases are dropped the retained complete cases retain their original names. Thus, any per-case series of computed values (such as $fit) has the names of the retained cases the values correspond to. These can be discovered by names(GLM$fit) but you don't want them as character strings, so convert them to integers: as.integer(names(GLM$fit)) Done! I hope this helps some people. Ted. E-Mail: (Ted Harding) [EMAIL PROTECTED] Fax-to-email: +44 (0)870 094 0861 Date: 27-Aug-08 Time: 00:45:47 -- XFMail -- __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] A Tip: lm, glm, and retained cases
On Tue, Aug 26, 2008 at 6:45 PM, Ted Harding [EMAIL PROTECTED] wrote: Hi Folks, This tip is probably lurking somewhere already, but I've just discovered it the hard way, so it is probably worth passing on for the benefit of those who might otherwise hack their way along the same path. Say (for example) you want to do a logistic regression of a binary response Y on variables X1, X2, X3, X4: GLM - glm(Y ~ X1 + X2 + X3 + X4) Say there are 1000 cases in the data. Because of missing values (NAs) in the variables, the number of complete cases retained for the regression is, say, 600. glm() does this automatically. QUESTION: Which cases are they? You can of course find out by hand on the lines of ix - which( (!is.na(Y))(!is.na(X1))...(!is.na(X4)) ) but one feels that GLM already knows -- so how to get it to talk? ANSWER: (e.g.) ix - as.integer(names(GLM$fit)) Alternatively, you can use: attr(GLM$model, na.action) Hadley -- http://had.co.nz/ __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Problem with Integrate for NEF-HS distribution
If you look at your sech(pi*x/2) you can write it as
sech(pi*x/2) = 2*exp(pi*x/2)/(1 + exp(pi*x))
For x -15, exp(pi*x) 10^-20, so for this interval you can replace
sech(pi*x/2) by 2*exp(pi*x/2) and so the integral from -Inf to -15 (or even -10
- depends on your accuracy requirements) can be computed analytically, so you
are left with integral from -10 (or -15) to your upper bound and this should be
all right for numerical integration.
--- On Wed, 27/8/08, xia wang [EMAIL PROTECTED] wrote:
From: xia wang [EMAIL PROTECTED]
Subject: RE: [R] Problem with Integrate for NEF-HS distribution
To: [EMAIL PROTECTED], [EMAIL PROTECTED]
Received: Wednesday, 27 August, 2008, 12:26 AM
Thanks. I revised the code but it doesn't make
difference. The cause of the error is just as stated in the
?integrate If the function is approximately constant
(in particular, zero) over nearly all
its range it is possible that the result and error estimate
may be seriously
wrong. I have searched R-archive. It may not be
feasible to solve the problem by rectangle
approximation as some postings suggested because the
integration is in fact within MCMC samplings as follows.
The lower bound is always -Inf. The upper bound and the
value of theta changes for each sampling in MCMC.
I tried to multiple the integrand by a large number ( like
10^16) and changes the rel.tol. It does help for some range
of theta but not all.
Xia
like.fun - function(beta,theta) {
integrand - function(X,theta)
{0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
upper-x%*%beta
prob-matrix(NA, nrow(covariate),1)
for (i in 1:nrow(covariate))
{prob[i]-integrate(integrand,lower=-Inf,
upper=upper[i],theta=theta)$value
}
likelihood - sum(log(dbinom(y,n,prob)))
return(likelihood)
}
Date: Tue, 26 Aug 2008 00:49:02 -0700
From: [EMAIL PROTECTED]
Subject: Re: [R] Problem with Integrate for NEF-HS
distribution
To: [EMAIL PROTECTED]
Shouldn't you define
integrand - function(X,theta) ..
and not
integrand - function(X) .
--- On Tue, 26/8/08, xia wang
[EMAIL PROTECTED] wrote:
From: xia wang [EMAIL PROTECTED]
Subject: [R] Problem with Integrate for NEF-HS
distribution
To: r-help@r-project.org
Received: Tuesday, 26 August, 2008, 3:00 PM
I need to calcuate the cumulative probability for
the
Natural Exponential Family - Hyperbolic secant
distribution
with a parameter theta between -pi/2 and pi/2.
The
integration should be between 0 and 1 as it is a
probability.
The function integrate works fine
when the
absolute value of theta is not too large. That
is, the
NEF-HS distribution is not too skewed. However,
once the
theta gets large in absolute value, such as -1 as
shown in
the example below, integrate keeps
give me error
message for non-finite function when
I put the
lower bound as -Inf. I suspect that it is caused
by the
very heavy tail of the distribution.
Is there any way that I can get around of this
and let
integrate work for the skewed
distribution? Or
is there any other function for integrating in
R-package?
Thanks a lot for your advice in advance!
theta--1
sech-function(X) 2/(exp(-X)+exp(X))
integrand - function(X)
{0.5*cos(theta)*exp(X*theta)*sech(pi*X/2)}
integrate(integrand, -3,1)
0.8134389 with absolute error 7e-09
integrate(integrand, -10,1)
0.9810894 with absolute error 5.9e-06
integrate(integrand, -15,1)
0.9840505 with absolute error 7e-09
integrate(integrand, -50,1)
0.9842315 with absolute error 4.4e-05
integrate(integrand, -100,1)
0.9842315 with absolute error 3.2e-05
integrate(integrand, -Inf,1)
Error in integrate(integrand, -Inf, 1) :
non-finite
function value
Xia
_
Be the filmmaker you always wanted to be—learn
how to
burn a DVD with Windows®.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained,
reproducible code.
_
Be the filmmaker you always wanted to be—learn how to
burn a DVD with Windows®.
http://clk.atdmt.com/MRT/go/108588797/direct/01/
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] A Tip: lm, glm, and retained cases
On 26-Aug-08 23:49:37, hadley wickham wrote: On Tue, Aug 26, 2008 at 6:45 PM, Ted Harding [EMAIL PROTECTED] wrote: Hi Folks, This tip is probably lurking somewhere already, but I've just discovered it the hard way, so it is probably worth passing on for the benefit of those who might otherwise hack their way along the same path. Say (for example) you want to do a logistic regression of a binary response Y on variables X1, X2, X3, X4: GLM - glm(Y ~ X1 + X2 + X3 + X4) Say there are 1000 cases in the data. Because of missing values (NAs) in the variables, the number of complete cases retained for the regression is, say, 600. glm() does this automatically. QUESTION: Which cases are they? You can of course find out by hand on the lines of ix - which( (!is.na(Y))(!is.na(X1))...(!is.na(X4)) ) but one feels that GLM already knows -- so how to get it to talk? ANSWER: (e.g.) ix - as.integer(names(GLM$fit)) Alternatively, you can use: attr(GLM$model, na.action) Hadley Thanks! I can see that it works -- though understanding how requires a deeper knowledge of R internals. However, since you've approached it from that direction, simply GLM$model is a dataframe of the retained cases (with corresponding row-names), all variables at once, and that is possibly an even simpler approach! Ted. E-Mail: (Ted Harding) [EMAIL PROTECTED] Fax-to-email: +44 (0)870 094 0861 Date: 27-Aug-08 Time: 01:31:46 -- XFMail -- __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Do I need a special package to run LSD tests?
Thanks for the answers, this one worked: library(Agricolae) michal33 wrote: Hi there, I am trying to run LSD.test(model) I used the following commands: attach(model) m1- glm(ttl.m ~ site+year, family=quasipoisson, data= model) df-df.residual(m1) MSerror-deviance(m1)/df The following command did not work: comparison- LSD.test (ttl.m, site+year, df, MSerror, alpha = 0.05, group=FALSE) I get an error message: Error: could not find function LSD.test Do I need to download a special package for that? Any idea which one? or how to do it without getting the error message? All comment will be highly appreciated! Thanks :) Michal -- View this message in context: http://www.nabble.com/Do-I-need-a-special-package-to-run-LSD-tests--tp19154567p19172843.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Do I need a special package to run LSD tests?
On 27/08/2008, at 12:32 PM, michal33 wrote:
Thanks for the answers, this one worked:
library(Agricolae)
No, it couldn't have. There is no such package as ``Agricolae''.
There ***is*** a package called ``agricolae''. :-)
cheers,
Rolf Turner
##
Attention:\ This e-mail message is privileged and confid...{{dropped:9}}
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Re: [R] A Tip: lm, glm, and retained cases
on 08/26/2008 07:31 PM (Ted Harding) wrote: On 26-Aug-08 23:49:37, hadley wickham wrote: On Tue, Aug 26, 2008 at 6:45 PM, Ted Harding [EMAIL PROTECTED] wrote: Hi Folks, This tip is probably lurking somewhere already, but I've just discovered it the hard way, so it is probably worth passing on for the benefit of those who might otherwise hack their way along the same path. Say (for example) you want to do a logistic regression of a binary response Y on variables X1, X2, X3, X4: GLM - glm(Y ~ X1 + X2 + X3 + X4) Say there are 1000 cases in the data. Because of missing values (NAs) in the variables, the number of complete cases retained for the regression is, say, 600. glm() does this automatically. QUESTION: Which cases are they? You can of course find out by hand on the lines of ix - which( (!is.na(Y))(!is.na(X1))...(!is.na(X4)) ) but one feels that GLM already knows -- so how to get it to talk? ANSWER: (e.g.) ix - as.integer(names(GLM$fit)) Alternatively, you can use: attr(GLM$model, na.action) Hadley Thanks! I can see that it works -- though understanding how requires a deeper knowledge of R internals. However, since you've approached it from that direction, simply GLM$model is a dataframe of the retained cases (with corresponding row-names), all variables at once, and that is possibly an even simpler approach! Or just use: model.frame(ModelObject) as the extractor function... :-) Another 'a priori' approach would be to use na.omit() or one of its brethren on the data frame before creating the model. Which function is used depends upon how 'na.action' is set. The returned value, or more specifically the 'na.action' attribute as appropriate, would yield information similar to Hadley's approach relative to which records were excluded. For example, using the simple data frame in ?na.omit: DF - data.frame(x = c(1, 2, 3), y = c(0, 10, NA)) DF x y 1 1 0 2 2 10 3 3 NA DF.na - na.omit(DF) DF.na x y 1 1 0 2 2 10 attr(DF.na, na.action) 3 3 attr(,class) [1] omit So you can see that record 3 was removed from the original data frame due to the NA for 'y'. HTH, Marc Schwartz __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
