[R] Splom plot:how to plot with different symbols?

[d. sarthi maheshwari]

Hi,

Kindly let me know if I posted a wrong question in the forum.

I want to draw a splom plot with different symbols in plot. My command is as
follows:

splom(~ log10(splomData[2:3]), groups = programs, data = splomData,
panel = panel.superpose,
key = list(title = paste(splomLoop,"Programs of Hog Analysis (Sorted by
LR(GB))"),
   panel = panel.superpose,
 columns = splomLoop,
   points = list(pch = c(1:splomLoop), #this will display different
symbols in legend but not in actual graph
 col = rainbow(splomLoop, s=1, v=1),
   text = trim(as.character(ProgramNameList[1:splomLoop,1])

Kindly help.

Thanks & Regards
Sarthi M.

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Different gridlines per panel in xyplot

[Rene Braeckman]

I solved my own problem. Suddenly remembered the list of very useful
functions under "Accessing Auxiliary Information During Plotting" in the
help pages.

Here is the line that did the trick:

panel.abline(h=yScale[[panel.number()]], v=xScale, col.line="gray")

Rene 

-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Rene Braeckman
Sent: Tuesday, February 20, 2007 8:52 PM
To: r-help@stat.math.ethz.ch
Subject: [R] Different gridlines per panel in xyplot

In the example R script below, horizontal gray gridlines are drawn at y
coordinates where the points are drawn with the code:
 
panel.abline(h=y, v=xScale, col.line="gray")

How do I change this so that the horizontal gray gridlines are drawn at y
coordinates where the y labels are drawn? The challenge is that each panel
has different y-ranges (in my real example the y-ranges and y-intervals are
even more different). For example, I wish I could use the yScale list as the
h parameter in abline, but it does not work with a list.
 
Thanks for any help.
Rene
 
library(lattice)
Subj <- rep(1:4,each=3)
Time <- rep(1:3,4) + 0.1
Conc <- (1:12) + 0.1
df <- data.frame(Subj,Time,Conc)
xScale <- 1:3
yScale <- list(1:3,4:6,7:9,10:12)
xyplot(Conc ~ Time | Subj,
   data = df,
   layout = c(2,2),
   type="b",
   scales=list(
  x=list(at=xScale),
  y=list(at=yScale,relation="free")
   ),
   panel = function(x,y,...) {
   panel.abline(h=y, v=xScale, col.line="gray")
   panel.xyplot(x,y,...)
  }
)

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Different gridlines per panel in xyplot

[Rene Braeckman]

Thanks for the reply. Trunc would only work when the truncated y values
result in the desired y coordinates for the grid lines as in this simple
data set. Since my real dataset contains many more points, it does not give
a general solution.

Rene  

-Original Message-
From: Bernd Weiss [mailto:[EMAIL PROTECTED] 
Sent: Tuesday, February 20, 2007 9:44 PM
To: Rene Braeckman; r-help@stat.math.ethz.ch
Subject: Re: [R] Different gridlines per panel in xyplot

Am 20 Feb 2007 um 20:52 hat Rene Braeckman geschrieben:

From:   "Rene Braeckman" <[EMAIL PROTECTED]>
To: 
Date sent:  Tue, 20 Feb 2007 20:52:29 -0800
Subject:[R] Different gridlines per panel in xyplot

> In the example R script below, horizontal gray gridlines are drawn at 
> y coordinates where the points are drawn with the code:
> 
> panel.abline(h=y, v=xScale, col.line="gray")
> 
> How do I change this so that the horizontal gray gridlines are drawn 
> at y coordinates where the y labels are drawn? The challenge is that 
> each panel has different y-ranges (in my real example the y-ranges and 
> y-intervals are even more different). For example, I wish I could use 
> the yScale list as the h parameter in abline, but it does not work 
> with a list.
> 
> Thanks for any help.
> Rene
> 
> library(lattice)
> Subj <- rep(1:4,each=3)
> Time <- rep(1:3,4) + 0.1
> Conc <- (1:12) + 0.1
> df <- data.frame(Subj,Time,Conc)
> xScale <- 1:3
> yScale <- list(1:3,4:6,7:9,10:12)
> xyplot(Conc ~ Time | Subj,
>data = df,
>layout = c(2,2),
>type="b",
>scales=list(
>   x=list(at=xScale),
>   y=list(at=yScale,relation="free")
>),
>panel = function(x,y,...) {
>panel.abline(h=y, v=xScale, col.line="gray")
>panel.xyplot(x,y,...)
>   }
> )
> 

Dear Rene,

I am not quite sure whether this is the most elegant/general solution, but
one option might be to use trunc(), see the full code below.

panel.abline(h=trunc(y), v=xScale, col.line="gray")


library(lattice)
Subj <- rep(1:4,each=3)
Time <- rep(1:3,4) + 0.1
Conc <- (1:12) + 0.1
df <- data.frame(Subj,Time,Conc)
xScale <- 1:3
yScale <- list(1:3,4:6,7:9,10:12)
xyplot(Conc ~ Time | Subj,
   data = df,
   layout = c(2,2),
   type="b",
   scales=list(
  x=list(at=xScale),
  y=list(at=yScale,relation="free")
   ),
   panel = function(x,y,...) {
## use trunc(y)
   panel.abline(h=trunc(y), v=xScale, col.line="gray")
   panel.xyplot(x,y,...)
  }
)


HTH,

Bernd

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Different gridlines per panel in xyplot

[Bernd Weiss]

Am 20 Feb 2007 um 20:52 hat Rene Braeckman geschrieben:

From:   "Rene Braeckman" <[EMAIL PROTECTED]>
To: 
Date sent:  Tue, 20 Feb 2007 20:52:29 -0800
Subject:[R] Different gridlines per panel in xyplot

> In the example R script below, horizontal gray gridlines are drawn at
> y coordinates where the points are drawn with the code:
> 
> panel.abline(h=y, v=xScale, col.line="gray")
> 
> How do I change this so that the horizontal gray gridlines are drawn
> at y coordinates where the y labels are drawn? The challenge is that
> each panel has different y-ranges (in my real example the y-ranges and
> y-intervals are even more different). For example, I wish I could use
> the yScale list as the h parameter in abline, but it does not work
> with a list.
> 
> Thanks for any help.
> Rene
> 
> library(lattice)
> Subj <- rep(1:4,each=3)
> Time <- rep(1:3,4) + 0.1
> Conc <- (1:12) + 0.1
> df <- data.frame(Subj,Time,Conc)
> xScale <- 1:3
> yScale <- list(1:3,4:6,7:9,10:12)
> xyplot(Conc ~ Time | Subj,
>data = df,
>layout = c(2,2),
>type="b",
>scales=list(
>   x=list(at=xScale),
>   y=list(at=yScale,relation="free")
>),
>panel = function(x,y,...) {
>panel.abline(h=y, v=xScale, col.line="gray")
>panel.xyplot(x,y,...)
>   }
> )
> 

Dear Rene,

I am not quite sure whether this is the most elegant/general 
solution, but one option might be to use trunc(), see the full code 
below.

panel.abline(h=trunc(y), v=xScale, col.line="gray")


library(lattice)
Subj <- rep(1:4,each=3)
Time <- rep(1:3,4) + 0.1
Conc <- (1:12) + 0.1
df <- data.frame(Subj,Time,Conc)
xScale <- 1:3
yScale <- list(1:3,4:6,7:9,10:12)
xyplot(Conc ~ Time | Subj,
   data = df,
   layout = c(2,2),
   type="b",
   scales=list(
  x=list(at=xScale),
  y=list(at=yScale,relation="free")
   ),
   panel = function(x,y,...) {
## use trunc(y)
   panel.abline(h=trunc(y), v=xScale, col.line="gray")
   panel.xyplot(x,y,...)
  }
)


HTH,

Bernd

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Different gridlines per panel in xyplot

[Rene Braeckman]

In the example R script below, horizontal gray gridlines are drawn at y
coordinates where the points are drawn with the code:
 
panel.abline(h=y, v=xScale, col.line="gray")

How do I change this so that the horizontal gray gridlines are drawn at y
coordinates where the y labels are drawn? The challenge is that each panel
has different y-ranges (in my real example the y-ranges and y-intervals are
even more different). For example, I wish I could use the yScale list as the
h parameter in abline, but it does not work with a list.
 
Thanks for any help.
Rene
 
library(lattice)
Subj <- rep(1:4,each=3)
Time <- rep(1:3,4) + 0.1
Conc <- (1:12) + 0.1
df <- data.frame(Subj,Time,Conc)
xScale <- 1:3
yScale <- list(1:3,4:6,7:9,10:12)
xyplot(Conc ~ Time | Subj,
   data = df,
   layout = c(2,2),
   type="b",
   scales=list(
  x=list(at=xScale),
  y=list(at=yScale,relation="free")
   ),
   panel = function(x,y,...) {
   panel.abline(h=y, v=xScale, col.line="gray")
   panel.xyplot(x,y,...)
  }
)

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Installing Package rgl - Compilation Fails - FreeBSD

[Duncan Murdoch]

On 2/20/2007 4:45 PM, Prof Brian Ripley wrote:
> The problem is that rgl is apparently written for GNU make, and has (as 
> shipped)
> 
> ifdef MAKINGAGL
> [EMAIL PROTECTED]@ -Iext
> [EMAIL PROTECTED]@
> else
> PKG_CPPFLAGS= -If:/R/R-2.4.1/src/extra/zlib -DHAVE_PNG_H 
> -If:/R/R-2.4.1/src/gnuwin32/bitmap/libpng  -Iext
> PKG_LIBS=-lgdi32 -lopengl32 -lglu32 
> -Lf:/R/R-2.4.1/src/gnuwin32/bitmap/libpng -
> lpng -Lf:/R/R-2.4.1/src/extra/zlib -lz
> endif
> 
> and similar for BUILDAGL.
> 
> That seems to have been written to make it workable on MacOS X. Given that 
> configure knows (or could know) the OS, it seems better to write (via 
> configure) a separate Makevars for MacOS X and remove all the 
> ifdef...endif stuff for everyone else.  (If you do that in Makevars.in it 
> should work.)

On MacOS X we currently build two libraries, one for X11 and one for 
AGL, and these tests choose between those two targets (which need 
different compilation of a few files, and separate linking of all files 
against different libraries).

Fixing this to work more portably isn't something that I know how to do. 
   If someone wants to fix it and send me a patch I'll incorporate it, 
but otherwise it likely won't get fixed.

By the way, the file you quote is Makevars, which shouldn't have been 
shipped:  it gets produced from Makevars.in.  I'll have to be careful to 
do the packaging from a clean checkout next time.

Duncan Murdoch

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Simplification of Generalised Linear mixed effects models using glmmPQL

[Andrew Robinson]

Hello Tom,

the problem is because R has assumed that pop and rep are integers,
not factor levels.  Try:

test <- read.table("test.txt",header=T)

sapply(test, class)

test$pop <- factor(test$pop)
test$rep <- factor(test$rep)

then try fitting the models.  Also, there has been substantial
discussion about the production of p-values for mixed-effects models
in R; it's now a FAQ:

http://cran.r-project.org/doc/FAQ/R-FAQ.html#Why-are-p_002dvalues-not-displayed-when-using-lmer_0028_0029_003f

The package writer (Professor Bates) recommends the use of mcmcsamp to
obtain inferential information about your model - see

?mcmcsamp

in the lme4 package for examples of its use.

Cheers,

Andrew

ps it's a bad idea to call things "data" :)

On Tue, Feb 20, 2007 at 01:13:25PM +, T.C. Cameron wrote:
> Dear R users I have built several glmm models using glmmPQL in the
> following structure:
>  
> m1<-glmmPQL(dev~env*har*treat+dens, random = ~1|pop/rep,  family =
> Gamma)
>  
> (full script below, data attached)
>  
> I have tried all the methods I can find to obtain some sort of model fit
> score or to compare between models using following the deletion of terms
> (i.e. AIC, logLik, anova.lme(m1,m2)), but I cannot get any of them to
> work.
> Yet I see on several R help pages that others have with similar models?
>  
> I have tried the functions in lme4 as well and lmer or lmer2 will not
> accept my random terms of "rep" (replicate) nested within "pop"
> population.
>  
> I have read the appropriate sections of the available books and R help
> pages but I am at a loss of where to move from here
>  
>  
> 
> data<-read.table("D:\\bgytcc\\MITES\\Data\\Analysis\\test.txt",header=T)
> attach(data)
> names(data)
>  
>  
>  
> m1<-glmmPQL(dev~env*har*treat+dens, random = ~1|pop/rep, family = Gamma)
> summary(m1)
> anova.lme(m1)
> m2<-update(m1,~.-env:har:treat)
> anova.lme(m1,m2)###this does not work
> AIC(m1)##this does not work
> logLik(m1)##this does not work?
>  
>  
>  
> ##this does not work
> class(m1) <- "lme"
> class(m2) <- "lme"
> anova.lme(m1,m2)
> #
>  
> m3<-lmer(dev~env*har*treat+dens + (1|pop/rep), family = Gamma)
>  
> ## this generates an error
> Error in lmerFactorList(formula, mf, fltype) : 
> number of levels in grouping factor(s) 'rep:pop', 'pop' is too
> large
> In addition: Warning messages:
> 1: numerical expression has 1851 elements: only the first used in:
> rep:pop 
> 2: numerical expression has 1851 elements: only the first used in:
> rep:pop 
>  
> 
> m4<-lmer(dev~env*har*treat + dens + (1|rep) +(1|pop), family = Gamma,
> method = "Laplace")
> ## this works but it does not give me an anova output with p values
> anova(m4)
> Analysis of Variance Table
>   Df  Sum Sq Mean Sq
> env1   17858   17858
> har2 879 439
> treat  226131306
> dens   1 1016476 1016476
> env:har2 870 435
> env:treat  21188 594
> har:treat  4 313  78
> env:har:treat  41188 297
> 
>  
> 
> 
> 
> Dr Tom C Cameron
> Genetics, Ecology and Evolution
> IICB, University of Leeds
> Leeds, UK
> Office: +44 (0)113 343 2837
> Lab:+44 (0)113 343 2854
> Fax:+44 (0)113 343 2835
> 
> 
> Email: [EMAIL PROTECTED]
> Webpage: click here
>  
> 
>  


-- 
Andrew Robinson  
Department of Mathematics and StatisticsTel: +61-3-8344-9763
University of Melbourne, VIC 3010 Australia Fax: +61-3-8344-4599
http://www.ms.unimelb.edu.au/~andrewpr
http://blogs.mbs.edu/fishing-in-the-bay/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] analysis of correlation matrices

[ablukacz]

Hello,

I'm looking for a package in R that performs, analysis of correlation 
matrices: cross-classified by 2 factors.

The orginal reference to this method is by CJ Brien Biometrica (1998) 
75(3):469-76.

THank you,

Agnes


-
E. Agnes Richards, Ph.D. Candidate
Department of Zoology
University of Toronto at Mississauga
3359 Mississauga Rd. North
Mississauga, Ont.
Canada
L5L 1C6

lab 905-828-5452
fax 905-828-3792

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Installing Package rgl - Compilation Fails - FreeBSD

[Prof Brian Ripley]

The problem is that rgl is apparently written for GNU make, and has (as 
shipped)

ifdef MAKINGAGL
[EMAIL PROTECTED]@ -Iext
[EMAIL PROTECTED]@
else
PKG_CPPFLAGS= -If:/R/R-2.4.1/src/extra/zlib -DHAVE_PNG_H 
-If:/R/R-2.4.1/src/gnuwin32/bitmap/libpng  -Iext
PKG_LIBS=-lgdi32 -lopengl32 -lglu32 
-Lf:/R/R-2.4.1/src/gnuwin32/bitmap/libpng -
lpng -Lf:/R/R-2.4.1/src/extra/zlib -lz
endif

and similar for BUILDAGL.

That seems to have been written to make it workable on MacOS X. Given that 
configure knows (or could know) the OS, it seems better to write (via 
configure) a separate Makevars for MacOS X and remove all the 
ifdef...endif stuff for everyone else.  (If you do that in Makevars.in it 
should work.)


On Tue, 20 Feb 2007, Rainer Hurling wrote:

> Duncan Murdoch schrieb:
>> On 2/20/2007 1:11 PM, Ranjan Maitra wrote:
>>> Hi Duncan,
>>> I don't know if this will list all the dependencies for your documentation, 
>>> since Rick's error messages did not involve libraries and header files 
>>> already installed by him for something else, perhaps.
>>
>> No, but it's a start:  I'm thinking of something more like a FAQ than
>> comprehensive documentation.  Comprehensive docs are not feasible to
>> maintain (there are so many systems that Daniel and I don't use), but
>> hints that two non-standard packages solved one person's problems might
>> be enough of a hint to get someone else going.
>
> Duncan,
>
> thank you for this purpose. I am such a person who could need some help
> with installing rgl and hope it is ok to place it in this thread.
>
> My trial to install rgl_0.70.tar.gz on R-2.4.1 on FreeBSD 7.0-CURRENT
> ends up with errors. Obiously something is wrong with Makevars. I have
> no idea what to do next. 'rgl' is one of the rare cases that do not
> install under R on FreeBSD.
>
> The install messages are short:
>
> -
> #R CMD INSTALL rgl_0.70.tar.gz
> * Installing to library '/usr/local/lib/R/library'
> * Installing *source* package 'rgl' ...
> checking for gcc... gcc
> checking for C compiler default output file name... a.out
> checking whether the C compiler works... yes
> checking whether we are cross compiling... no
> checking for suffix of executables...
> checking for suffix of object files... o
> checking whether we are using the GNU C compiler... yes
> checking whether gcc accepts -g... yes
> checking for gcc option to accept ANSI C... none needed
> checking how to run the C preprocessor... gcc -E
> checking for X... libraries /usr/X11R6/lib, headers /usr/X11R6/include
> checking for libpng-config... yes
> configure: using libpng-config
> configure: using libpng dynamic linkage
> configure: creating ./config.status
> config.status: creating src/Makevars
> ** libs
> "Makevars", line 9: Need an operator
> "Makevars", line 12: Need an operator
> "Makevars", line 15: Need an operator
> "Makevars", line 21: Need an operator
> "Makevars", line 23: Need an operator
> "Makevars", line 36: Need an operator
> "Makevars", line 38: Need an operator
> make: fatal errors encountered -- cannot continue
> chmod: /usr/local/lib/R/library/rgl/libs/*: No such file or directory
> ERROR: compilation failed for package 'rgl'
> ** Removing '/usr/local/lib/R/library/rgl'
> -
>
> Are there any experiences with installing rgl on FreeBSD? Do you know if
> it is at all possible to get rgl to work on FreeBSD? If I can do
> anything like testing or giving more information let me know.
>
> I appreciate any help. Thank you in advance.
>
> Rainer Hurling
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Installing Package rgl - Compilation Fails - FreeBSD

[Rainer Hurling]

Duncan Murdoch schrieb:
> On 2/20/2007 1:11 PM, Ranjan Maitra wrote:
>> Hi Duncan,
>> I don't know if this will list all the dependencies for your documentation, 
>> since Rick's error messages did not involve libraries and header files 
>> already installed by him for something else, perhaps.
> 
> No, but it's a start:  I'm thinking of something more like a FAQ than 
> comprehensive documentation.  Comprehensive docs are not feasible to 
> maintain (there are so many systems that Daniel and I don't use), but 
> hints that two non-standard packages solved one person's problems might 
> be enough of a hint to get someone else going.

Duncan,

thank you for this purpose. I am such a person who could need some help 
with installing rgl and hope it is ok to place it in this thread.

My trial to install rgl_0.70.tar.gz on R-2.4.1 on FreeBSD 7.0-CURRENT 
ends up with errors. Obiously something is wrong with Makevars. I have 
no idea what to do next. 'rgl' is one of the rare cases that do not 
install under R on FreeBSD.

The install messages are short:

-
#R CMD INSTALL rgl_0.70.tar.gz
* Installing to library '/usr/local/lib/R/library'
* Installing *source* package 'rgl' ...
checking for gcc... gcc
checking for C compiler default output file name... a.out
checking whether the C compiler works... yes
checking whether we are cross compiling... no
checking for suffix of executables...
checking for suffix of object files... o
checking whether we are using the GNU C compiler... yes
checking whether gcc accepts -g... yes
checking for gcc option to accept ANSI C... none needed
checking how to run the C preprocessor... gcc -E
checking for X... libraries /usr/X11R6/lib, headers /usr/X11R6/include
checking for libpng-config... yes
configure: using libpng-config
configure: using libpng dynamic linkage
configure: creating ./config.status
config.status: creating src/Makevars
** libs
"Makevars", line 9: Need an operator
"Makevars", line 12: Need an operator
"Makevars", line 15: Need an operator
"Makevars", line 21: Need an operator
"Makevars", line 23: Need an operator
"Makevars", line 36: Need an operator
"Makevars", line 38: Need an operator
make: fatal errors encountered -- cannot continue
chmod: /usr/local/lib/R/library/rgl/libs/*: No such file or directory
ERROR: compilation failed for package 'rgl'
** Removing '/usr/local/lib/R/library/rgl'
-

Are there any experiences with installing rgl on FreeBSD? Do you know if 
it is at all possible to get rgl to work on FreeBSD? If I can do 
anything like testing or giving more information let me know.

I appreciate any help. Thank you in advance.

Rainer Hurling

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Calculating the Sharpe ratio

[Bernd Dittmann]

Hi Mark,

thanks for your email.
I used your formula for cumul. returns and plugged them into sharpe:

 > mysharpe <- function(x){
+ return(sharpe(cret(x), r=0, scale=1))
+ }

whereby "cret" is my cumul. returns function as defined by:

 > cret
function(x){
cumprod(diff(log(x))+1)-1
}

For the index series "Index" I obtain a sharpe ratio (r=0 and scale=1) of:

 > mysharpe(Index)
[1] 0.8836429

Do you reckon this result and the method above are correct?

Many thanks in advance!

Bernd




Leeds, Mark (IED) schrieb:
> If the doc says to use cumulated and you didn't, then I supsect the call
> to shaprp in
> Tseries is not correct.  Also, to get PercetnREturns, I hope you did
> diff(log(series))
> Where series is an object containing prices. It's not so clear
> Form your email.
>
> If you want to send in cumulative returns ( which
> You should do if the doc says to ) you just take the returns ( by doing
> above )
> and then , add 1 to each element, do  a cumprod and then subtract 1 so
> something like :
>
> rtns<-diff(log(priceseries)
> oneplusrtns<-1+rtns
> cumprodrtns<-cumprod(oneplusreturns)
>
> cumrtns<-cumprodrtns-1.
>
> Then, the elements in cumrtns represent the cumulative reeturn upto that
> point.
>
> But, test it out with an easy example to make sure because I didn't.
>
>
>
>
> -Original Message-
> From: [EMAIL PROTECTED]
> [mailto:[EMAIL PROTECTED] On Behalf Of Bernd Dittmann
> Sent: Monday, February 19, 2007 8:39 AM
> To: r-help@stat.math.ethz.ch
> Subject: [R] Calculating the Sharpe ratio
>
> Hi useRs,
>
> I am trying to calculate the Sharpe ratio with "sharpe" of the library
> "tseries".
>
> The documentation requires the univariate time series to be a
> portfolio's cumulated returns. In this case, the example given
>
> data(EuStockMarkets)
> dax <- log(EuStockMarkets[,"FTSE"])
>
> is however not the cumulated returns but rather the daily returns of the
> FTSE stock index.
>
> Is this way of calculating the Sharpe ratio correct?
>
> Here are my own data:
>
> yearIndexPercentReturns
> 19851170.091
> 1986129.90.11
> 1987149.90.154
> 1988184.80.233
> 1989223.10.208
> 1990223.20
> 1991220.5-0.012
> 1992208.1-0.056
> 1993202.1-0.029
> 1994203.10.005
> 1995199.6-0.017
> 1996208.60.045
> 1997221.70.063
> 1998233.70.054
> 1999250.50.072
> 2000275.10.098
> 2001298.60.085
> 2002350.60.174
> 2003429.10.224
> 2004507.60.183
> 2005536.60.057
> 2006581.30.083
>
>
> I calculated the Sharpe ratio in two different ways:
> (1) using natural logs as approximation of % returns, using "sharpe" of
> "tseries".
> (2) using the % returns using a variation the "sharpe" function.
>
> In both cases I used the risk free rate r=0 and scale=1 since I am using
> annual data already.
>
> My results:
>
> METHOD 1: "sharpe":
>
>  > index <- log(Index)
>  > sharpe(index, scale=1)
> [1] 0.9614212
>
>
>
> METHOD 2: my own %-based formula:
>
>  > mysharp
> function(x, r=0, scale=sqrt(250))
> {
> if (NCOL(x) > 1)
> stop("x is not a vector or univariate time series") if (any(is.na(x)))
> stop("NAs in x") if (NROW(x) ==1)
> return(NA)
> else{
> return(scale * (mean(x) - r)/sd(x))
> }
> }
>
>
>
>  > mysharp(PercentReturns, scale=1)
> [1] 0.982531
>
>
> Both Sharp ratios differ only slightly since logs approximate percentage
> changes (returns).
>
>
> Are both methods correct, esp. since I am NOT using cumulated returns as
>
> the manual says?
>
> If cumulated returns were supposed to be used, could I cumulate the 
> %-returns with "cumsum(PercentReturns)"?
>
> Many thanks in advance!
>
> Bernd
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
>
> This is not an offer (or solicitation of an offer) to buy/sell the 
> securities/instruments mentioned or an official confirmation.  Morgan Stanley 
> may deal as principal in or own or act as market maker for 
> securities/instruments mentioned or may advise the issuers.  This is not 
> research and is not from MS Research but it may refer to a research 
> analyst/research report.  Unless indicated, these views are the author's and 
> may differ from those of Morgan Stanley research or others in the Firm.  We 
> do not represent this is accurate or complete and we may not update this.  
> Past performance is not indicative of future returns.  For additional 
> information, research reports and important disclosures, contact me or see 
> https://secure.ms.com/servlet/cls.  You should not use e-mail to request, 
> authorize or effect the purchase or sale of any security or instrument, to 
> send transfer instructions

Re: [R] linux gplots install unhappy

[Marc Schwartz]

On Tue, 2007-02-20 at 18:59 +, Prof Brian Ripley wrote:
> On Tue, 20 Feb 2007, Marc Schwartz wrote:
> 
> > On Tue, 2007-02-20 at 10:17 -0800, Randy Zelick wrote:
> >> Hello all,
> >>
> >> I use R on both windows and a "mainframe" linux installation (RedHat
> >> enterprise 3.0, which they tell me is soon to be upgraded to 4.0). On
> >> windows I installed the package gplots without trouble, and it works fine.
> >> When I attempted to do the same on the unix computer, the following error
> >> message was forthcoming:
> >>
> >>
> >>
> >>
> >> downloaded 216Kb
> >>
> >> * Installing *source* package 'gplots' ...
> >> ** R
> >> ** data
> >> ** inst
> >> ** preparing package for lazy loading
> >> Loading required package: gtools
> >> Warning in library(pkg, character.only = TRUE, logical = TRUE, lib.loc =
> >> lib.loc) :
> >>   there is no package called 'gtools'
> >> Error: package 'gtools' could not be loaded
> >> Execution halted
> >> ERROR: lazy loading failed for package 'gplots'
> >> ** Removing '/n/fs/disk/resuser02/u/zelickr/R/library/gplots'
> >>
> >> The downloaded packages are in
> >>  /tmp/RtmpikM2JW/downloaded_packages
> >> Warning messages:
> >> 1: installation of package 'gplots' had non-zero exit status in:
> >> install.packages("gplots", lib = "~/R/library")
> >> 2: cannot create HTML package index in:
> >> tools:::unix.packages.html(.Library)
> >>
> >>
> >>
> >> Can someone provide the bit of information I need to progress with this?
> >>
> >> Thanks very much,
> >>
> >> =Randy=
> >
> > gplots has a dependency on other packages (gtools and gdata).
> >
> > Thus, when you install it use:
> >
> >  install.packages("gplots", dependencies = TRUE, ...)
> >
> > That will download and install the other packages as well.
> >
> > There should have been a similar requirement under Windows, so not sure
> > what may have been different there, unless there is some confounding due
> 
> You don't need gtools to install a binary package, just to use it!
> I suspect the menu was used on Windows, where dependencies = TRUE is the 
> default.
> 
> Come 2.5.0 this sort of query will go away, as the essential dependencies 
> become the default on all platforms.

Thanks for the clarifications.

Regards,

Marc

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] contextstack overflow

[Prof Brian Ripley]

On Tue, 20 Feb 2007, Steven Finch wrote:


Hello!

I written several implementations in R of Rémy's algorithm for
generating random ordered strictly binary trees on 2n+1 vertices.

One implementation involves manufacturing character strings like:

"X <- list(0,list(0,0))"

for the case n=2.  If I perform the following two steps:

cmd <- "X <- list(0,list(0,0))"
eval(parse(text=cmd))

then X becomes a true nested list in R.  This works fine for n=2,
but often for n=200, an error message:

Error in parse(text = cmd) : contextstack overflow

appears and execution stops.  Clearly there exists an upper bound
on the allowable depth of nestings in R!  Can this upper bound be
easily increased?


It is hardcoded (as 50) in 6 places in gram.c, so fairly easily changed.


Other implementations avoid this problem, so this issue is not
crucial to me.  I do wish, however, to understand the limits of
this particular approach.  Thank you!

Steve Finch
http://algo.inria.fr/bsolve/

P.S.  If anyone else has written R code for generating random
trees (on a fixed number of vertices), I would enjoy seeing this!

_
Don’t miss your chance to WIN 10 hours of private jet travel from Microsoft® 
Office Live http://clk.atdmt.com/MRT/go/mcrssaub0540002499mrt/direct/01/





--
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] printing intermediate lines while in a function

[Bart Joosen]

Thanks for al the tips, it was the readline() function who did the trick for 
me, so thanks for all of your input

Kind Regards

Bart


 Original message:
>
>> 
>>
>> But unfortunately R will do first the calculations and then
>> afterwards return the strings.
>>
>> Is there a way around?
>

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] baseline fitters

[Tuszynski, Jaroslaw W.]

I am not surprised at slowness of runquantile, since it is trying to
perform n=4500 partial sorts of k=225 elements. Here are some thoughts
at speeding it up:
1) playing with different endrule settings can save some time, but
usually results with undesirable effects at first and last 112 values.
All rum* functions in caTools use low level C code for inner elements
between k/2 and n-k/2. However the elements at the edge are calculated
using R functions. In case of runquantile with endrule="func" that means
k calls of R quantile function. One option for endrule not available at
present would be to pad both sides of the array with k/2 numbers and
than use endrule="trim". The trick would be to pick good value for the
padding number.

2) you mentioned that you "jimmied something together
with runmin and runmedian". I would try something like runmean with
window of size 5, 15, 25 and than runmin with window size k. The first
one should get rid of your 'reverse-spikes' and second would take
running min of your smoothed function.

Best,
Jarek Tuszynski

-Original Message-
From: Thaden, John J [mailto:[EMAIL PROTECTED] 
Sent: Tuesday, February 20, 2007 1:23 PM
To: r-help@stat.math.ethz.ch
Cc: Tuszynski, Jaroslaw W.
Subject: baseline fitters

I am pretty pleased with baselines I fit to chromatograms using the
runquantile() function in caTools(v1.6) when its probs parameter is 
set to 0.2 and its k parameter to ~1/20th of n (e.g., k ~ 225 for n ~ 
4500, where n is time series length).  This ignores occasional low-
side outliers, and, after baseline subtraction, I can re-adjust any
negative values to zero.
  
But runquantile's computation time proves exceedingly long for my large
datasets, particularly if I set the endrule parameter to 'func'.  Here
is
what caTools author Jarek Tuszynski says about relative speeds of
various
running-window functions:

   - runmin, runmax, runmean run at O(n) 
   - runmean(..., alg="exact") can have worst case speed of O(n^2) for 
 some small data vectors, but average case is still close to O(n). 
   - runquantile and runmad run at O(n*k) 
   - runmed - related R function runs at O(n*log(k))

The obvious alternative runmin() performs poorly due to dropout (zero-
and low-value 'reverse-spikes') in the data. And runmed fits a baseline
that,
upon subtraction, obviously will send half the values into the negative,
not
suitable for my application. I jimmied something together
with runmin and runmedian that is considerably faster; unfortunately,
the fit seems less good, at least by eye, due still to the bad runmin
behavior.

I'd be interested in other baseline fitting suggestions implemented
or implementable in R (I'm using v. 2.4.1pat under WinXP).  Why, for
instance, can I not find a running trimmean function? Because it 
offers no computational savings over runquantile?

Also, the 'func' setting for the caTools endrule parameter--which
adjusts the
value of k as ends are approached--is said not to be optimized (using
this option does add further time to computations).  Is there an alter-
native that would be speedier, e.g., setting endrule = "keep" and then
subsequently treating ends with R's smoothEnds() function?

-John Thaden
Little Rock, Arkansas USA

Confidentiality Notice: This e-mail message, including any a...{{dropped}}

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Installing Package rgl - Compilation Fails - Summary

[Prof Brian Ripley]

On Tue, 20 Feb 2007, Ranjan Maitra wrote:

> Hi Duncan,
>
> I don't know if this will list all the dependencies for your 
> documentation, since Rick's error messages did not involve libraries and 
> header files already installed by him for something else, perhaps.

It will not, and most of this is already in the README.  The problem is 
the penchant for linux distros to break standard pieces of software into 
ever smaller pieces, so naming the pieces is often no great help 6 months 
later.

On  FC5/6 I think

yum install mesa-libGL-devel mesa-libGLU-devel libXext-devel libpng-devel

should do it.  (Note that mesa-libGL-devel appears to have been already 
there, and libXext-devel pulls in the rest of the X11 stuff.)

A quick check shows -lXext is actually not needed, so it could be removed 
from configure.ac and the requirements: that would make a lot more sense.

>
> Just a thought.
>
> Ranjan
>
>
>
> On Tue, 20 Feb 2007 11:59:24 -0500 Rick Bilonick <[EMAIL PROTECTED]> wrote:
>
>> Summarizing:
>>
>> I'm running R 2.4.1 on a current FC6 32-bit system. In order to have the
>> rgl R package install, I needed to install both mesa-libGLU-devel (FC6
>> version is 6.5.1-9) and libXext-devel (FC6) rpm packages. Thanks to
>> everyone who commented.
>>
>> Rick B.


-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] linux gplots install unhappy

[Prof Brian Ripley]

On Tue, 20 Feb 2007, Marc Schwartz wrote:

> On Tue, 2007-02-20 at 10:17 -0800, Randy Zelick wrote:
>> Hello all,
>>
>> I use R on both windows and a "mainframe" linux installation (RedHat
>> enterprise 3.0, which they tell me is soon to be upgraded to 4.0). On
>> windows I installed the package gplots without trouble, and it works fine.
>> When I attempted to do the same on the unix computer, the following error
>> message was forthcoming:
>>
>>
>>
>>
>> downloaded 216Kb
>>
>> * Installing *source* package 'gplots' ...
>> ** R
>> ** data
>> ** inst
>> ** preparing package for lazy loading
>> Loading required package: gtools
>> Warning in library(pkg, character.only = TRUE, logical = TRUE, lib.loc =
>> lib.loc) :
>>   there is no package called 'gtools'
>> Error: package 'gtools' could not be loaded
>> Execution halted
>> ERROR: lazy loading failed for package 'gplots'
>> ** Removing '/n/fs/disk/resuser02/u/zelickr/R/library/gplots'
>>
>> The downloaded packages are in
>>  /tmp/RtmpikM2JW/downloaded_packages
>> Warning messages:
>> 1: installation of package 'gplots' had non-zero exit status in:
>> install.packages("gplots", lib = "~/R/library")
>> 2: cannot create HTML package index in:
>> tools:::unix.packages.html(.Library)
>>
>>
>>
>> Can someone provide the bit of information I need to progress with this?
>>
>> Thanks very much,
>>
>> =Randy=
>
> gplots has a dependency on other packages (gtools and gdata).
>
> Thus, when you install it use:
>
>  install.packages("gplots", dependencies = TRUE, ...)
>
> That will download and install the other packages as well.
>
> There should have been a similar requirement under Windows, so not sure
> what may have been different there, unless there is some confounding due

You don't need gtools to install a binary package, just to use it!
I suspect the menu was used on Windows, where dependencies = TRUE is the 
default.

Come 2.5.0 this sort of query will go away, as the essential dependencies 
become the default on all platforms.

> to your not installing the packages in standard locations, given some of
> the output above. I presume that this is because you don't have root
> access on the Linux server and you are installing to your local user
> path.
>
> HTH,
>
> Marc Schwartz
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] User defined split function in rpart

[Tobias Guennel]

I have made some progress with the user defined splitting function and I got
a lot of the things I needed to work. However, I am still stuck on accessing
the node data. It would probably be enough if somebody could tell me, how I
can access the original data frame of the call to rpart. 
So if the call is: fit0 <- rpart(Sat ~Infl +Cont+ Type,
 housing, control=rpart.control(minsplit=10, xval=0),
 method=alist)
how can I access the housing data frame within the user defined splitting
function?

Any input would be highly appreciated!

Thank you
Tobias Guennel

-Original Message-
From: Tobias Guennel [mailto:[EMAIL PROTECTED] 
Sent: Monday, February 19, 2007 3:40 PM
To: 'r-help@stat.math.ethz.ch'
Subject: [R] User defined split function in rpart

Maybe I should explain my Problem a little bit more detailed.
The rpart package allows for user defined split functions. An example is
given in the source/test directory of the package as usersplits.R.
The comments say that three functions have to be supplied:
1. "The 'evaluation' function.  Called once per node.
  Produce a label (1 or more elements long) for labeling each node,
  and a deviance." 
2. The split function, where most of the work occurs.
   Called once per split variable per node.
3. The init function:
   fix up y to deal with offsets
   return a dummy parms list
   numresp is the number of values produced by the eval routine's "label".

I have altered the evaluation function and the split function for my needs.
Within those functions, I need to fit a proportional odds model to the data
of the current node. I am using the polr() routine from the MASS package to
fit the model. 
Now my problem is, how can I call the polr() function only with the data of
the current node. That's what I tried so far:

evalfunc <- function(y,x,parms,data) {
   
pomnode<-polr(data$y~data$x,data,weights=data$Freq)
parprobs<-predict(pomnode,type="probs")
dev<-0
K<-dim(parprobs)[2]
N<-dim(parprobs)[1]/K
for(i in 1:N){
tempsum<-0
Ni<-0
for(l in 1:K){
Ni<-Ni+data$Freq[K*(i-1)+l]
}
for(j in 1:K){
tempsum<-tempsum+data$Freq[K*(i-1)+j]/Ni*log(parprobs[i,j]*Ni/data$Freq[K*(i
-1)+j])
}
dev=dev+Ni*tempsum
}
dev=-2*dev
wmean<-1
list(label= wmean, deviance=dev)

} 

I get the error: Error in eval(expr, envir, enclos) : argument "data" is
missing, with no default

How can I use the data of the current node?

Thank you
Tobias Guennel

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] linux gplots install unhappy

[Marc Schwartz]

On Tue, 2007-02-20 at 10:17 -0800, Randy Zelick wrote:
> Hello all,
> 
> I use R on both windows and a "mainframe" linux installation (RedHat 
> enterprise 3.0, which they tell me is soon to be upgraded to 4.0). On 
> windows I installed the package gplots without trouble, and it works fine. 
> When I attempted to do the same on the unix computer, the following error 
> message was forthcoming:
> 
> 
> 
> 
> downloaded 216Kb
> 
> * Installing *source* package 'gplots' ...
> ** R
> ** data
> ** inst
> ** preparing package for lazy loading
> Loading required package: gtools
> Warning in library(pkg, character.only = TRUE, logical = TRUE, lib.loc = 
> lib.loc) :
>   there is no package called 'gtools'
> Error: package 'gtools' could not be loaded
> Execution halted
> ERROR: lazy loading failed for package 'gplots'
> ** Removing '/n/fs/disk/resuser02/u/zelickr/R/library/gplots'
> 
> The downloaded packages are in
>  /tmp/RtmpikM2JW/downloaded_packages
> Warning messages:
> 1: installation of package 'gplots' had non-zero exit status in: 
> install.packages("gplots", lib = "~/R/library")
> 2: cannot create HTML package index in: 
> tools:::unix.packages.html(.Library)
> 
> 
> 
> Can someone provide the bit of information I need to progress with this?
> 
> Thanks very much,
> 
> =Randy=

gplots has a dependency on other packages (gtools and gdata).

Thus, when you install it use:

  install.packages("gplots", dependencies = TRUE, ...)

That will download and install the other packages as well.

There should have been a similar requirement under Windows, so not sure
what may have been different there, unless there is some confounding due
to your not installing the packages in standard locations, given some of
the output above. I presume that this is because you don't have root
access on the Linux server and you are installing to your local user
path.

HTH,

Marc Schwartz

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] baseline fitters

[Thaden, John J]

I am pretty pleased with baselines I fit to chromatograms using the
runquantile() function in caTools(v1.6) when its probs parameter is 
set to 0.2 and its k parameter to ~1/20th of n (e.g., k ~ 225 for n ~ 
4500, where n is time series length).  This ignores occasional low-
side outliers, and, after baseline subtraction, I can re-adjust any
negative values to zero.
  
But runquantile's computation time proves exceedingly long for my large
datasets, particularly if I set the endrule parameter to 'func'.  Here is
what caTools author Jarek Tuszynski says about relative speeds of various
running-window functions:

   - runmin, runmax, runmean run at O(n) 
   - runmean(..., alg="exact") can have worst case speed of O(n^2) for 
 some small data vectors, but average case is still close to O(n). 
   - runquantile and runmad run at O(n*k) 
   - runmed - related R function runs at O(n*log(k))

The obvious alternative runmin() performs poorly due to dropout (zero-
and low-value 'reverse-spikes') in the data. And runmed fits a baseline that,
upon subtraction, obviously will send half the values into the negative, not
suitable for my application. I jimmied something together
with runmin and runmedian that is considerably faster; unfortunately,
the fit seems less good, at least by eye, due still to the bad runmin
behavior.

I'd be interested in other baseline fitting suggestions implemented
or implementable in R (I'm using v. 2.4.1pat under WinXP).  Why, for
instance, can I not find a running trimmean function? Because it 
offers no computational savings over runquantile?

Also, the 'func' setting for the caTools endrule parameter--which adjusts the
value of k as ends are approached--is said not to be optimized (using
this option does add further time to computations).  Is there an alter-
native that would be speedier, e.g., setting endrule = "keep" and then
subsequently treating ends with R's smoothEnds() function?

-John Thaden
Little Rock, Arkansas USA

Confidentiality Notice: This e-mail message, including any a...{{dropped}}

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] linux gplots install unhappy

[Benilton Carvalho]

well, it's complaining because you don't have gtools installed.

how about:

install.packages("gplots", dep=T)

?

b

On Feb 20, 2007, at 1:17 PM, Randy Zelick wrote:

> Hello all,
>
> I use R on both windows and a "mainframe" linux installation (RedHat
> enterprise 3.0, which they tell me is soon to be upgraded to 4.0). On
> windows I installed the package gplots without trouble, and it  
> works fine.
> When I attempted to do the same on the unix computer, the following  
> error
> message was forthcoming:
>
>
>
>
> downloaded 216Kb
>
> * Installing *source* package 'gplots' ...
> ** R
> ** data
> ** inst
> ** preparing package for lazy loading
> Loading required package: gtools
> Warning in library(pkg, character.only = TRUE, logical = TRUE,  
> lib.loc =
> lib.loc) :
>   there is no package called 'gtools'
> Error: package 'gtools' could not be loaded
> Execution halted
> ERROR: lazy loading failed for package 'gplots'
> ** Removing '/n/fs/disk/resuser02/u/zelickr/R/library/gplots'
>
> The downloaded packages are in
>  /tmp/RtmpikM2JW/downloaded_packages
> Warning messages:
> 1: installation of package 'gplots' had non-zero exit status in:
> install.packages("gplots", lib = "~/R/library")
> 2: cannot create HTML package index in:
> tools:::unix.packages.html(.Library)
>
>
>
> Can someone provide the bit of information I need to progress with  
> this?
>
> Thanks very much,
>
> =Randy=
>
> R. Zelick email: [EMAIL PROTECTED]
> Department of Biology voice: 503-725-3086
> Portland State University fax:   503-725-3888
>
> mailing:
> P.O. Box 751
> Portland, OR 97207
>
> shipping:
> 1719 SW 10th Ave, Room 246
> Portland, OR 97201
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting- 
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] linux gplots install unhappy

[Aimin Yan]

install gtools package firstly

Aimin

At 12:17 PM 2/20/2007, Randy Zelick wrote:
>gtools'

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Installing Package rgl - Compilation Fails - Summary

[Duncan Murdoch]

On 2/20/2007 1:11 PM, Ranjan Maitra wrote:
> Hi Duncan,
> 
> I don't know if this will list all the dependencies for your documentation, 
> since Rick's error messages did not involve libraries and header files 
> already installed by him for something else, perhaps.

No, but it's a start:  I'm thinking of something more like a FAQ than 
comprehensive documentation.  Comprehensive docs are not feasible to 
maintain (there are so many systems that Daniel and I don't use), but 
hints that two non-standard packages solved one person's problems might 
be enough of a hint to get someone else going.

Duncan Murdoch

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] linux gplots install unhappy

[Randy Zelick]

Hello all,

I use R on both windows and a "mainframe" linux installation (RedHat 
enterprise 3.0, which they tell me is soon to be upgraded to 4.0). On 
windows I installed the package gplots without trouble, and it works fine. 
When I attempted to do the same on the unix computer, the following error 
message was forthcoming:




downloaded 216Kb

* Installing *source* package 'gplots' ...
** R
** data
** inst
** preparing package for lazy loading
Loading required package: gtools
Warning in library(pkg, character.only = TRUE, logical = TRUE, lib.loc = 
lib.loc) :
  there is no package called 'gtools'
Error: package 'gtools' could not be loaded
Execution halted
ERROR: lazy loading failed for package 'gplots'
** Removing '/n/fs/disk/resuser02/u/zelickr/R/library/gplots'

The downloaded packages are in
 /tmp/RtmpikM2JW/downloaded_packages
Warning messages:
1: installation of package 'gplots' had non-zero exit status in: 
install.packages("gplots", lib = "~/R/library")
2: cannot create HTML package index in: 
tools:::unix.packages.html(.Library)



Can someone provide the bit of information I need to progress with this?

Thanks very much,

=Randy=

R. Zelick   email: [EMAIL PROTECTED]
Department of Biology   voice: 503-725-3086
Portland State University   fax:   503-725-3888

mailing:
P.O. Box 751
Portland, OR 97207

shipping:
1719 SW 10th Ave, Room 246
Portland, OR 97201

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Installing Package rgl - Compilation Fails - Summary

[Ranjan Maitra]

Hi Duncan,

I don't know if this will list all the dependencies for your documentation, 
since Rick's error messages did not involve libraries and header files already 
installed by him for something else, perhaps.

Just a thought.

Ranjan



On Tue, 20 Feb 2007 11:59:24 -0500 Rick Bilonick <[EMAIL PROTECTED]> wrote:

> Summarizing:
> 
> I'm running R 2.4.1 on a current FC6 32-bit system. In order to have the
> rgl R package install, I needed to install both mesa-libGLU-devel (FC6
> version is 6.5.1-9) and libXext-devel (FC6) rpm packages. Thanks to
> everyone who commented.
> 
> Rick B.
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] interaction term and scatterplot3d

[Patrick Giraudoux]

Sorry to answer to myself. A solution was trivial with lattice... (as 
often !)

library(lattice)

modbusetmp<-lm(IKA_buse ~ Ct *Cc,data=dtbuse)

G1<-cloud(IKA_buse~Ct*Cc,type="h",data=dtbuse)
G2<-cloud(IKA_buse~Ct*Cc,data=dtbuse)

seqCc<-seq(min(dtbuse$Cc),max(dtbuse$Cc),l=20)
seqCt<-seq(min(dtbuse$Ct),max(dtbuse$Ct),l=20)
grille<-expand.grid(Cc=seqCc,Ct=seqCt)
zbuse<-predict(modbusetmp,newdata=grille,type=response)
G3<-wireframe(zbuse~grille$Ct*grille$Cc)

print(G3,more=T)
print(G1,more=T)
print(G2)

Some obvious improvements can be done with labels and differencial 
colors (above or below the surface) can be attributed comparing 
observations to predicted values at the same x y values...


Patrick Giraudoux a écrit :
> Dear Listers,
>
> I would be interested in representing a trend surface including an 
> interaction term z = f(x,y,x.y) - eg the type of chart obtained with 
> persp() or wireframe(), then adding datapoints as a cloud, ideally 
> with dots which are under the surface in a color, and those who are 
> above in another color. An other option would be to draw segments 
> between the dots and the ground of the chart.
>
> scatterplot3d looks like being close to do such things except it does 
> not to include (to my knowledge) a coefficient for the interaction 
> term (thus just model z = f(x,y).
>
> Does anybody has an idea where to start with this and the main steps? 
> Or a place/website where some script examples can be found?
>
> Patrick
>
>

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Calculating the Sharpe ratio

[AA]

Hi Bernd,

It seems to me that the sharpe function uses the diff to calculate the first
differnce of the input which
is a cumalative return series. This is consistant. see sharpe
Basically you need the rate of returns ((priceFinal - Price
Initial)/priceInitial) which you used. In your case r = 0, so you are
calculating the simple signal to noise ratio (mean/std). It looks correct to
me.
The example given in the sharpe function does not seem to be correct though.

data(EuStockMarkets)
dax <- log(EuStockMarkets[,"DAX"])

The input is a time series of closing prices. the simple log does not give a
rate of returns.
you need to use the log Difference of prices to get the returns. and then
divide by std of the return series.

I Hope this helps.

A.
On 2/19/07, Bernd Dittmann <[EMAIL PROTECTED]> wrote:
>
> Hi useRs,
>
> I am trying to calculate the Sharpe ratio with "sharpe" of the library
> "tseries".
>
> The documentation requires the univariate time series to be a
> portfolio's cumulated returns. In this case, the example given
>
> data(EuStockMarkets)
> dax <- log(EuStockMarkets[,"FTSE"])
>
> is however not the cumulated returns but rather the daily returns of the
> FTSE stock index.
>
> Is this way of calculating the Sharpe ratio correct?
>
> Here are my own data:
>
> yearIndexPercentReturns
> 19851170.091
> 1986129.90.11
> 1987149.90.154
> 1988184.80.233
> 1989223.10.208
> 1990223.20
> 1991220.5-0.012
> 1992208.1-0.056
> 1993202.1-0.029
> 1994203.10.005
> 1995199.6-0.017
> 1996208.60.045
> 1997221.70.063
> 1998233.70.054
> 1999250.50.072
> 2000275.10.098
> 2001298.60.085
> 2002350.60.174
> 2003429.10.224
> 2004507.60.183
> 2005536.60.057
> 2006581.30.083
>
>
> I calculated the Sharpe ratio in two different ways:
> (1) using natural logs as approximation of % returns, using "sharpe" of
> "tseries".
> (2) using the % returns using a variation the "sharpe" function.
>
> In both cases I used the risk free rate r=0 and scale=1 since I am using
> annual data already.
>
> My results:
>
> METHOD 1: "sharpe":
>
> > index <- log(Index)
> > sharpe(index, scale=1)
> [1] 0.9614212
>
>
>
> METHOD 2: my own %-based formula:
>
> > mysharp
> function(x, r=0, scale=sqrt(250))
> {
> if (NCOL(x) > 1)
> stop("x is not a vector or univariate time series")
> if (any(is.na(x)))
> stop("NAs in x")
> if (NROW(x) ==1)
> return(NA)
> else{
> return(scale * (mean(x) - r)/sd(x))
> }
> }
>
>
>
> > mysharp(PercentReturns, scale=1)
> [1] 0.982531
>
>
> Both Sharp ratios differ only slightly since logs approximate percentage
> changes (returns).
>
>
> Are both methods correct, esp. since I am NOT using cumulated returns as
> the manual says?
>
> If cumulated returns were supposed to be used, could I cumulate the
> %-returns with "cumsum(PercentReturns)"?
>
> Many thanks in advance!
>
> Bernd
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] memory management uestion [Broadcast]

[Federico Calboli]

Charles C. Berry wrote:

> 
> 
> This is a bit different than your original post (where it appeared that 
> you were manipulating one row of a matrix at a time), but the issue is 
> the same.
> 
> As suggested in my earlier email this looks like a caching issue, and 
> this is not peculiar to R.
> 
> Viz.
> 
> "Most modern CPUs are so fast that for most program workloads the 
> locality of reference of memory accesses, and the efficiency of the 
> caching and memory transfer between different levels of the hierarchy, 
> is the practical limitation on processing speed. As a result, the CPU 
> spends much of its time idling, waiting for memory I/O to complete."
> 
> (from http://en.wikipedia.org/wiki/Memory_hierarchy)
> 
> 
> The computation you have is challenging to your cache, and the effect of 
> dropping unused columns of your 'data' object by assiging the columns 
> used  to 'srow' and 'drow' has lightened the load.
> 
> If you do not know why SAXPY and friends are written as they are, a 
> little bit of study will be rewarded by a much better understanding of 
> these issues. I think Golub and Van Loan's 'Matrix Computations' touches 
> on this (but I do not have my copy close to hand to check).

Thanks for the clarifications. My bottom line is, I prefer dummy variables 
because they allow me to write cleaner code, with a shorter line for the same 
instruction i.e. less chances of creeping errors (+ turning into -, etc).

I've been challenged that that's memory inefficent, and I wanted to have the 
opinion of people with more experience than mine on the matter.

Best,

Fede

-- 
Federico C. F. Calboli
Department of Epidemiology and Public Health
Imperial College, St Mary's Campus
Norfolk Place, London W2 1PG

Tel  +44 (0)20 7594 1602 Fax (+44) 020 7594 3193

f.calboli [.a.t] imperial.ac.uk
f.calboli [.a.t] gmail.com

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] memory management uestion [Broadcast]

[Charles C. Berry]

On Tue, 20 Feb 2007, Federico Calboli wrote:

> Liaw, Andy wrote:
>>  I don't see why making copies of the columns you need inside the loop is
>>  "better" memory management.  If the data are in a matrix, accessing
>>  elements is quite fast.  If you're worrying about speed of that, do what
>>  Charles suggest: work with the transpose so that you are accessing
>>  elements in the same column in each iteration of the loop.
>
> As I said, this is pretty academic, I am not looking for how to do something 
> differetly.
>
> Having said that, let me present this code:
>
> for(i in gp){
>new[i,1] = ifelse(srow[i]>0, new[srow[i],zippo[i]], sav[i])
>new[i,2] = ifelse(drow[i]>0, new[drow[i],zappo[i]], sav[i])
>  }
>
> where gp is large vector and srow and drow are the dummy variables for:
>
> srow = data[,2]
> drow = data[,4]
>
> If instead of the dummy variable I access the array directly (and its' a 
> 60 x 6 array) the loop takes 2/3 days --not sure here, I killed it after 
> 48 hours.
>
> If I use dummy variables the code runs in 10 minutes-ish.
>
> Comments?


This is a bit different than your original post (where it appeared that 
you were manipulating one row of a matrix at a time), but the issue is the 
same.

As suggested in my earlier email this looks like a caching issue, and this 
is not peculiar to R.

Viz.

"Most modern CPUs are so fast that for most program workloads the locality 
of reference of memory accesses, and the efficiency of the caching and 
memory transfer between different levels of the hierarchy, is the 
practical limitation on processing speed. As a result, the CPU spends much 
of its time idling, waiting for memory I/O to complete."

(from http://en.wikipedia.org/wiki/Memory_hierarchy)


The computation you have is challenging to your cache, and the effect of 
dropping unused columns of your 'data' object by assiging the 
columns used  to 'srow' and 'drow' has lightened the load.

If you do not know why SAXPY and friends are written as they are, a little 
bit of study will be rewarded by a much better understanding of these 
issues. I think Golub and Van Loan's 'Matrix Computations' touches on this 
(but I do not have my copy close to hand to check).


>
> Best,
>
> Fede
>
> -- 
> Federico C. F. Calboli
> Department of Epidemiology and Public Health
> Imperial College, St Mary's Campus
> Norfolk Place, London W2 1PG
>
> Tel  +44 (0)20 7594 1602 Fax (+44) 020 7594 3193
>
> f.calboli [.a.t] imperial.ac.uk
> f.calboli [.a.t] gmail.com
>

Charles C. Berry(858) 534-2098
  Dept of Family/Preventive Medicine
E mailto:[EMAIL PROTECTED]   UC San Diego
http://biostat.ucsd.edu/~cberry/ La Jolla, San Diego 92093-0901

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] testing slopes

[Dimitris Rizopoulos]

two options are to use an offset term or the linear.hypothesis()  
function from package car, e.g.,

y <- rnorm(100, 2 - 1 * (x <- runif(100, -3, 3)), 3)

fit0 <- lm(y ~ 1 + offset(-x))
fit1 <- lm(y ~ x)
anova(fit0, fit1)

library(car)
linear.hypothesis(fit1, c("x = -1"))


I hope it helps.

Best,
Dimitris


Dimitris Rizopoulos
Ph.D. Student
Biostatistical Centre
School of Public Health
Catholic University of Leuven

Address: Kapucijnenvoer 35, Leuven, Belgium
Tel: +32/(0)16/336899
Fax: +32/(0)16/337015
Web: http://med.kuleuven.be/biostat/
  http://www.student.kuleuven.be/~m0390867/dimitris.htm


Quoting Indermaur Lukas <[EMAIL PROTECTED]>:

> Hello
> Instead of testing against 0 i would like to test regression slopes   
> against -1. Any idea if there's an R script (package?) available.
>
> Thanks for any hint.
> Cheers
> Lukas
>
>
>
>
> °°°
> Lukas Indermaur, PhD student
> eawag / Swiss Federal Institute of Aquatic Science and Technology
> ECO - Department of Aquatic Ecology
> Überlandstrasse 133
> CH-8600 Dübendorf
> Switzerland
>
> Phone: +41 (0) 71 220 38 25
> Fax: +41 (0) 44 823 53 15
> Email: [EMAIL PROTECTED]
> www.lukasindermaur.ch
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>



Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Installing Package rgl - Compilation Fails - Summary

[Rick Bilonick]

Summarizing:

I'm running R 2.4.1 on a current FC6 32-bit system. In order to have the
rgl R package install, I needed to install both mesa-libGLU-devel (FC6
version is 6.5.1-9) and libXext-devel (FC6) rpm packages. Thanks to
everyone who commented.

Rick B.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Help with xlab and ylab distance from axes

[Marc Schwartz]

On Tue, 2007-02-20 at 11:31 -0500, Michael Kubovy wrote:
> Dear r-helpers,
> 
> In basic graphics, I have a figure with x and y axes suppresed. I  
> would like to move the xlab and the ylab closer to the axes. Do I  
> have to use mtext()?

Michael,

You could set the first value in par("mgp") in the plot() call.  

Compare:

# Default values
plot(1, mgp = c(3, 1, 0), axes = FALSE)

# Move the labels closer
plot(1, mgp = c(1, 1, 0), axes = FALSE)


Changing that value is comparable to changing the 'line' argument in
mtext().

See ?par

HTH,

Marc Schwartz

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Standardized residual variances in SEM

[Mathieu d'Acremont]

Hello,

I'm using the "sem" package to do a confirmatory factor analysis on data 
collected with a questionnaire. In the model, there is a unique factor G 
and 23 items. I would like to calculate the standardized residual 
variance of the observed variables. "Sem" only gives the residual 
variance with the "summary" function, or the standardized loadings with 
the "standardized.coefficients" function (see below). Does anybody know 
how to standardized the residual variance ?

Sincerely yours,

 > summary(semModif.tmp, digits=2)

  Model Chisquare =  601   Df =  229 Pr(>Chisq) = 0
  Chisquare (null model) =  2936   Df =  253
  Goodness-of-fit index =  0.81
  Adjusted goodness-of-fit index =  0.78
  RMSEA index =  0.08   90% CI: (0.072, 0.088)
  Bentler-Bonnett NFI =  0.8
  Tucker-Lewis NNFI =  0.85
  Bentler CFI =  0.86
  BIC =  -667

  Normalized Residuals
Min. 1st Qu.  MedianMean 3rd Qu.Max.
-2.6300 -0.5640 -0.0728 -0.0067  0.5530  3.5500

  Parameter Estimates
 Estimate Std Error z value Pr(>|z|)
param1  0.78 0.073 10.70.0e+00  Q1 <--- G
param2  0.79 0.065 12.10.0e+00  Q2 <--- G
param3  0.63 0.073  8.50.0e+00  Q3 <--- G
param4  0.74 0.066 11.20.0e+00  Q4 <--- G
param6  0.81 0.068 11.90.0e+00  Q6 <--- G
param7  0.65 0.060 11.00.0e+00  Q7 <--- G
param9  0.60 0.059 10.10.0e+00  Q9 <--- G
param10 0.64 0.065  9.90.0e+00  Q10 <--- G
param11 0.72 0.054 13.30.0e+00  Q11 <--- G
param12 0.59 0.063  9.30.0e+00  Q12 <--- G
param13 0.61 0.069  8.70.0e+00  Q13 <--- G
param14 0.70 0.074  9.60.0e+00  Q14 <--- G
param15 0.68 0.066 10.40.0e+00  Q15 <--- G
param16 0.75 0.056 13.30.0e+00  Q16 <--- G
param17 0.86 0.060 14.30.0e+00  Q17 <--- G
param18 0.63 0.059 10.70.0e+00  Q18 <--- G
param19 0.75 0.062 12.20.0e+00  Q19 <--- G
param20 0.68 0.060 11.40.0e+00  Q20 <--- G
param21 0.64 0.068  9.30.0e+00  Q21 <--- G
param22 0.63 0.065  9.70.0e+00  Q22 <--- G
param23 0.71 0.065 10.90.0e+00  Q23 <--- G
param24 0.70 0.052 13.70.0e+00  Q24 <--- G
param25 0.41 0.066  6.33.4e-10  Q25 <--- G
param26 0.98 0.091 10.80.0e+00  Q1 <--> Q1
param27 0.72 0.068 10.60.0e+00  Q2 <--> Q2
param28 1.09 0.099 11.00.0e+00  Q3 <--> Q3
param29 0.77 0.072 10.70.0e+00  Q4 <--> Q4
param31 0.79 0.075 10.60.0e+00  Q6 <--> Q6
param32 0.64 0.059 10.70.0e+00  Q7 <--> Q7
param34 0.66 0.061 10.80.0e+00  Q9 <--> Q9
param35 0.79 0.073 10.80.0e+00  Q10 <--> Q10
param36 0.45 0.043 10.40.0e+00  Q11 <--> Q11
param37 0.79 0.072 10.90.0e+00  Q12 <--> Q12
param38 0.96 0.088 11.00.0e+00  Q13 <--> Q13
param39 1.05 0.096 10.90.0e+00  Q14 <--> Q14
param40 0.80 0.074 10.80.0e+00  Q15 <--> Q15
param41 0.49 0.047 10.40.0e+00  Q16 <--> Q16
param42 0.51 0.050 10.10.0e+00  Q17 <--> Q17
param43 0.63 0.059 10.80.0e+00  Q18 <--> Q18
param44 0.64 0.060 10.60.0e+00  Q19 <--> Q19
param45 0.63 0.059 10.70.0e+00  Q20 <--> Q20
param46 0.91 0.084 10.90.0e+00  Q21 <--> Q21
param47 0.82 0.076 10.90.0e+00  Q22 <--> Q22
param48 0.77 0.071 10.80.0e+00  Q23 <--> Q23
param49 0.39 0.038 10.30.0e+00  Q24 <--> Q24
param50 0.95 0.086 11.10.0e+00  Q25 <--> Q25
param51 0.29 0.047  6.34.0e-10  Q9 <--> Q7

  Iterations =  16
 > standardized.coefficients(semModif.tmp, digits=2)
 Std. Estimate
param1  param1  0.62  Q1 <--- G
param2  param2  0.68  Q2 <--- G
param3  param3  0.51  Q3 <--- G
param4  param4  0.64  Q4 <--- G
param6  param6  0.67  Q6 <--- G
param7  param7  0.63  Q7 <--- G
param9  param9  0.59  Q9 <--- G
param10 param10 0.58  Q10 <--- G
param11 param11 0.73  Q11 <--- G
param12 param12 0.55  Q12 <--- G
param13 param13 0.53  Q13 <--- G
param14 param14 0.57  Q14 <--- G
param15 param15 0.61  Q15 <--- G
param16 param16 0.73  Q16 <--- G
param17 param17 0.77  Q17 <--- G
param18 param18 0.62  Q18 <--- G
param19 param19 0.69  Q19 <--- G
param20 param20 0.65  Q20 <--- G
param21 param21 0.55  Q21 <--- G
param22 param22 0.57  Q22 <--- G
param23 param23 0.63  Q23 <--- G
param24 param24 0.75  Q24 <--- G
param25 param25 0.39  Q25 <--- G
 >

-- 
Mathieu d'Acremont, PhD [EMAIL PROTECTED]
Maître-Assistanttel/fax +4122 379 98 20/44

Pôle de Recherche National en Sciences Affectives
CISA - Université de Genève
Rue des Battoirs 7
CH-1205 Genève
http://affect.unige.ch/

__
R-help@stat.

Re: [R] summary polr

[Paolo Accadia]

Thank you very much
Paolo

>>> Prof Brian Ripley <[EMAIL PROTECTED]> 20/02/07 4:00 PM >>>
That's not it (the function is 'coef' not 'coeff', and R can tell 
functions and lists apart).

If you read the help page for polr you will see you could have used 
Hess=TRUE.  It works then.  THAT is why we needed an example, to see how 
you used the function.

On Tue, 20 Feb 2007, Michael Dewey wrote:

> At 14:41 20/02/2007, you wrote:
> Please do not just reply to me,
> 1 - I might not know
> 2 - it breaks the threading
>
>> Hi
>>
>> here there is an example extracted from polr help in MASS:
>>
>> The function could be:
>>
>>temp <- function(form, dat) {
>> house.plr <- polr(formula =
>> form, weights = Freq, data = dat)
>> coeff <- summary(house.plr)
>>return(coeff)}
>
> Why do you try to redefine the coeff extractor function?
> Try calling the results of summary summ or some
> other obvious name and I think you will find the problem goes away.
>
> See also
> > library(fortunes)
> > fortune("dog")
>
> Firstly, don't call your matrix 'matrix'. Would you call your dog 'dog'?
> Anyway, it might clash with the function 'matrix'.
>-- Barry Rowlingson
>   R-help (October 2004)
>
>
>
>> the function can be called by:
>>
>>   temp(Sat ~ Infl + Type + Cont, housing)
>>
>> where all data is available from MASS, as it is
>> an example in R Help on 'polr'.
>>
>> Results are:
>>
>>Re-fitting to get Hessian
>>
>>Error in eval(expr, envir, enclos) : object "dat" not found
>>
>> Paolo Accadia
>>
>>
> Michael Dewey <[EMAIL PROTECTED]> 20/02/07 1:43 PM >>>
>> At 15:21 19/02/2007, Paolo Accadia wrote:
>>> Hi all,
>>>
>>> I have a problem to estimate Std. Error and
>>> t-value by “polr† in library Mass.
>> s.
>>> They result from the summary of a polr object.
>>>
>>> I can obtain them working in the R environment
>> with the following statements:
>>>
>>>  temp <- polr(formula = formula1,  data = data1)
>>>  coeff <- summary(temp),
>>>
>>> but when the above statements are enclosed in a
>>> function, summary reports the following error:
>>>
>>> Error in eval(expr, envir, enclos) : object "dat" not found
>>>
>>> Someone knows how I can solve the problem?
>>
>> By giving us a short reproducible example?
>>
>> Specifically we do not know:
>> 1 - what formula1 is
>> 2 - what the structure of data1 is
>> 3 - what the enclosing function looks like
>> 4 - what dat is
>>
>>
>>> Thanks for any help.
>>> Paolo
>>
>> Michael Dewey
>> http://www.aghmed.fsnet.co.uk
>
> Michael Dewey
> http://www.aghmed.fsnet.co.uk
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] testing slopes

[Indermaur Lukas]

Hello
Instead of testing against 0 i would like to test regression slopes against -1. 
Any idea if there's an R script (package?) available.
 
Thanks for any hint.
Cheers
Lukas
 
 
 
 
°°° 
Lukas Indermaur, PhD student 
eawag / Swiss Federal Institute of Aquatic Science and Technology 
ECO - Department of Aquatic Ecology
Überlandstrasse 133
CH-8600 Dübendorf
Switzerland
 
Phone: +41 (0) 71 220 38 25
Fax: +41 (0) 44 823 53 15 
Email: [EMAIL PROTECTED]
www.lukasindermaur.ch

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Help with xlab and ylab distance from axes

[Michael Kubovy]

Dear r-helpers,

In basic graphics, I have a figure with x and y axes suppresed. I  
would like to move the xlab and the ylab closer to the axes. Do I  
have to use mtext()?
_
Professor Michael Kubovy
University of Virginia
Department of Psychology
USPS: P.O.Box 400400Charlottesville, VA 22904-4400
Parcels:Room 102Gilmer Hall
 McCormick RoadCharlottesville, VA 22903
Office:B011+1-434-982-4729
Lab:B019+1-434-982-4751
Fax:+1-434-982-4766
WWW:http://www.people.virginia.edu/~mk9y/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Mahalanobis distance and probability of group membership using Hotelling's T2 distribution

[Mike White]

I want to calculate the probability that a group will include a particular
point using the squared Mahalanobis distance to the centroid. I understand
that the squared Mahalanobis distance is distributed as chi-squared but that
for a small number of random samples from a multivariate normal population
the Hotellings T2 (T squared) distribution should be used.
I cannot find a function for Hotelling's T2 distribution in R (although from
a previous post I have been provided with functions for the Hotelling Test).
My understanding is that the Hotelling's T2 distribution is related to the F
distribution using the equation:
 T2(u,v) = F(u, v-u+1)*vu/(v-u+1)
where u is the number of variables and v the number of group members.

I have written the R code below to compare the results from the chi-squared
distribution with the Hotelling's T2 distribution for probability of a
member being included within a group.
Please can anyone confirm whether or not this is the correct way to use
Hotelling's T2 distribution for probability of group membership. Also, when
testing a particular group member, is it preferable to leave that member out
when calculating the centre and covariance of the group for the Mahalanobis
distances?

Thanks
Mike White



## Hotelling T^2 distribution function
ph<-function(q, u, v, ...){
# q vector of quantiles as in function pf
# u number of independent variables
# v number of observations
if (!v > u+1) stop("n must be greater than p+1")
df1 <- u
df2 <- v-u+1
pf(q*df2/(v*u), df1, df2, ...)
}

# compare Chi-squared and Hotelling T^2 distributions for a group member
u<-3
v<-10
set.seed(1)
mat<-matrix(rnorm(v*u), nrow=v, ncol=u)
MD2<-mahalanobis(mat, center=colMeans(mat), cov=cov(mat))
d<-MD2[order(MD2)]
# select a point midway between nearest and furthest from centroid
dm<-d[length(d)/2]
1-ph(dm,u,v)# probability using Hotelling T^2 distribution
# [1] 0.6577069
1-pchisq(dm, u) # probability using Chi-squared distribution
# [1] 0.5538466

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] summary polr

[Prof Brian Ripley]

That's not it (the function is 'coef' not 'coeff', and R can tell 
functions and lists apart).


If you read the help page for polr you will see you could have used 
Hess=TRUE.  It works then.  THAT is why we needed an example, to see how 
you used the function.


On Tue, 20 Feb 2007, Michael Dewey wrote:


At 14:41 20/02/2007, you wrote:
Please do not just reply to me,
1 - I might not know
2 - it breaks the threading


Hi

here there is an example extracted from polr help in MASS:

The function could be:

   temp <- function(form, dat) {
house.plr <- polr(formula =
form, weights = Freq, data = dat)
coeff <- summary(house.plr)
   return(coeff)}


Why do you try to redefine the coeff extractor function?
Try calling the results of summary summ or some
other obvious name and I think you will find the problem goes away.

See also
> library(fortunes)
> fortune("dog")

Firstly, don't call your matrix 'matrix'. Would you call your dog 'dog'?
Anyway, it might clash with the function 'matrix'.
   -- Barry Rowlingson
  R-help (October 2004)




the function can be called by:

  temp(Sat ~ Infl + Type + Cont, housing)

where all data is available from MASS, as it is
an example in R Help on 'polr'.

Results are:

   Re-fitting to get Hessian

   Error in eval(expr, envir, enclos) : object "dat" not found

Paolo Accadia



Michael Dewey <[EMAIL PROTECTED]> 20/02/07 1:43 PM >>>

At 15:21 19/02/2007, Paolo Accadia wrote:

Hi all,

I have a problem to estimate Std. Error and
t-value by “polr† in library Mass.

s.

They result from the summary of a polr object.

I can obtain them working in the R environment

with the following statements:


 temp <- polr(formula = formula1,  data = data1)
 coeff <- summary(temp),

but when the above statements are enclosed in a
function, summary reports the following error:

Error in eval(expr, envir, enclos) : object "dat" not found

Someone knows how I can solve the problem?


By giving us a short reproducible example?

Specifically we do not know:
1 - what formula1 is
2 - what the structure of data1 is
3 - what the enclosing function looks like
4 - what dat is



Thanks for any help.
Paolo


Michael Dewey
http://www.aghmed.fsnet.co.uk


Michael Dewey
http://www.aghmed.fsnet.co.uk

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.



--
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] R: Re: summary polr

[Paolo Accadia]

Hi all,

The problem is that when you try to use the function summary of a polr object 
in a function, it does not work.
The problem is not related to the formula or the structure of data involved.
It is probably related to the use of the function "vcov" in the code of summary 
for polr, and the iterative procedure to estimate the Hessian.

Anyway, here there is an example extracted from polr help in MASS:

The function could be:

  temp <- function(form, dat) {
   house.plr <- polr(formula = form, weights = Freq, data = dat)
   summ <- summary(house.plr)
  return(summ)}

the function can be called by:

 temp(Sat ~ Infl + Type + Cont, housing)

where all data is available from MASS, as it is an example in R Help on 'polr'.

Results are:

  Re-fitting to get Hessian

  Error in eval(expr, envir, enclos) : object "dat" not found

Paolo Accadia

>>> Michael Dewey <[EMAIL PROTECTED]> 20/02/07 13.43 >>>
At 15:21 19/02/2007, Paolo Accadia wrote:
>Hi all,
>
>I have a problem to estimate Std. Error and 
>t-value by “polr† in library Mass.
>They result from the summary of a polr object.
>
>I can obtain them working in the R environment with the following statements:
>
>  temp <- polr(formula = formula1,  data = data1)
>  coeff <- summary(temp),
>
>but when the above statements are enclosed in a 
>function, summary reports the following error:
>
>Error in eval(expr, envir, enclos) : object "dat" not found
>
>Someone knows how I can solve the problem?

By giving us a short reproducible example?

Specifically we do not know:
1 - what formula1 is
2 - what the structure of data1 is
3 - what the enclosing function looks like
4 - what dat is


>Thanks for any help.
>Paolo

Michael Dewey
http://www.aghmed.fsnet.co.uk

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Use R source in .net (R Source in .net einbinden)

[Martin]

Hi!

 

I would like to use an existing R-code in an .net project. Now I realize
that there exists a DCOM interface which enables me to use the Statconnector
object in my c# project. That way I can get the results of terms using the
.evaluate() method. 

My question is: Is it somehow possible to take a whole file of R code and
have it evaluated, or is it possible to sort of compile the R source in a
way that it could be called as a method from c# .net (for example compile to
dll). It does not seem so, but I thought it couldn't hurt to ask ;-)

Maybe someone had a similar problem and could offer me advice.

Thanks!

 

Hallo!

Ich würde gerne einen existierenden R code in ein .net projekt (C#)
einbinden. Übers DCOM Interface und die Evaluate Funktion kann ich ja
einzelne Ausdrücke auswerten lassen, aber ist es vielleicht möglich
irgendwie ein ganzes code file zu übergeben oder vielleicht aus einem R code
eine dll oder ähnliches zu kompilieren, die man dann in C# .net einbinden
könnte?

Hab dazu nichts gefunden, aber dachte es könnte nicht schaden zu fragen. 

Vielleicht hat ja jemand schon was in der Art gemacht und kann mir einen Rat
geben.

Jedenfalls, danke im Voraus!


[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] summary polr

[Michael Dewey]

At 14:41 20/02/2007, you wrote:
Please do not just reply to me,
1 - I might not know
2 - it breaks the threading

>Hi
>
>here there is an example extracted from polr help in MASS:
>
>The function could be:
>
>temp <- function(form, dat) {
> house.plr <- polr(formula = 
> form, weights = Freq, data = dat)
> coeff <- summary(house.plr)
>return(coeff)}

Why do you try to redefine the coeff extractor function?
Try calling the results of summary summ or some 
other obvious name and I think you will find the problem goes away.

See also
 > library(fortunes)
 > fortune("dog")

Firstly, don't call your matrix 'matrix'. Would you call your dog 'dog'?
Anyway, it might clash with the function 'matrix'.
-- Barry Rowlingson
   R-help (October 2004)



>the function can be called by:
>
>   temp(Sat ~ Infl + Type + Cont, housing)
>
>where all data is available from MASS, as it is 
>an example in R Help on 'polr'.
>
>Results are:
>
>Re-fitting to get Hessian
>
>Error in eval(expr, envir, enclos) : object "dat" not found
>
>Paolo Accadia
>
>
> >>> Michael Dewey <[EMAIL PROTECTED]> 20/02/07 1:43 PM >>>
>At 15:21 19/02/2007, Paolo Accadia wrote:
> >Hi all,
> >
> >I have a problem to estimate Std. Error and
> >t-value by “polr† in library Mass.
>s.
> >They result from the summary of a polr object.
> >
> >I can obtain them working in the R environment 
> with the following statements:
> >
> >  temp <- polr(formula = formula1,  data = data1)
> >  coeff <- summary(temp),
> >
> >but when the above statements are enclosed in a
> >function, summary reports the following error:
> >
> >Error in eval(expr, envir, enclos) : object "dat" not found
> >
> >Someone knows how I can solve the problem?
>
>By giving us a short reproducible example?
>
>Specifically we do not know:
>1 - what formula1 is
>2 - what the structure of data1 is
>3 - what the enclosing function looks like
>4 - what dat is
>
>
> >Thanks for any help.
> >Paolo
>
>Michael Dewey
>http://www.aghmed.fsnet.co.uk

Michael Dewey
http://www.aghmed.fsnet.co.uk

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Reading Post-Script files

[Ted Harding]

On 20-Feb-07 Ralf Finne wrote:
> Hi everybody!
> Is there any way to read a postscrit file into R?
> 
> All the best to you
> Ralf Finne
> SYH University of Applied Sciences
> Vasa Finland

Well, yes ... since a PostScript file is ASCII text you could
use readline() ... but what you'd do with it after that I cannot
imagine! [So I'm joking here]

So, to come to the point: Why do you want to do that?

If, for example, the PS file displays tabular data which you
want to read into R as data, then perhaps one way would be to
use suitable utility software to convert the PS file to a PDF
file. Then you can display the PDF file using Acrobat Reader,
and use the mouse to copy the data into a file which you can
then edit up so that it can be read nicely into a dataframe
in R.

There are also utilities called ps2ascii (part of ghostscript;
often fails to do a clean job) and pstotext which you can get
from

  http://www.cs.wisc.edu/~ghost/doc/pstotext.htm

(which claims to do a better job, though I've not tested it)
which can convert a PS file into an ASCII text file which
gives essentially the same text layout as seen when you
display/print the PS file.

You may then be able to edit this into a form which R can read
as you want.

But, for better targeted advice, it would be useful to know
why you want to "read a PDF file into R"!

Best wishes,
Ted.


E-Mail: (Ted Harding) <[EMAIL PROTECTED]>
Fax-to-email: +44 (0)870 094 0861
Date: 20-Feb-07   Time: 15:11:34
-- XFMail --

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Reading Post-Script files

[Roger Bivand]

On Tue, 20 Feb 2007, Ralf Finne wrote:

> Hi everybody!
> Is there any way to read a postscrit file into R?

See http://www.r-project.org/useR-2006/Slides/Murrell.pdf

page 4, the grImport package, now on CRAN, with further notes on Paul 
Murrell's home page:

http://www.stat.auckland.ac.nz/~paul/

> 
> All the best to you
> Ralf Finne
> SYH University of Applied Sciences
> Vasa Finland
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: [EMAIL PROTECTED]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Installing Package rgl - Compilation Fails

[Duncan Murdoch]

On 2/20/2007 9:44 AM, Rick Bilonick wrote:
> On Mon, 2007-02-19 at 15:11 -0600, Ranjan Maitra wrote:
>> The error is different now. It now cannot find Xext library. Do a yum search 
>> on this and install that.
>> 
>> yum provides libXext
>> 
>> which will give you the package Xext which needs to be installed.
>> 
>> You may also need to install the libXext-devel.i386 package.
>> 
>> HTH.
>> Ranjan
>> 
>> 
> 
> Thanks. Installing libXext-devel allowed rgl to install.

Would you mind summarizing what you started from, and what extras you 
needed to install?  I'd like to add this information to the rgl docs.

Duncan Murdoch

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] "contingency table" for several variables

[Frank E Harrell Jr]

David LEVY wrote:
> DearList,
> 
> 
> 
> I ‘m trying to draw ONE table that summarize SEVERAL categorical variables
> according to one classification variable, say “sex”. The result would look
> like several contingency tables appended one to the other. All the variables
> belong to a data frame.
> 
> The summary.formula in Hmisc package does something pretty close and is
> ready for a Latex export  but I need either to get rid off the percentage
> (or put the count prior to the percentage )in the “reverse” option or to add
> a chisquare test in the “response” method.

Not a good idea.  The % is normalized for sample size differences so 
should be emphasized.

Frank

> 
> 
> 
> The result would look like the one of
> 
>> summary(sex~treatment+Symptoms, fun = table, method = “response”)
> 
> in the help of  summary.formula but with chisquare tests attached.
> 
> Or :
> 
>> summary(sex~treatment+Symptoms, fun = table, method = “reverse”, test= T)
> 
>  gives all the information, but I can’t use it for its form is not
> appropriate.
> 
> 
> 
> Is there any package where I could find a solution ? Any way to create a
> function that would add a test in the “response” method ?
> 
> Otherwise, ftable should help but I don’t know how to use it with a
> data.frame.
> 
> 
> 
> Thanks for your help.
> 
> Regards,
> 
> David
> 
>   [[alternative HTML version deleted]]
> 
> 
> 
> 
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.


-- 
Frank E Harrell Jr   Professor and Chair   School of Medicine
  Department of Biostatistics   Vanderbilt University

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Installing Package rgl - Compilation Fails

[Rick Bilonick]

On Mon, 2007-02-19 at 15:11 -0600, Ranjan Maitra wrote:
> The error is different now. It now cannot find Xext library. Do a yum search 
> on this and install that.
> 
> yum provides libXext
> 
> which will give you the package Xext which needs to be installed.
> 
> You may also need to install the libXext-devel.i386 package.
> 
> HTH.
> Ranjan
> 
> 

Thanks. Installing libXext-devel allowed rgl to install.

Rick B.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] "contingency table" for several variables

[David LEVY]

DearList,



I ‘m trying to draw ONE table that summarize SEVERAL categorical variables
according to one classification variable, say “sex”. The result would look
like several contingency tables appended one to the other. All the variables
belong to a data frame.

The summary.formula in Hmisc package does something pretty close and is
ready for a Latex export  but I need either to get rid off the percentage
(or put the count prior to the percentage )in the “reverse” option or to add
a chisquare test in the “response” method.



The result would look like the one of

> summary(sex~treatment+Symptoms, fun = table, method = “response”)

in the help of  summary.formula but with chisquare tests attached.

Or :

> summary(sex~treatment+Symptoms, fun = table, method = “reverse”, test= T)

 gives all the information, but I can’t use it for its form is not
appropriate.



Is there any package where I could find a solution ? Any way to create a
function that would add a test in the “response” method ?

Otherwise, ftable should help but I don’t know how to use it with a
data.frame.



Thanks for your help.

Regards,

David

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] contextstack overflow

[Steven Finch]

Hello!

I written several implementations in R of Rémy's algorithm for
generating random ordered strictly binary trees on 2n+1 vertices.

One implementation involves manufacturing character strings like:

"X <- list(0,list(0,0))"

for the case n=2.  If I perform the following two steps:

cmd <- "X <- list(0,list(0,0))"
eval(parse(text=cmd))

then X becomes a true nested list in R.  This works fine for n=2,
but often for n=200, an error message:

Error in parse(text = cmd) : contextstack overflow

appears and execution stops.  Clearly there exists an upper bound
on the allowable depth of nestings in R!  Can this upper bound be
easily increased?

Other implementations avoid this problem, so this issue is not
crucial to me.  I do wish, however, to understand the limits of
this particular approach.  Thank you!

Steve Finch
http://algo.inria.fr/bsolve/

P.S.  If anyone else has written R code for generating random
trees (on a fixed number of vertices), I would enjoy seeing this!

_
Don’t miss your chance to WIN 10 hours of private jet travel from Microsoft® 
Office Live http://clk.atdmt.com/MRT/go/mcrssaub0540002499mrt/direct/01/


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Reshape (pivot) question

[Lauri Nikkinen]

Thanks guys! Jim got it right what I was looking for. I think that the
reason why I don't get this right with cast-function is that I don't know
how to write a proper stats-function. SQL is still much easier for me.

Lauri

2007/2/20, hadley wickham <[EMAIL PROTECTED]>:
>
> > Haven't quite learned to 'cast' yet, but I have always used the 'apply'
> > functions for this type of processing:
> >
> > > x <- "id patient_id date code class eala
> > + ID1564262 1562 6.4.200612:00  1 NA
> > + ID1564262 1562 6.4.200612:00  1 NA
> > + ID1564264 1365 14.2.200614:35  1 50
> > + ID1564265 1342 7.4.200614:30  2 50
> > + ID1564266 1648 7.4.200614:30  2 50
> > + ID1564267 1263 10.2.200615:45  2 10
> > + ID1564267 1263 10.2.200615:45  3 10
> > + ID1564269 5646 13.5.200617:02  3 10
> > + ID1564270 7561 13.5.200617:02  1 10
> > + ID1564271 1676 15.5.200620:41  2 20"
> > >
> > > x.in <- read.table(textConnection(x), header=TRUE)
> > > # 'by' seems to drop NAs so convert to a character string for
> processing
> > > x.in$eala <- ifelse(is.na(x.in$eala), "NA", as.character(x.in$eala))
> > > # convert date to POSIXlt so we can access the year and month
> > > myDate <- strptime(x.in$date, "%d.%m.%Y%H:%M")
> > > x.in$year <- myDate$year + 1900
> > > x.in$month <- myDate$mon+1
>
> To do this with the reshape package, all you need is:
>
> > x.in$patient_id <- factor(x.in$patient_id)
> > dfm <- melt(x.in, id=c("code", "month", "year"), m=c("id","patient_id"))
> > cast(dfm, code + month + year ~ variable, stats)
> code month year id_n id_uniikit patient_id_n patient_id_uniikit
> 1  2 20061  11  1
> 2  4 20062  22  2
> 3  5 20061  11  1
> 4  2 20061  11  1
> 5  5 20061  11  1
> 6  2 20061  11  1
> 7  4 20062  12  1
> 8  5 20061  11  1
>
> Which looks like what you want from your R code, but not your SQL.
>
> Hadley
>

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Reading Post-Script files

[Ralf Finne]

Hi everybody!
Is there any way to read a postscrit file into R?

All the best to you
Ralf Finne
SYH University of Applied Sciences
Vasa Finland

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Sample size

[Frank E Harrell Jr]

SAULEAU Erik-André wrote:
> Dear R-list,
> 
> I have to design the validation of a score (ordinal values between 0 and 6) 
> reputed to separate 4 groups of patients with known frequencies in the 
> population. I think the more accurate is to calculate sample size under 
> median test. Is there a function for that in R (not in the pwr package)?
> 
> Thanks in advance, with all my best, erik S.
> 
> 
> 
> 
> =
> Erik-André SAULEAU
>  
> SEAIM
> Centre Hospitalier
> 87, Avenue d'Altkirch
> BP 1070
> 68051 Mulhouse Cédex
> Tel: 03 89 64 67 53
> Mel: [EMAIL PROTECTED]
> Web: www.ch-mulhouse.fr

The median test has low power and should be avoided.  For a validation 
study a hypothesis test is seldom of interest and you might base the 
sample size on the needed precision (margin of error; confidence 
interval width) of some estimate.

-- 
Frank E Harrell Jr   Professor and Chair   School of Medicine
  Department of Biostatistics   Vanderbilt University

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] memory management uestion [Broadcast]

[Federico Calboli]

Liaw, Andy wrote:
> I don't see why making copies of the columns you need inside the loop is
> "better" memory management.  If the data are in a matrix, accessing
> elements is quite fast.  If you're worrying about speed of that, do what
> Charles suggest: work with the transpose so that you are accessing
> elements in the same column in each iteration of the loop.

As I said, this is pretty academic, I am not looking for how to do something 
differetly.

Having said that, let me present this code:

for(i in gp){
new[i,1] = ifelse(srow[i]>0, new[srow[i],zippo[i]], sav[i])
new[i,2] = ifelse(drow[i]>0, new[drow[i],zappo[i]], sav[i])
  }

where gp is large vector and srow and drow are the dummy variables for:

srow = data[,2]
drow = data[,4]

If instead of the dummy variable I access the array directly (and its' a 60 
x 6 array) the loop takes 2/3 days --not sure here, I killed it after 48 hours.

If I use dummy variables the code runs in 10 minutes-ish.

Comments?

Best,

Fede

-- 
Federico C. F. Calboli
Department of Epidemiology and Public Health
Imperial College, St Mary's Campus
Norfolk Place, London W2 1PG

Tel  +44 (0)20 7594 1602 Fax (+44) 020 7594 3193

f.calboli [.a.t] imperial.ac.uk
f.calboli [.a.t] gmail.com

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] summary polr

[Michael Dewey]

At 15:21 19/02/2007, Paolo Accadia wrote:
>Hi all,
>
>I have a problem to estimate Std. Error and 
>t-value by “polr” in library Mass.
>They result from the summary of a polr object.
>
>I can obtain them working in the R environment with the following statements:
>
>  temp <- polr(formula = formula1,  data = data1)
>  coeff <- summary(temp),
>
>but when the above statements are enclosed in a 
>function, summary reports the following error:
>
>Error in eval(expr, envir, enclos) : object "dat" not found
>
>Someone knows how I can solve the problem?

By giving us a short reproducible example?

Specifically we do not know:
1 - what formula1 is
2 - what the structure of data1 is
3 - what the enclosing function looks like
4 - what dat is


>Thanks for any help.
>Paolo

Michael Dewey
http://www.aghmed.fsnet.co.uk

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Reshape (pivot) question

[hadley wickham]

> Haven't quite learned to 'cast' yet, but I have always used the 'apply'
> functions for this type of processing:
>
> > x <- "id patient_id date code class eala
> + ID1564262 1562 6.4.200612:00  1 NA
> + ID1564262 1562 6.4.200612:00  1 NA
> + ID1564264 1365 14.2.200614:35  1 50
> + ID1564265 1342 7.4.200614:30  2 50
> + ID1564266 1648 7.4.200614:30  2 50
> + ID1564267 1263 10.2.200615:45  2 10
> + ID1564267 1263 10.2.200615:45  3 10
> + ID1564269 5646 13.5.200617:02  3 10
> + ID1564270 7561 13.5.200617:02  1 10
> + ID1564271 1676 15.5.200620:41  2 20"
> >
> > x.in <- read.table(textConnection(x), header=TRUE)
> > # 'by' seems to drop NAs so convert to a character string for processing
> > x.in$eala <- ifelse(is.na(x.in$eala), "NA", as.character(x.in$eala))
> > # convert date to POSIXlt so we can access the year and month
> > myDate <- strptime(x.in$date, "%d.%m.%Y%H:%M")
> > x.in$year <- myDate$year + 1900
> > x.in$month <- myDate$mon+1

To do this with the reshape package, all you need is:

> x.in$patient_id <- factor(x.in$patient_id)
> dfm <- melt(x.in, id=c("code", "month", "year"), m=c("id","patient_id"))
> cast(dfm, code + month + year ~ variable, stats)
  code month year id_n id_uniikit patient_id_n patient_id_uniikit
1  2 20061  11  1
2  4 20062  22  2
3  5 20061  11  1
4  2 20061  11  1
5  5 20061  11  1
6  2 20061  11  1
7  4 20062  12  1
8  5 20061  11  1

Which looks like what you want from your R code, but not your SQL.

Hadley

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] memory management uestion [Broadcast]

[Liaw, Andy]

I don't see why making copies of the columns you need inside the loop is
"better" memory management.  If the data are in a matrix, accessing
elements is quite fast.  If you're worrying about speed of that, do what
Charles suggest: work with the transpose so that you are accessing
elements in the same column in each iteration of the loop.

Andy 

From: Federico Calboli
> 
> Charles C. Berry wrote:
> 
> > Whoa! You are accessing one ROW at a time.
> > 
> > Either way this will tangle up your cache if you have many rows and 
> > columns in your orignal data.
> > 
> > You might do better to do
> > 
> > Y <- t( X ) ### use '<-' !
> > 
> > for (i in whatever ){
> > do something using Y[ , i ]
> > }
> 
> My question is NOT how to write the fastest code, it is 
> whether dummy variables (for lack of better words) make the 
> memory management better, i.e. faster, or not.
> 
> Best,
> 
> Fede
> 
> --
> Federico C. F. Calboli
> Department of Epidemiology and Public Health Imperial 
> College, St Mary's Campus Norfolk Place, London W2 1PG
> 
> Tel  +44 (0)20 7594 1602 Fax (+44) 020 7594 3193
> 
> f.calboli [.a.t] imperial.ac.uk
> f.calboli [.a.t] gmail.com
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 
> 


--
Notice:  This e-mail message, together with any attachments,...{{dropped}}

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Sample size

[SAULEAU Erik-André]

Dear R-list,

I have to design the validation of a score (ordinal values between 0 and 6) 
reputed to separate 4 groups of patients with known frequencies in the 
population. I think the more accurate is to calculate sample size under median 
test. Is there a function for that in R (not in the pwr package)?

Thanks in advance, with all my best, erik S.




=
Erik-André SAULEAU
 
SEAIM
Centre Hospitalier
87, Avenue d'Altkirch
BP 1070
68051 Mulhouse Cédex
Tel: 03 89 64 67 53
Mel: [EMAIL PROTECTED]
Web: www.ch-mulhouse.fr
=
Registre des Cancers du Haut-Rhin
9, Rue du Dr Mangeney
BP 1370
68070 Mulhouse Cédex
Tel: 03 89 64 62 51
Web: www.arer68.org
=


**
Afin d'eviter toute propagation de virus informatique, et en complement 
des dispositifs en place, ce message (et ses pieces jointes s'il y en a) 
a ete automatiquement analyse par un antivirus de messagerie. 
**


[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Reshape (pivot) question

[jim holtman]

Haven't quite learned to 'cast' yet, but I have always used the 'apply'
functions for this type of processing:

> x <- "id patient_id date code class eala
+ ID1564262 1562 6.4.200612:00  1 NA
+ ID1564262 1562 6.4.200612:00  1 NA
+ ID1564264 1365 14.2.200614:35  1 50
+ ID1564265 1342 7.4.200614:30  2 50
+ ID1564266 1648 7.4.200614:30  2 50
+ ID1564267 1263 10.2.200615:45  2 10
+ ID1564267 1263 10.2.200615:45  3 10
+ ID1564269 5646 13.5.200617:02  3 10
+ ID1564270 7561 13.5.200617:02  1 10
+ ID1564271 1676 15.5.200620:41  2 20"
>
> x.in <- read.table(textConnection(x), header=TRUE)
> # 'by' seems to drop NAs so convert to a character string for processing
> x.in$eala <- ifelse(is.na(x.in$eala), "NA", as.character(x.in$eala))
> # convert date to POSIXlt so we can access the year and month
> myDate <- strptime(x.in$date, "%d.%m.%Y%H:%M")
> x.in$year <- myDate$year + 1900
> x.in$month <- myDate$mon+1
> # split the data by eala, year, month and summarize
> x.by <- by(x.in, list(x.in$eala, x.in$year, x.in$month), function(x){
+ data.frame(eala=x$eala[1], month=x$month[1], year=x$year[1],
+ icount=length(unique(x$id)), pcount=length(unique(x$patient_id)),
+ count1=sum(x$class == 1), count2=sum(x$class == 2),
count3=sum(x$class == 3))
+ })
> # convert back to a data frame
> do.call(rbind, x.by)
  eala month year icount pcount count1 count2 count3
1   10 2 2006  1  1  0  1  1
2   50 2 2006  1  1  1  0  0
3   50 4 2006  2  2  0  2  0
4   NA 4 2006  1  1  2  0  0
5   10 5 2006  2  2  1  0  1
6   20 5 2006  1  1  0  1  0
>
>



On 2/20/07, Lauri Nikkinen <[EMAIL PROTECTED]> wrote:
>
> Hi R-users,
>
> I have a data set like this (first ten rows):
>
>id patient_id date code class eala ID1564262 1562 6.4.2006 12:00  1
> NA ID1564262 1562 6.4.2006 12:00  1 NA ID1564264 1365 14.2.2006 14:35
>  1 50 ID1564265 1342 7.4.2006 14:30  2 50 ID1564266 1648
> 7.4.200614:30
>  2 50 ID1564267 1263 10.2.2006 15:45  2 10 ID1564267 1263
> 10.2.200615:45
>  3 10 ID1564269 5646 13.5.2006 17:02  3 10 ID1564270 7561
> 13.5.200617:02
>  1 10 ID1564271 1676 15.5.2006 20:41  2 20
>
> How can I do a new (pivot?) data.frame in R which I can achieve by MS SQL:
>
> select  eala,
> datepart(month, date) as month,
> datepart(year, date) as year,
> count(distinct id) as id_count,
> count(distinct patient_id) as patient_count,
> count(distinct(case when class = 1 then code else null end)) as count_1,
> count(distinct(case when class = 2 then code else null end)) as count_2,
> count(distinct(case when class = 3 then code else null end)) as count_3,
> into temp2
> from temp1
> group by datepart(month, date), datepart(year, date), eala
> order by datepart(month, date), datepart(year, date), eala
>
> I tried something like this but could not go further:
>
> stats <- function(x) {
>count <- function(x) length(na.omit(x))
>c(
>  n = count(x),
>  uniikit = length(unique(x))
>  )
> }
> library(reshape)
> attach(dframe)
> dfm <- melt(dframe, measure.var=c("id","patient_id"), id.var=c
> ("code",""this
> should be month"",""this should be year), variable_name="variable")
>
> cast(dfm, code + month + year ~ variable, stats)
>
> Regards,
>
> Lauri
>
>[[alternative HTML version deleted]]
>
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
Jim Holtman
Cincinnati, OH
+1 513 646 9390

What is the problem you are trying to solve?

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Problems with obtaining t-tests of regression coefficients applying consistent standard errors after run 2SLS estimation. Clearer !!!!!

[Guillermo Julián San Martín]

First I have to say I am sorry because I have not been so clear in my
previous e-mails. I will try to explain clearer what it is my problem.

I have the following model:



lnP=Sc+Ag+Ag2+Var+R+D



In this model the variable Sc is endogenous and the rest are all objective
exogenous variables. I verified that Sc is endogenous through a standard
Hausman test. To determine this I defined before a new instrumental
variable, I2. Also I detected through a Breusch Pagan Test a problem of
heteroskedasticity.

With the intention to avoid the problem of the endogenous variable and the
heteroskedasticity I want to apply first the technique 2SLS and then based
in these results I want to obtain the t-tests of the coefficients applying
Heteroskedasticity Consistent Standard Errors (HCSE) or Huber-White errors.

Like I showed above I have just one structural equation in the model. In
this situation, to apply 2SLS in R until I know there two possible ways:
First to use the function tsls() from package sem, or second, to use the
function systemfit() from package systemfit. I thought that systemfit  was
for situations when there are more than one structural equation in the
model. Anyway I probed with the two ways and I obtained similar results.
Below, I show the program lines:



*Program lines 1:*

* *

*> First option: applying tsls *

*>library (sem)*

*>Reg1 <-tsls (LnP~Sc+Ag+Ag2+Var+R+D,~I2+Ag+Ag2+Var+R+D)  *

*>summary (Reg1)*

* *

*> Second option: applying systemfit *

*>library (systemfit)*

*>RS <- LnP~Sc+Ag+Ag2+Var+R+D   # structural equation*

*>Inst <- ~I2+Ag+Ag2+Var+R+D# instrumental variables*

*>labels <- list("RS")*

*>system <- list(RS)*

*>Reg2 <- systemfit("2SLS", system, labels, Inst, saveMemory=TRUE) *

*>summary (Reg2)*



Now I want to obtain the t-tests of the coefficients but applying the HCSE.
I know two different ways to obtain them: First, applying the function
robcov() from package Design, or Second, applying coeftest() and vcovHV()
from packages lmtest and sandwich. The program lines are the following:



*Program lines 2:*



*> First option: using robcov() *

*>library (Design)*

*>options(scipen=20)*
*>robcov(Reg1) ### I tried with Reg 2 too*
**
**
**
***>Second option: using coeftest and vcovHC *

*>library (lmtest)*

*>library (sandwich)*

*>coeftest (Reg1, vcov=vcovHC(Reg1 or Reg2, type="HC0"))  ### I tried with
Reg 2 too*



In the two cases after trying to apply robcov or coeftest I obtained a
message of error:



*With robcov:*

*> Error in rep.default(1, p) : rep() incorrect type for second argument *

* *

*With coeftest:*

*>** **Error in terms.default(object) : no terms component*



In the following lines I copy all the sequence together from R. If I try
with systemfit instead of tsls I obtain the same result:



*Program lines 3:*

*>ibrary (sem)*

*> Reg1 <-tsls(LnP~Sc+Ag+Ag2+Var+R+D,~I2+Ag+Ag2+Var+R+D)*

*> library (Design)*

*Loading required package: Hmisc*

*Loading required package: chron*

*Hmisc library by Frank E Harrell Jr*

*Type library(help='Hmisc'), ?Overview, or ?Hmisc.Overview')*

*to see overall documentation.*

*NOTE:Hmisc no longer redefines [.factor to drop unused levels when*

*subsetting.  To get the old behavior of Hmisc type dropUnusedLevels().*

*Loading required package: survival*

*Loading required package: splines*

*Attaching package: 'survival'*

*The following object(s) are masked from package:Hmisc :*

* untangle.specials *

*Design library by Frank E Harrell Jr*

*Type library(help='Design'), ?Overview, or ?Design.Overview')*

*to see overall documentation.*

*Attaching package: 'Design'*

*The following object(s) are masked from package:survival :*

* cox.zph Surv survfit *

*The following object(s) are masked from package:Hmisc :*

* .noGenenerics .R. .SV4. *

*> options(scipen=20)*

*> robcov(Reg1)*

*Error in rep.default(1, p) : rep() incorrect type for second argument*

*> library (lmtest)*

*Loading required package: zoo*

*Attaching package: 'lmtest'*

*The following object(s) are masked from package:Design :*

* lrtest *

*> library (sandwich)*

*> coeftest (Reg1, vcov=vcovHC(Reg1, type="HC0"))*

*Error in terms.default(object) : no terms component*

*>*

I would like to solve these problems and to understand what the meanings of
the errors messages are. Also I want to know if there is another way to
obtain the HCSE through a function, because if I can not fix the problem my
last alternative is to do it manually.

I will thank a lot if somebody can help me

Best wishes

Guillermo

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Lengend function and moving average plot

[John Kane]

--- amna khan <[EMAIL PROTECTED]> wrote:

> Sir I am very new user of R. I am not understanding
> the how to write the
> plot description in box outside the plot. 

Have a look at 
http://finzi.psych.upenn.edu/R/Rhelp02a/archive/68585.html

You should read up on par 
?par 

I hope this helps

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Mahalanobis distance and probability of group membership using Hotelling's T2 distribution

[Mike White]

I want to calculate the probability that a group will include a particular
point using the squared Mahalanobis distance to the centroid. I understand
that the squared Mahalanobis distance is distributed as chi-squared but that
for a small number of random samples from a multivariate normal population
the Hotellings T2 (T squared) distribution should be used.
I cannot find a function for Hotelling's T2 distribution in R (although from
a previous post I have been provided with functions for the Hotelling Test).
My understanding is that the Hotelling's T2 distribution is related to the F
distribution using the equation:
 T2(u,v) = F(u, v-u+1)*vu/(v-u+1)
where u is the number of variables and v the number of group members.

I have written the R code below to compare the results from the chi-squared
distribution with the Hotelling's T2 distribution for probability of a
member being included within a group.
Please can anyone confirm whether or not this is the correct way to use
Hotelling's T2 distribution for probability of group membership. Also, when
testing a particular group member, is it preferable to leave that member out
when calculating the centre and covariance of the group for the Mahalanobis
distances?

Thanks
Mike White



## Hotelling T^2 distribution function
ph<-function(q, u, v, ...){
# q vector of quantiles as in function pf
# u number of independent variables
# v number of observations
if (!v > u+1) stop("n must be greater than p+1")
df1 <- u
df2 <- v-u+1
pf(q*df2/(v*u), df1, df2, ...)
}

# compare Chi-squared and Hotelling T^2 distributions for a group member
u<-3
v<-10
set.seed(1)
mat<-matrix(rnorm(v*u), nrow=v, ncol=u)
MD2<-mahalanobis(mat, center=colMeans(mat), cov=cov(mat))
d<-MD2[order(MD2)]
# select a point midway between nearest and furthest from centroid
dm<-d[length(d)/2]
1-ph(dm,u,v)# probability using Hotelling T^2 distribution
# [1] 0.6577069
1-pchisq(dm, u) # probability using Chi-squared distribution
# [1] 0.5538466

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] tree()

[Prof Brian Ripley]

This is a function of your data and the tuning parameters you chose to 
use.  See ?tree.control.

On Tue, 20 Feb 2007, stephenc wrote:

> Hi
>
> I am trying to use tree() to classify movements in a futures contract.  My
> data is like this:
>
> diff  dip  dim adx
> 1  0100.08650.100.0
> 2  0 93.185402044.5455 93.18540
> 3  0 90.309951549.1169 90.30995
> 4  1 85.22030 927.0419 85.22030
> 5  1 85.36084 785.6480 85.36084
> 6  0 85.72627 663.3814 85.72627
> 7  0 78.06721 500.1113 78.06721
> 8  1 69.59398 376.7558 69.59398
> 9  1 71.15429 307.4533 71.15429
> 10 1 71.81023 280.6238 71.81023
>
> plus another 6000 lines
>
> The cpus example works fine and I am trying this:
>
> tree.model <- tree(as.factor(indi$diff) ~ indi$dim + indi$dip + indi$adx,
> indi[1:4000,])

Oh, please!  use the data= argument properly.

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] Difficulties with dataframe filter using elements from an array created using a for loop or seq()

[Prof Brian Ripley]

FAQ Q7.31

On Tue, 20 Feb 2007, Todd A. Johnson wrote:

> Hi All-
>
> This seems like such a pathetic problem to be posting about, but I have no
> idea why this testcase does not work.  I have tried this using R 2.4.1,
> 2.4.0, 2.3.0, and 2.0.0 on several different computers (Mac OS 10.4.8,
> Windows XP, Linux).  Below the signature, you will find my test case R code.
>
> My point in this folly is to take a dataframe of 300,000 rows, create a
> filter based on two of the rows, and count the number of rows in the
> filtered and unfiltered dataframe.  One column in the filter only has the
> numbers 0.05, 0.15, 0.25, 0.35, 0.45, 0.55, 0.65, 0.75, 0.85, 0.95, so I
> thought that I could just iterate in a for loop and get the job done. Just
> the simple single column filter case is presented here. Obviously, there are
> only ten numbers, so the "manual" method is easy, but I would like to have a
> more flexible program. (Plus it worries me if the simple things don't do
> what I expect... :-) )
>
>> From the output, you can see that the loop using the "handmadevector" that
> creates a filter and counts the elements, correctly finds one match for each
> element in the vector, but the seq() and for loop produced vectors each give
> a mixture of true and false matches.
>
> Can anyone tell me why the "loopvector" and "seqvector" do not provide the
> same output as the "handmadevector".
>
>
> Thank you for your assistance!
>
> Todd
>
>

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

Re: [R] bootstrapping Levene's test

[Chuck Cleland]

[EMAIL PROTECTED] wrote:
> Hello all,
> 
> I am low down on the learning curve of R but so far I have had little  
> trouble using most of the packages. However, recently I have run into  
> a wall when it comes to bootstrapping a Levene's test (from the car  
> package) and thought you might be able to help. I have not been able  
> to find R examples for the "boot" package where the test statistic  
> specifically uses a grouping variable (or at least a simple example  
> with this condition). I would like to do a  non-parametric bootstrap  
> to eventually get 95% confidence intervals using the boot.ci command.  
> I have included the coding I have tried on a simple data set below. If  
> anyone could provide some help, specifically with regards to how the  
> "statistic" arguement should be set up in the boot package, it would  
> be greatly appreciated.
> 
>> library(boot)
>> library(car)
>> data<-c(2,45,555,1,77,1,2,1,2,1)
>> group<-c(1,1,1,1,1,2,2,2,2,2)
>> levene.test(data,group)
> Levene's Test for Homogeneity of Variance
>Df F value Pr(>F)
> group  1  1.6929 0.2294
> 8
>> stat<-function(a){levene.test(a,group)}
>> trial1<-boot(data,statistic,100)
> Error in statistic(data, original, ...) : unused argument(s) ( ...)

  One problem is that levene.test() returns an ANOVA table, which is not
a statistic.  Looking inside levene.test() to see what might be used as
a statistic, for the two-group case it seems like you could do something
like this:

library(boot)
library(car)
set.seed(671969)

df <- data.frame(Y = c(rnorm(100,0,1), rnorm(100,0,1)),
 G = rep(c(1,2), each = 100))

boot.levene <- function(data,indices){
  levene.diff <- diff(tapply(df$Y[indices],
list(df$G[indices]), mad))
  return(levene.diff)
  }

first.try <- boot(df, boot.levene, 1000)

first.try

ORDINARY NONPARAMETRIC BOOTSTRAP

Call:
boot(data = df, statistic = boot.levene, R = 1000)


Bootstrap Statistics :
  original biasstd. error
t1* -0.1944924 0.02439172   0.1753512

boot.ci(first.try, index=1, type="bca")

BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
Based on 1000 bootstrap replicates

CALL :
boot.ci(boot.out = first.try, type = "bca", index = 1)

Intervals :
Level   BCa
95%   (-0.5772,  0.1159 )
Calculations and Intervals on Original Scale

> Best regards,
> Kevin
> 
> __
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

-- 
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Difficulties with dataframe filter using elements from an array created using a for loop or seq()

[Todd A. Johnson]

Hi All-
 
This seems like such a pathetic problem to be posting about, but I have no
idea why this testcase does not work.  I have tried this using R 2.4.1,
2.4.0, 2.3.0, and 2.0.0 on several different computers (Mac OS 10.4.8,
Windows XP, Linux).  Below the signature, you will find my test case R code.
 
My point in this folly is to take a dataframe of 300,000 rows, create a
filter based on two of the rows, and count the number of rows in the
filtered and unfiltered dataframe.  One column in the filter only has the
numbers 0.05, 0.15, 0.25, 0.35, 0.45, 0.55, 0.65, 0.75, 0.85, 0.95, so I
thought that I could just iterate in a for loop and get the job done. Just
the simple single column filter case is presented here. Obviously, there are
only ten numbers, so the "manual" method is easy, but I would like to have a
more flexible program. (Plus it worries me if the simple things don't do
what I expect... :-) )

>From the output, you can see that the loop using the "handmadevector" that
creates a filter and counts the elements, correctly finds one match for each
element in the vector, but the seq() and for loop produced vectors each give
a mixture of true and false matches.

Can anyone tell me why the "loopvector" and "seqvector" do not provide the
same output as the "handmadevector".
 

Thank you for your assistance!

Todd

-- 
Todd A. Johnson
Research Associate, Laboratory for Medical Informatics
SNP Research Center,RIKEN
1-7-22Suehiro,Tsurumi-ku,Yokohama
Kanagawa 230-0045,Japan

Cellphone: 090-5309-5867

E-mail: [EMAIL PROTECTED]



Here's the testcase, with the sample code between the lines and the output
following:
 
_
## Set up three different vectors, each with the numbers 0.05, 0.15, 0.25,
0.35, 0.45, 0.55, 0.65, 0.75, 0.85, 0.95
## each of which is used to select records from a dataframe based on
equality to a particular column
## The first vector is created by using a for loop
loopvector <- c()
for (i in 0:9){
loopvector <- c(loopvector, (i*0.10)+0.05);
}
## The second vector is made "by hand"
handmadevector <- c(0.05, 0.15, 0.25, 0.35, 0.45, 0.55, 0.65, 0.75, 0.85,
0.95)
## The third vector is made using seq()
seqvector <- seq(0.05, 0.95, 0.10)
## Are the vectors the same?
all.equal(loopvector, handmadevector)
all.equal(loopvector, seqvector)
print(handmadevector)
print(loopvector)
print(seqvector)
## As a simple testcase, I create a dataframe with two variables, a varA of
dummy data, and bBins
## which is the column on which I was trying to filter.
a <- c(0,1,2,0,1,3,4,5,3,5)
b <- c(0.05,0.15,0.25,0.35,0.45,0.55,0.65,0.75,0.85,0.95)
testdf <- data.frame(varA = a, bBins = b)
attach(testdf)
## Loop through each of the vectors, create a filter on the dataframe based
on equality with the current iteration,
## and print that number and the count of records in the dataframe that
match that number.
for (i in loopvector){
aqs_filt <- bBins==i;
print(i);
print(length(testdf$varA[aqs_filt]));
}
for (i in handmadevector){
aqs_filt <- bBins==i;
print(i);
print(length(testdf$varA[aqs_filt]));
}
for (i in seqvector){
aqs_filt <- bBins==i;
print(i);
print(length(testdf$varA[aqs_filt]));
}
 
_
 
Here's the output from R 2.4.1 running on an Apple 12" Powerbook.
 
> ## Set up three different vectors, each with the numbers 0.05, 0.15, 0.25,
0.35, 0.45, 0.55, 0.65, 0.75, 0.85, 0.95
> ## each of which is used to select records from a dataframe based on equality
to a particular column
> ## The first vector is created by using a for loop
> loopvector <- c()
> for (i in 0:9){
+ loopvector <- c(loopvector, (i*0.10)+0.05);
+ }
> ## The second vector is made "by hand"
> handmadevector <- c(0.05, 0.15, 0.25, 0.35, 0.45, 0.55, 0.65, 0.75, 0.85,
0.95)
> ## The thirs vector is made using seq()
> seqvector <- seq(0.05, 0.95, 0.10)
> ## Are the vectors the same?
> all.equal(loopvector, handmadevector)
[1] TRUE
> all.equal(loopvector, seqvector)
[1] TRUE
> 
> print(handmadevector)
 [1] 0.05 0.15 0.25 0.35 0.45 0.55 0.65 0.75 0.85 0.95
> print(loopvector)
 [1] 0.05 0.15 0.25 0.35 0.45 0.55 0.65 0.75 0.85 0.95
> print(seqvector)
 [1] 0.05 0.15 0.25 0.35 0.45 0.55 0.65 0.75 0.85 0.95
> ## As a simple testcase, I create a dataframe with two variables, a varA of
dummy data, and bBins
> ## which is the column on which I was trying to filter.
> a <- c(0,1,2,0,1,3,4,5,3,5)
> b <- c(0.05,0.15,0.25,0.35,0.45,0.55,0.65,0.75,0.85,0.95)
> testdf <- data.frame(varA = a, bBins = b)
> attach(testdf)
> ## Loop through each of the vectors, create a filter on the dataframe based on
equality with the current iteration,
> ## and print that number and the count of records in the dataframe that match
that number.
> for (i in loopvector){
+ aqs_filt <- bBins==i;
+ print(i);
+ print(length(testdf$varA[aqs_filt]));
+ }
[1] 0.05
[1] 1
[1] 0.15
[1] 0
[1] 0.25
[1] 1
[1] 0.35
[1] 0
[1] 0.45
[1] 1
[1] 0.55
[1] 1
[1] 0.65
[1] 0
[1

Re: [R] Another subsetting enigma

[Johannes Graumann]

jim holtman wrote:

>> matches <- sapply(result, function(x){
> +     .comma <- strsplit(x, ',')[[1]]  # get the fields
> +     any(my.list %in% .comma)
> + })

Thanks for that!

Joh

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] interaction term and scatterplot3d

[Patrick Giraudoux]

Dear Listers,

I would be interested in representing a trend surface including an 
interaction term z = f(x,y,x.y) - eg the type of chart obtained with 
persp() or wireframe(), then adding datapoints as a cloud, ideally with 
dots which are under the surface in a color, and those who are above in 
another color. An other option would be to draw segments between the 
dots and the ground of the chart.

scatterplot3d looks like being close to do such things except it does 
not to include (to my knowledge) a coefficient for the interaction term 
(thus just model z = f(x,y).

Does anybody has an idea where to start with this and the main steps? Or 
a place/website where some script examples can be found?

Patrick

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.

[R] Reshape (pivot) question

[Lauri Nikkinen]

Hi R-users,

I have a data set like this (first ten rows):

id patient_id date code class eala ID1564262 1562 6.4.2006 12:00  1
NA ID1564262 1562 6.4.2006 12:00  1 NA ID1564264 1365 14.2.2006 14:35
 1 50 ID1564265 1342 7.4.2006 14:30  2 50 ID1564266 1648 7.4.200614:30
 2 50 ID1564267 1263 10.2.2006 15:45  2 10 ID1564267 1263
10.2.200615:45
 3 10 ID1564269 5646 13.5.2006 17:02  3 10 ID1564270 7561
13.5.200617:02
 1 10 ID1564271 1676 15.5.2006 20:41  2 20

How can I do a new (pivot?) data.frame in R which I can achieve by MS SQL:

select  eala,
 datepart(month, date) as month,
 datepart(year, date) as year,
 count(distinct id) as id_count,
 count(distinct patient_id) as patient_count,
 count(distinct(case when class = 1 then code else null end)) as count_1,
 count(distinct(case when class = 2 then code else null end)) as count_2,
 count(distinct(case when class = 3 then code else null end)) as count_3,
into temp2
from temp1
group by datepart(month, date), datepart(year, date), eala
order by datepart(month, date), datepart(year, date), eala

I tried something like this but could not go further:

stats <- function(x) {
count <- function(x) length(na.omit(x))
c(
  n = count(x),
  uniikit = length(unique(x))
  )
}
library(reshape)
attach(dframe)
dfm <- melt(dframe, measure.var=c("id","patient_id"), id.var=c("code",""this
should be month"",""this should be year), variable_name="variable")

cast(dfm, code + month + year ~ variable, stats)

Regards,

Lauri

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.