Re: [R] randomForest

2005-07-11 Thread Martin Maechler
 Duncan == Duncan Murdoch [EMAIL PROTECTED]
 on Thu, 07 Jul 2005 15:44:38 -0400 writes:

Duncan On 7/7/2005 3:38 PM, Weiwei Shi wrote:
 Hi there:
 I have a question on random foresst:
 
 recently i helped a friend with her random forest and i came with this 
problem:
 her dataset has 6 classes and since the sample size is pretty small:
 264 and the class distr is like this (Diag is the response variable)

 sample.size - lapply(1:6, function(i) sum(Diag==i))
 sample.size
 [[1]]
 [1] 36



and later you get problems because you didn't know that a *list*
such as 'sample.size' should be made into a so called
*atomic vector* {and there's a function  is.atomic(.) ! to test for it}
and Duncan and others told you about unlist().

Now there are two things I'd want to add:

1) If you had used 

  s.size - table(Diag)

   you had used a faster and simpler expression with the same result.
   Though in general (when there can be zero counts), to give the
   same result, you'd need

  s.size - table(factor(Diag, levels = 1:6))

   Still a bit preferable to the lapply(.) IMO


2)  You should get into the habit of using
 sapply(.)   rather than  lapply(.).

sapply() originally was exactly devised for the above task:
and stands for ``[s]implified lapply''.

It always returns an ``unlisted'' result when appropriate.

Regards,
Martin Maechler, ETH Zurich

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] choice of graph

2005-07-11 Thread Navarre Sabine
Hi,
 
It's about 2 weeks that I think about a graph to translate my datas. But I 
don't have an really idea.
I 'm going to expose you my problem:

I have a questionnaire with 15 questions and you have more possibilties to 
answer to these.
For example:
 The trainer is competent:   Yes   No
 I have learn a lot at the training:   Bad Quit bad   Medium 
Good  Quit good
 etc

I would like to represent all my question on the same plot and differentiate 
the different type of answer possible.

On my first impression, I had doing a classification of the different type of 
answer and doing a barplot, but my responsable don't want a classification but 
want to see all the question with their own type of answer.

I have attached my first idea!

Can you please help me?


Thanks a lot
 
Sabine



-

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] validation, calibration and Design

2005-07-11 Thread Williams Scott
 

Hi R experts,

 

I am trying to do a prognostic model validation study, using cancer
survival data. There are 2 data sets - 1500 cases used to develop a
nomogram, and another of 800 cases used as an independent validation
cohort.  I have validated the nomogram in the original data (easy with
the Design tools), and then want to show that it also has good results
with the independent data using 60 month survival. I would also like to
show that the nomogram is significantly different to an existing model
based on 60 month survival data generated by it (eg by McNemar's test).

Hence, somewhat shortened:   

 

#using R 2.01 on Windows

library(Hmisc)

library(Design)

 

data1 #dataframe with predictor variables A and B, cens and time 

  columns (months)

ddist1 - datadist(data1) 

options(datadist='ddist1') 

 

s1 - Surv(data1$time, data1$cens)

 

cph.nomo - cph(s1 ~ A+B, surv=T, x=T, y=T, time.inc=60)

 

survcph - Survival(cph.nomo, x=T, y=T, time.inc=60, surv=T)

surv5 - function(lp) survcph(60, lp)

nomogram(cph.nomo, lp=T, conf.int=F, fun=list(surv5, surv7), 

funlabel=c(5 yr DFS))

 

# now have a useful nomogram model, with good discrimination and

#calibration when checked with validate and calibrate (not shown)

#move on to validation cohort of n=800

 

Data2 #Validation data with same predictor variables A, B, cens, time

# do I need to put data2 into datadist??

 

s2 - Surv(data2$time, data2$cens)

 

#able to derive 60 month estimates of survival using

data2.est5 - survest(cph.nomo, expand.grid(A=data2$A, B=data2$B), 

times=c(60), conf.int=0)

 

rcorr.cens(data2.est5$surv, s2) # tests discrimination of the model 

#against the validation data observed censored data

 

# I cant find a way to use calibrate in this setting though??

# Also, if I have the 5 year estimates for 2 different models, I can 

# use rcorr.cens to show discrimination, but which values are 

# suitable for a test of difference (eg with McNemars)?

# I have tried predict / newdata function a number of ways but it 

# typically returns an error relating to unequal vector lengths

 

 

What I cant work out is where to go now to derive a calibration curve of
the predicted 5 year result (val.data5) and the observed  (s2). Or can I
do it another way? For example, could I merge the 2 data frames and use
lines1:1500 to build the model and the last 800 lines to validate?

 

Obviously I am a novice, and sure to be missing something simple. I have
spent countless hours pouring over Prof Harrell's text (which is great
but doesn't have a specific example of this) and Design Help plus the R
news archive with no success, so any help is very much appreciated. 

 

Scott Williams MD

Peter MacCallum Cancer Centre

Melbourne Australia

 


[[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


Re: [R] time series regression

2005-07-11 Thread Achim Zeileis
On Fri, 8 Jul 2005, yyan liu wrote:

 Hi:
   I have two time series y(t) and x(t). I want to
 regress Y on X. Because Y is a time series and may
 have autocorrelation such as AR(p), so it is not
 efficient to use OLS directly. The model I am trying
 to fit is like
 Y(t)=beta0+beta1*X(t)+rho*Y(t-1)+e(t)

 e(t) is iid normal random error. Anybody know whether
 there is a function in R can fit such models? The
 function can also let me specify how many beta's and
 rho's I can have in the model.

If you want to estimate the model by ML, you can use arima() and specify
further regressors via the `xreg' argument.
Estimation by OLS can be done via lm(), but that typically requires
setting up the lags yourself. More convenient interfaces are provided in
the `dyn' package by Gabor Grothendieck and my `dynlm' package.
Z

 Thx a lot!

 liu

 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] choice of graph

2005-07-11 Thread TEMPL Matthias
 Hi,
  
 It's about 2 weeks that I think about a graph to translate my 
 datas. But I don't have an really idea. I 'm going to expose 
 you my problem:
 
 I have a questionnaire with 15 questions and you have more 
 possibilties to answer to these. For example:
  The trainer is competent:   Yes   No
  I have learn a lot at the training:   Bad Quit bad   
 Medium Good  Quit good
  etc
 
 I would like to represent all my question on the same plot 
 and differentiate the different type of answer possible.
 
 On my first impression, I had doing a classification of the 
 different type of answer and doing a barplot, 

Did you try a mosaicplot too?
?mosaicplot

Best,
Matthias

 but my 
 responsable don't want a classification but want to see all 
 the question with their own type of answer.
 
 I have attached my first idea!
 
 Can you please help me?
 
 
 Thanks a lot
  
 Sabine
 
 
   
 -
 


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] (no subject)

2005-07-11 Thread Smit, R. \(Robin\)
Hello,
 
The estimate of glm dispersion can be based on the deviance or on the
Pearson statistic.
I have compared output from R glm() to another statastical package and
it appears that R uses the Pearson statistic.
I was wondering if it is possible to make use R the deviance instead by
modifying the glm(...) function?
Thanks for your attention.
 
Kind regards,
Robin Smit

 

This e-mail and its contents are subject to the DISCLAIMER at 
http://www.tno.nl/disclaimer/email.html
[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Problems with R on OS X

2005-07-11 Thread Heinz Schild
I used R on OS X 10.3x quite some time with no serious problems.  
Sometimes R stopped when I tried to execute a bigger program. After  
updating to OS X to version 10.4 R worked but I still had the problem  
with bigger programs. Therefore I re-installed R on top of the  
existing R version. The installation finished properly but suddenly R  
did not work. Then I reinstalled OS X 10.4 because I thought some  
left over registry information from R may caused the problem. R still  
crashed. I hoped the problem would be overcome with the new R version  
1.12. Unfortunately this is not the case. The error report (see  
appendix) says that Thread 0 crashed.

What can I do?
Heinz Schild

Date/Time:  2005-07-08 18:13:16.748 +0200
OS Version: 10.4.1 (Build 8B15)
Report Version: 3

Command: R
Path:/Applications/R.app/Contents/MacOS/R
Parent:  WindowServer [62]

Version: 1.12 (1622)

PID:269
Thread: 0

Exception:  EXC_BREAKPOINT (0x0006)
Code[0]:0x0001
Code[1]:0x92897200


Thread 0 Crashed:
0   com.apple.Foundation0x92897200 _NSRaiseError + 264
1   com.apple.Foundation0x92896f3c +[NSException raise:format:] + 40
2   com.apple.Foundation0x9286c120 -[NSString stringByAppendingString:] 
+ 104
3   org.R-project.R 0x6ce0 -[RController doLoadHistory:] + 812 
(crt.c:300)
4   org.R-project.R 0x3ce0 -[RController awakeFromNib] + 1980 
(crt.c:300)
5   com.apple.Foundation0x92886788 -[NSSet makeObjectsPerformSelector:] 
+ 164
6   com.apple.AppKit0x93636f94 -[NSIBObjectData 
nibInstantiateWithOwner:topLevelObjects:] + 864
7   com.apple.AppKit0x93623324 loadNib + 240
8   com.apple.AppKit0x93622d7c +[NSBundle(NSNibLoading) 
_loadNibFile:nameTable:withZone:ownerBundle:] + 716
9   com.apple.AppKit0x9367a05c +[NSBundle(NSNibLoading) 
loadNibFile:externalNameTable:withZone:] + 156
10  com.apple.AppKit0x93709e10 +[NSBundle(NSNibLoading) 
loadNibNamed:owner:] + 344
11  com.apple.AppKit0x93709bb0 NSApplicationMain + 344
12  org.R-project.R 0x2dd0 _start + 392 (crt.c:267)
13  dyld0x8fe01048 _dyld_start + 60

Thread 1:
0   libSystem.B.dylib   0x9000a778 mach_msg_trap + 8
1   libSystem.B.dylib   0x9000a6bc mach_msg + 60
2   com.apple.CoreFoundation0x9074a4d8 __CFRunLoopRun + 832
3   com.apple.CoreFoundation0x90749ddc CFRunLoopRunSpecific + 268
4   com.apple.Foundation0x928716b8 -[NSConnection sendInvocation:] + 
1468
5   com.apple.Foundation0x92870154 -[NSObject(NSForwardInvocation) 
forward::] + 408
6   libobjc.A.dylib 0x909b10d0 _objc_msgForward + 176
7   org.R-project.R 0x4244 -[RController readThread:] + 244 
(crt.c:300)
8   com.apple.Foundation0x9287c2b4 forkThreadForFunction + 108
9   libSystem.B.dylib   0x9002c3d4 _pthread_body + 96

Thread 0 crashed with PPC Thread State:
  srr0: 0x92897200 srr1: 0x0202f030vrsave: 0x
cr: 0x24002422  xer: 0x0004   lr: 0x928971d8  ctr: 0x9285903c
r0: 0x   r1: 0xb1f0   r2: 0xa2856508   r3: 0xbfffeda0
r4: 0x   r5: 0x92858964   r6: 0xbfffee24   r7: 0x00ff
r8: 0xbfffee10   r9: 0x03c09f20  r10: 0x909aa668  r11: 0x24002422
   r12: 0x9285903c  r13: 0x  r14: 0x  r15: 0x
   r16: 0x0002  r17: 0x00043524  r18: 0x00043524  r19: 0x00043524
   r20: 0x00043524  r21: 0x00043524  r22: 0x00043524  r23: 0x00043524
   r24: 0xa365e07c  r25: 0x000369b4  r26: 0x01175740  r27: 0x03c2e7d0
   r28: 0x03c37c30  r29: 0xa285b7ec  r30: 0x90a239c0  r31: 0x92897108

Binary Images Description:
0x1000 -0x38fff org.R-project.R 1.12 (1622) 
/Applications/R.app/Contents/MacOS/R
   0x52000 -0x84fff libreadline.5.0.dylib   
/Library/Frameworks/R.framework/Resources/lib/libreadline.5.0.dylib
   0x98000 -0xc3fff libncurses.5.dylib  /usr/lib/libncurses.5.dylib
  0x205000 -   0x392fff libR.dylib  
/Library/Frameworks/R.framework/Resources/lib/libR.dylib
0x8fe0 - 0x8fe50fff dyld 43 /usr/lib/dyld
0x9000 - 0x901a6fff libSystem.B.dylib   /usr/lib/libSystem.B.dylib
0x901fe000 - 0x90202fff libmathCommon.A.dylib   
/usr/lib/system/libmathCommon.A.dylib
0x90204000 - 0x90257fff com.apple.CoreText 1.0.0 (???)  
/System/Library/Frameworks/ApplicationServices.framework/Versions/A/Frameworks/CoreText.framework/Versions/A/CoreText
0x90284000 - 0x90335fff ATS 
/System/Library/Frameworks/ApplicationServices.framework/Versions/A/Frameworks/ATS.framework/Versions/A/ATS
0x90364000 - 0x9069cfff com.apple.CoreGraphics 1.256.5 (???)
/System/Library/Frameworks/ApplicationServices.framework/Versions/A/Frameworks/CoreGraphics.framework/Versions/A/CoreGraphics
0x90727000 - 0x90800fff com.apple.CoreFoundation 6.4.1 (368.1)  
/System/Library/Frameworks/CoreFoundation.framework/Versions/A/CoreFoundation
0x90849000 - 0x90849fff 

[R] R on kubuntu

2005-07-11 Thread Constant Depièreux
Hello all,

I am planning to redeploy my workstation under KUBUNTU.

Does any body has any r experience installing/using r on this platform?

Best regards.


-- 
Constant Depièreux
Managing Director
Applied QUality Technologies Europe sprl
Rue des Déportés 123, B-4800 Verviers
(Tel) +32 87 292175 - (Fax) +32 87 292171 - (Mobile) +32 475 555 818
(Web) http://www.aqte.be - (Courriel) [EMAIL PROTECTED]
(Skype) cdepiereux

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] R on kubuntu

2005-07-11 Thread Peter Dalgaard
Constant Depièreux [EMAIL PROTECTED] writes:

 Hello all,
 
 I am planning to redeploy my workstation under KUBUNTU.
 
 Does any body has any r experience installing/using r on this platform?

I don't think it's any different than any other Debian-derived
platform. 

-- 
   O__   Peter Dalgaard Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark  Ph:  (+45) 35327918
~~ - ([EMAIL PROTECTED])  FAX: (+45) 35327907

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] R on kubuntu

2005-07-11 Thread ernesto
Constant Depièreux wrote:

Hello all,

I am planning to redeploy my workstation under KUBUNTU.

Does any body has any r experience installing/using r on this platform?

Best regards.


  


Hi,

I've just moved from SuSE 9.1 to Ubuntu 5.04. The instalation is clean
and easy, although you have to waste some time selecting packages.

After I got all the packages installed I just compiled R like usual and
it is working very well.

Regards

EJ

ps: I have to tell you that Ubuntu was a very good surprise for me.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] plot(cox.zph()): customize xlab ylab

2005-07-11 Thread Dan Bebber
Hello,

plot(cox.zph(my.ph),var=1,xlab=Year)

gives the error:
Error in plot.default(range(xx), yr, type = n, xlab
= Time, ylab = ylab[i],  : formal argument xlab
matched by multiple actual arguments

How can I customize the xlab and ylab for plots of
cox.zph?

Thanks,
Dan Bebber

Department of Plant Sciences
University of Oxford
UK



___ 
How much free photo storage do you get? Store your holiday

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] estVARXar parameter significance

2005-07-11 Thread Kemp S E \(Comp\)
Hi,

Does anyone know how to check the significance (p-values/t-values) of the 
estimated parameters via the estVARXar function (in the dse package)?

Thanks in advance,

Sam.

[[alternative HTML version deleted]]

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] calculating dispersion formula using deviance ( was Re: (no subject) )

2005-07-11 Thread Adaikalavan Ramasamy
Please try to use a meaningful subject line. See below for comments.



On Mon, 2005-07-11 at 11:30 +0200, Smit, R. (Robin) wrote:
 Hello,
  
 The estimate of glm dispersion can be based on the deviance or on the
 Pearson statistic.
 I have compared output from R glm() to another statastical package and
 it appears that R uses the Pearson statistic.

A quick search would also highlight the following thread
http://www.r-project.org/nocvs/mail/r-help/2002/6938.html

 I was wondering if it is possible to make use R the deviance instead by
 modifying the glm(...) function?

I don't know what the formula for using the deviance is but _IF_ it is
the square root of ratio of null deviance by its degrees of freedom,
then sqrt( fit$deviance / fit$df.null ) should be useful.


 Thanks for your attention.
  
 Kind regards,
 Robin Smit
 
  
 
 This e-mail and its contents are subject to the DISCLAIMER at 
 http://www.tno.nl/disclaimer/email.html
   [[alternative HTML version deleted]]
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] plot(cox.zph()): customize xlab ylab

2005-07-11 Thread Adaikalavan Ramasamy
I am not sure if there is an easy way around this. An ugly hack is to
make a copy the function survival:::plot.cox.zph and make your
modified function. But there are others in the list who might know
neater solutions.

Regards, Adai


On Mon, 2005-07-11 at 11:10 +0100, Dan Bebber wrote:
 Hello,
 
 plot(cox.zph(my.ph),var=1,xlab=Year)
 
 gives the error:
 Error in plot.default(range(xx), yr, type = n, xlab
 = Time, ylab = ylab[i],  : formal argument xlab
 matched by multiple actual arguments
 
 How can I customize the xlab and ylab for plots of
 cox.zph?
 
 Thanks,
 Dan Bebber
 
 Department of Plant Sciences
 University of Oxford
 UK
 
 
   
 ___ 
 How much free photo storage do you get? Store your holiday
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] R on kubuntu

2005-07-11 Thread Lefebure Tristan
No problem at all.
If you allow universe packages, many binary R packages are available (from 
the GNU/Linux Debian sid). 


On Monday 11 July 2005 11:54, Constant Depièreux wrote:
 Hello all,

 I am planning to redeploy my workstation under KUBUNTU.

 Does any body has any r experience installing/using r on this platform?

 Best regards.

-- 

Tristan LEFEBURE
Laboratoire d'écologie des hydrosystèmes fluviaux (UMR 5023)
Université Lyon I - Campus de la Doua
Bat. Darwin C 69622 Villeurbanne - France

Phone: (33) (0)4 26 23 44 02
Fax: (33) (0)4 72 43 15 23

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Sweave and complex numbers

2005-07-11 Thread Robin Hankin
Hi

When using Sweave, most of my functions get called with complex  
arguments.

They get typeset in with additions that I don't want; 1+1i appears  
as 1 + (0 + 1i)
and I would rather have plain old 1+1i.

Example follows:




\documentclass[a4paper]{article}

\title{A Test File}
\author{Robin Hankin}

\usepackage{a4wide}

\begin{document}

\maketitle

A simple example:
print=TRUE=
f - function(x){
x+c(1,1+1i)}
f(6+7i)
@


Question: why does \verb=c(1,1+1i)= get printed as
\verb=c(1,1 + (0+1i))= ?

And how can I stop it?

\end{document}




can anyone help me here?







--
Robin Hankin
Uncertainty Analyst
National Oceanography Centre, Southampton
European Way, Southampton SO14 3ZH, UK
  tel  023-8059-7743

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] polr (MASS) link functions

2005-07-11 Thread I M S White
When we analyse ordered categorical data, with categories II,
should we expect an equivalent analysis if the ordering is changed to
II ?  Not with a cloglog link (see example below). However, I
suspect the cloglog analysis on the reordered categories is equivalent
to a loglog analysis on the original ordering.  Perhaps for consistency
polr should include a loglog link? Or warn that a different answer is
available, possibly more appropriate?

Cheese tasting example from McCullagh and Nelder:

 library(MASS)
 frqs - c(0,0,1,7,8,8,19,8,1,
+ 6,9,12,11,7,6,1,0,0,
+ 1,1,6,8,23,7,5,1,0,
+ 0,0,0,1,3,7,14,16,11)
 Y1 - ordered(rep(1:9,4))
 Y2 -  ordered(rep(9:1,4))
 trt - factor(rep(c(A,B,C,D),each=9))
 fit1 - polr(Y1 ~ trt, method=clog,weight=frqs)
 fit2 - polr(Y2 ~ trt, method=clog,weight=frqs)
 print(fit1)
Call:
polr(formula = Y1 ~ trt, weights = frqs, method = clog)

Coefficients:
 trtB  trtC  trtD
-1.809487 -0.879997  0.987389

Intercepts:
1|2 2|3 3|4 4|5 5|6 6|7
-2.50668784 -2.04055560 -1.42592187 -0.80166219 -0.01722569  0.55135471
7|8 8|9
 1.55688848  2.84302173

Residual Deviance: 725.2341
AIC: 747.2341

 print(fit2)
Call:
polr(formula = Y2 ~ trt, weights = frqs, method = clog)

Coefficients:
  trtB   trtC   trtD
 1.9513455  1.1033902 -0.9588354

Intercepts:
   1|22|33|44|55|66|77|8
-1.3832854 -0.5216350  0.4225419  0.9994969  1.8317538  2.5777980  3.4532818
   8|9
 4.4026602

Residual Deviance: 718.0052
AIC: 740.0052


==
I.White
University of Edinburgh
Ashworth Laboratories, West Mains Road
Edinburgh EH9 3JT
Tel: 0131 650 5490  Fax: 0131 650 6564
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


Re: [R] plot(cox.zph()): customize xlab ylab

2005-07-11 Thread Duncan Murdoch
Adaikalavan Ramasamy wrote:
 I am not sure if there is an easy way around this. An ugly hack is to
 make a copy the function survival:::plot.cox.zph and make your
 modified function. But there are others in the list who might know
 neater solutions.

This hack is uglier (and might not work properly on some devices, but 
you could try):

plot(cox.zph(my.ph),var=1,col.lab=white)
title(xlab=Year, ylab=Beta for blah blah blah)

Duncan Murdoch
 
 Regards, Adai
 
 
 On Mon, 2005-07-11 at 11:10 +0100, Dan Bebber wrote:
 
Hello,

plot(cox.zph(my.ph),var=1,xlab=Year)

gives the error:
Error in plot.default(range(xx), yr, type = n, xlab
= Time, ylab = ylab[i],  : formal argument xlab
matched by multiple actual arguments

How can I customize the xlab and ylab for plots of
cox.zph?

Thanks,
Dan Bebber

Department of Plant Sciences
University of Oxford
UK


  
___ 
How much free photo storage do you get? Store your holiday

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

 
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Sweave and complex numbers

2005-07-11 Thread Duncan Murdoch
Robin Hankin wrote:
 Hi
 
 When using Sweave, most of my functions get called with complex  
 arguments.
 
 They get typeset in with additions that I don't want; 1+1i appears  
 as 1 + (0 + 1i)
 and I would rather have plain old 1+1i.
 
 Example follows:
 
 
 
 
 \documentclass[a4paper]{article}
 
 \title{A Test File}
 \author{Robin Hankin}
 
 \usepackage{a4wide}
 
 \begin{document}
 
 \maketitle
 
 A simple example:
 print=TRUE=
 f - function(x){
 x+c(1,1+1i)}
 f(6+7i)
 @
 
 
 Question: why does \verb=c(1,1+1i)= get printed as
 \verb=c(1,1 + (0+1i))= ?
 
 And how can I stop it?
 
 \end{document}
 
 
 
 
 can anyone help me here?

The R parser only understands pure imaginary constants as complex 
numbers.  It parses 1+1i as the sum of the real constant 1 and the 
complex constant 0+1i.

This isn't easy to work around.  Sweave could special case these, or the 
parser or deparser could, but it looks messy to me.  I'd guess the best 
solution would be if it happened in the parser (i.e. a real constant 
plus an imaginary constant was folded into a complex constant), but 
right now our parser doesn't do any sorts of optimizations like this. 
It's not easy to add the first one, not least because the parser doesn't 
know for sure that 1+1i really is a complex number:  you might have 
redefined the meaning of +.

Maybe someone else has an answer to your 2nd question.

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


[R] Misbehaviour of DSE

2005-07-11 Thread Ajay Shah

Folks,

I am finding problems with using dse:

 library(dse1)
Loading required package: tframe
Error: c(package '%s' required by '%s' could not be found, setRNG, dse1)
 library(dse2)
Loading required package: setRNG
Error: package 'setRNG' could not be loaded
In addition: Warning message:
there is no package called 'setRNG' in: library(pkg, character.only = TRUE, 
logical = TRUE, lib.loc = lib.loc) 

This is on R 2.1 on an Apple ibook (OS X) panther.

Thanks,

-ans.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Boxplot philosophy {was Boxplot in R}

2005-07-11 Thread Martin Maechler
 AdaiR == Adaikalavan Ramasamy [EMAIL PROTECTED]
 on Mon, 11 Jul 2005 03:04:44 +0100 writes:

AdaiR Just an addendum on the philosophical aspect of doing
AdaiR this.  By selecting the 5% and 95% quantiles, you are
AdaiR always going to get 10% of the data as extreme and
AdaiR these points may not necessarily outliers.  So when
AdaiR you are comparing information from multiple columns
AdaiR (i.e.  boxplots), it is harder to say which column
AdaiR contains more extreme value compared to others etc.

Yes, indeed!

People {and software implementations} have several times provided
differing definitions of how the boxplot whiskers should be defined.

I strongly believe that this is very often a very bad idea!!

A boxplot should be a universal mean communication and so one
should be *VERY* reluctant redefining the outliers.

I just find that Matlab (in their statistics toolbox)
does *NOT* use such a silly 5% / 95% definition of the whiskers,
at least not according to their documentation.
That's very good (and I wonder where you, Larry, got the idea of
the 5 / 95 %).
Using such a fixed percentage is really a very inferior idea to
John Tukey's definition {the one in use in all implementations
of S (including R) probably for close to 20 years now}.

I see one flaw in Tukey's definition {which is shared of course
by any silly percentage based ``outlier'' definition}:

   The non-dependency on the sample size.

If you have a 1000 (or even many more) points,
you'll get more and more `outliers' even for perfectly normal data.

But then, I assume John Tukey would have told us to do more
sophisticated things {maybe things like the violin plots} than
boxplot  if you have really very many data points, you may want
to see more features -- or he would have agreed to use 
   boxplot(*,  range = monotone_slowly_growing(n) )
for largish sample sizes n.

Martin Maechler, ETH Zurich




AdaiR Regards, Adai

AdaiR On Sun, 2005-07-10 at 18:10 -0500, Larry Xie wrote:
 I am trying to draw a plot like Matlab does: 
 
 The upper extreme whisker represents 95% of the data;
 The upper hinge represents 75% of the data;
 The median represents 50% of the data;
 The lower hinge represents 25% of the data;
 The lower extreme whisker represents 5% of the data.
 
 It looks like:
 
 --- 95%
 |
 |
 ---   75%
 | |
 |-|   50%
 | |
 | |
 ---   25%
 |
 --- 5%
 
 Anyone can give me some hints as to how to draw a boxplot like that?
 What function does it? I tried boxplot() but couldn't figure it out.
 If it's boxplot(), what arguments should I pass to the function? Thank
 you for your help. I'd appreciate it.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] class- vs. as()

2005-07-11 Thread stefan . albrecht




Dear all,

I would appreciate a lot, if someone could explain to me in a simple
way, why the assignment class- is not always working and one has to
take as() like in the example below.

 (v - matrix(1:9, 3))
 [,1] [,2] [,3]
[1,]147
[2,]258
[3,]369
 class(v)
[1] matrix

 class(v) - integer
 class(v)
[1] matrix
 v
 [,1] [,2] [,3]
[1,]147
[2,]258
[3,]369

 (vi - as(v, integer))
[1] 1 2 3 4 5 6 7 8 9
 class(vi)
[1] integer


With many thanks,

Stefan
[[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


Re: [R] validation, calibration and Design

2005-07-11 Thread Frank E Harrell Jr
Williams Scott wrote:
  
 
 Hi R experts,
 
  
 
 I am trying to do a prognostic model validation study, using cancer
 survival data. There are 2 data sets - 1500 cases used to develop a
 nomogram, and another of 800 cases used as an independent validation
 cohort.  I have validated the nomogram in the original data (easy with
 the Design tools), and then want to show that it also has good results
 with the independent data using 60 month survival. I would also like to
 show that the nomogram is significantly different to an existing model
 based on 60 month survival data generated by it (eg by McNemar's test).

Scott,

A nomogram is a graphical device, not a model to validate.  It merely 
represents a model.

If the 800 subjects came from the same hospitals in roughly the same 
time era, you are doing an internal validation and this is an 
exceedingly inefficient way to do it.  Not only is this wasting 800 
subjects from developing the model, but the validation sample is not 
large enough to yield reliable accuracy estimates.  And I don't see how 
McNemar's test applies as nothing is binary about this problem.

If you have two models that have identical degrees of overfitting (e.g. 
were based on the same number of CANDIDATE degrees of freedom) you can 
use the rcorrp.cens function to test for differences in discrimination 
or paired predictions.

If you really have an external sample (say 800 subjects from another 
country) you can use the groupkm function in Design to get a 
Kaplan-Meier-based calibration curve.  Otherwise I would recombine the 
data, develop the model on all subjects you can get, and use the 
bootstrap to validate it.

Frank

 
 Hence, somewhat shortened:   
 
  
 
 #using R 2.01 on Windows
 
 library(Hmisc)
 
 library(Design)
 
  
 
 data1 #dataframe with predictor variables A and B, cens and time 
 
   columns (months)
 
 ddist1 - datadist(data1) 
 
 options(datadist='ddist1') 
 
  
 
 s1 - Surv(data1$time, data1$cens)
 
  
 
 cph.nomo - cph(s1 ~ A+B, surv=T, x=T, y=T, time.inc=60)
 
  
 
 survcph - Survival(cph.nomo, x=T, y=T, time.inc=60, surv=T)
 
 surv5 - function(lp) survcph(60, lp)
 
 nomogram(cph.nomo, lp=T, conf.int=F, fun=list(surv5, surv7), 
 
 funlabel=c(5 yr DFS))
 
  
 
 # now have a useful nomogram model, with good discrimination and
 
 #calibration when checked with validate and calibrate (not shown)
 
 #move on to validation cohort of n=800
 
  
 
 Data2 #Validation data with same predictor variables A, B, cens, time
 
 # do I need to put data2 into datadist??
 
  
 
 s2 - Surv(data2$time, data2$cens)
 
  
 
 #able to derive 60 month estimates of survival using
 
 data2.est5 - survest(cph.nomo, expand.grid(A=data2$A, B=data2$B), 
 
 times=c(60), conf.int=0)
 
  
 
 rcorr.cens(data2.est5$surv, s2) # tests discrimination of the model 
 
 #against the validation data observed censored data
 
  
 
 # I cant find a way to use calibrate in this setting though??
 
 # Also, if I have the 5 year estimates for 2 different models, I can 
 
 # use rcorr.cens to show discrimination, but which values are 
 
 # suitable for a test of difference (eg with McNemars)?
 
 # I have tried predict / newdata function a number of ways but it 
 
 # typically returns an error relating to unequal vector lengths
 
  
 
  
 
 What I cant work out is where to go now to derive a calibration curve of
 the predicted 5 year result (val.data5) and the observed  (s2). Or can I
 do it another way? For example, could I merge the 2 data frames and use
 lines1:1500 to build the model and the last 800 lines to validate?
 
  
 
 Obviously I am a novice, and sure to be missing something simple. I have
 spent countless hours pouring over Prof Harrell's text (which is great
 but doesn't have a specific example of this) and Design Help plus the R
 news archive with no success, so any help is very much appreciated. 
 
  
 
 Scott Williams MD
 
 Peter MacCallum Cancer Centre
 
 Melbourne Australia
 
  
 
 
   [[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
 


-- 
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


[R] Off topic -2 Ln Lambda and Chi square

2005-07-11 Thread Laura Holt
Dear R :

Sorry for the off topic question, but does anyone know the reference for
the -2 Ln Lambda following a Chi Square distribution, please?

Possibly one of Bartlett's?

Thanks in advance!

Sincerely,
Laura Holt
mailto: [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


Re: [R] Sweave and complex numbers

2005-07-11 Thread Peter Dalgaard
Duncan Murdoch [EMAIL PROTECTED] writes:


  Question: why does \verb=c(1,1+1i)= get printed as
  \verb=c(1,1 + (0+1i))= ?
 
 The R parser only understands pure imaginary constants as complex 
 numbers.  It parses 1+1i as the sum of the real constant 1 and the 
 complex constant 0+1i.
 
 This isn't easy to work around.  Sweave could special case these, or the 
 parser or deparser could, but it looks messy to me.  I'd guess the best 
 solution would be if it happened in the parser (i.e. a real constant 
 plus an imaginary constant was folded into a complex constant), but 
 right now our parser doesn't do any sorts of optimizations like this. 
 It's not easy to add the first one, not least because the parser doesn't 
 know for sure that 1+1i really is a complex number:  you might have 
 redefined the meaning of +.

H... I'd say that the user would deserve what she gets in that
case. I suspect that the real issue is that R's tokenizer is
hand-written and not smart enough to recognize complex constants.
Figuring out whether there is a good reason for that or whether we
might actually use an automatic tokenizer like flex is somewhere on
the far end of my want-to-do list...

One thing that could be done quite immediately is to let the deparser
forget about a zero real part and drop the parentheses in that case.
After all 1i is perfectly legal R, and

 identical(1i,0+1i)
[1] TRUE

(Those parse/deparse asymmetries are maddening. *Next* time we
redesign R we should get this straightened out)

-- 
   O__   Peter Dalgaard Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark  Ph:  (+45) 35327918
~~ - ([EMAIL PROTECTED])  FAX: (+45) 35327907

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] plot(cox.zph()): customize xlab ylab

2005-07-11 Thread Adaikalavan Ramasamy
Duncan, your solution could be simplified using ann=FALSE in the plot

 fit - coxph( Surv(futime, fustat) ~ age + rx, ovarian)
 plot( cox.zph(fit), ann=F )
 title( xlab=My own label, ylab=A new label, main=A clever title)

Now, why did I not think of this before ?

Regards, Adai



On Mon, 2005-07-11 at 07:50 -0400, Duncan Murdoch wrote:
 Adaikalavan Ramasamy wrote:
  I am not sure if there is an easy way around this. An ugly hack is to
  make a copy the function survival:::plot.cox.zph and make your
  modified function. But there are others in the list who might know
  neater solutions.
 
 This hack is uglier (and might not work properly on some devices, but 
 you could try):
 
 plot(cox.zph(my.ph),var=1,col.lab=white)
 title(xlab=Year, ylab=Beta for blah blah blah)
 
 Duncan Murdoch
  
  Regards, Adai
  
  
  On Mon, 2005-07-11 at 11:10 +0100, Dan Bebber wrote:
  
 Hello,
 
 plot(cox.zph(my.ph),var=1,xlab=Year)
 
 gives the error:
 Error in plot.default(range(xx), yr, type = n, xlab
 = Time, ylab = ylab[i],  : formal argument xlab
 matched by multiple actual arguments
 
 How can I customize the xlab and ylab for plots of
 cox.zph?
 
 Thanks,
 Dan Bebber
 
 Department of Plant Sciences
 University of Oxford
 UK
 
 
 
 ___ 
 How much free photo storage do you get? Store your holiday
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 
  
  
  __
  R-help@stat.math.ethz.ch mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide! 
  http://www.R-project.org/posting-guide.html
 


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] plot(cox.zph()): customize xlab ylab

2005-07-11 Thread Duncan Murdoch
On 7/11/2005 9:10 AM, Adaikalavan Ramasamy wrote:
 Duncan, your solution could be simplified using ann=FALSE in the plot
 
  fit - coxph( Surv(futime, fustat) ~ age + rx, ovarian)
  plot( cox.zph(fit), ann=F )
  title( xlab=My own label, ylab=A new label, main=A clever title)
 
 Now, why did I not think of this before ?

That's nearly perfect.  I think the only problem is that plot.cox.zph 
can do multiple frames, and this will only allow annotation on the last 
one.  But if you put it in a loop and do them one at a time, this should 
be fine.

Duncan Murdoch

 
 Regards, Adai
 
 
 
 On Mon, 2005-07-11 at 07:50 -0400, Duncan Murdoch wrote:
 Adaikalavan Ramasamy wrote:
  I am not sure if there is an easy way around this. An ugly hack is to
  make a copy the function survival:::plot.cox.zph and make your
  modified function. But there are others in the list who might know
  neater solutions.
 
 This hack is uglier (and might not work properly on some devices, but 
 you could try):
 
 plot(cox.zph(my.ph),var=1,col.lab=white)
 title(xlab=Year, ylab=Beta for blah blah blah)
 
 Duncan Murdoch
  
  Regards, Adai
  
  
  On Mon, 2005-07-11 at 11:10 +0100, Dan Bebber wrote:
  
 Hello,
 
 plot(cox.zph(my.ph),var=1,xlab=Year)
 
 gives the error:
 Error in plot.default(range(xx), yr, type = n, xlab
 = Time, ylab = ylab[i],  : formal argument xlab
 matched by multiple actual arguments
 
 How can I customize the xlab and ylab for plots of
 cox.zph?
 
 Thanks,
 Dan Bebber
 
 Department of Plant Sciences
 University of Oxford
 UK
 
 

 ___ 
 How much free photo storage do you get? Store your holiday
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 
  
  
  __
  R-help@stat.math.ethz.ch mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide! 
  http://www.R-project.org/posting-guide.html
 


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] misc3d package

2005-07-11 Thread Mark Edmondson-Jones
Hi,

I am trying to install the misc3d package on a Windows (XP) installation of R 
2.0.1 using install.packages(misc3d) but with no success.  I have used this 
approach with other packages OK, but for misc3d I get the following output...

trying URL `http://cran.r-project.org/bin/windows/contrib/2.0/PACKAGES'
Content type `text/plain; charset=iso-8859-1' length 27996 bytes
opened URL
downloaded 27Kb

Warning message: 
No package misc3d on CRAN. in: download.packages(pkgs, destdir = tmpd, 
available = available,  

Any help would be much appreciated.

Regards,
Mark


This message has been checked for viruses but the contents of an attachment
may still contain software viruses, which could damage your computer system:
you are advised to perform your own checks. Email communications with the
University of Nottingham may be monitored as permitted by UK legislation.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Problems with corARMA

2005-07-11 Thread Constantinos Antoniou
Dear All,

I just came across the same error message (running a gnls{nlme}  
model). [R-2.1.1 on Mac OS 10.4]

For a reproducible example, please download the file:

http://mit.edu/costas/www/GR.txt

and run the following:


GR - read.table(GR.txt,header=T)
attach(GR)
myyear - year-1969

library(nlme)
gnls4a - gnls(fatalities/vehicles1000~a0*(vehicles1000/ 
population1000)^a1,start=c 
(a0=1,a1=1),data=GR,na.action=na.omit,correlation=corARMA(form=~myyear))


Thanks a lot,

Costas



On 10 Ιουν 2005, at 8:43 ΠΜ, Pamela McCaskie wrote:

 Dear all
 I am tryiing to fit the following lme with an ARMA correlation  
 structure:

 test - lme(fixed=fev1f~year, random=~1|id2, data=pheno2,
 correlation=corARMA(value=0.2, form=~year|id2), na.action=na.omit)

 But I get the following error message:

 Error in getGroupsFormula.default(correlation, asList = TRUE) :
 Form argument must be a formula

 I have used this same form argument with differerent correlation
 structures and it has worked fine. Does anyone know why it won't
 recognise ~year | id2 (or even ~ 1 | id2) as a formula?

 Any help would be great
 Pam

 -- 
 Pamela A McCaskie
 BSc(Mathematical Sciences)(Hons)

 Western Australian Institute for Medical Research
 University of Western Australia
 SCGH Campus
 Ground Floor, B Block
 QE-II Medical Centre
 Hospital Avenue, Nedlands
 Western Australia  6009
 AUSTRALIA
 Email:[EMAIL PROTECTED]
 Phone:+61-8-9346 1612
 Mob:  0417 926 607

 __
 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





-- 
Constantinos Antoniou, Ph.D.
Massachusetts Institute of Technology
Intelligent Transportation Systems Program
77 Massachusetts Ave., Rm. 1-249, Cambridge, MA 02139

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

Re: [R] misc3d package

2005-07-11 Thread Achim Zeileis
On Mon, 11 Jul 2005, Mark Edmondson-Jones wrote:

 Hi,

 I am trying to install the misc3d package on a Windows (XP) installation
 of R 2.0.1 using install.packages(misc3d) but with no success.  I have
 used this approach with other packages OK, but for misc3d I get the
 following output...

Try using a current version of R for which there is a pre-compiled
windows binary available or install the source package.
Best,
Z

 trying URL `http://cran.r-project.org/bin/windows/contrib/2.0/PACKAGES'
 Content type `text/plain; charset=iso-8859-1' length 27996 bytes
 opened URL
 downloaded 27Kb

 Warning message:
 No package misc3d on CRAN. in: download.packages(pkgs, destdir = tmpd, 
 available = available,

 Any help would be much appreciated.

 Regards,
 Mark


 This message has been checked for viruses but the contents of an attachment
 may still contain software viruses, which could damage your computer system:
 you are advised to perform your own checks. Email communications with the
 University of Nottingham may be monitored as permitted by UK legislation.

 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] misc3d package

2005-07-11 Thread Duncan Murdoch
On 7/11/2005 9:57 AM, Mark Edmondson-Jones wrote:
 Hi,
 
 I am trying to install the misc3d package on a Windows (XP) installation of R 
 2.0.1 using install.packages(misc3d) but with no success.  I have used this 
 approach with other packages OK, but for misc3d I get the following output...
 
 trying URL `http://cran.r-project.org/bin/windows/contrib/2.0/PACKAGES'
 Content type `text/plain; charset=iso-8859-1' length 27996 bytes
 opened URL
 downloaded 27Kb
 
 Warning message: 
 No package misc3d on CRAN. in: download.packages(pkgs, destdir = tmpd, 
 available = available,  
 
 Any help would be much appreciated.

The message should be read at face value.  There is no version of misc3d 
for R 2.0.1.  I think misc3d was released in June, well after 2.1.0 came 
out, and we don't do binary builds for old releases.

If you want to use that package, you'll need to upgrade (or try to 
compile it yourself).

Duncan Murdoch

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] misc3d package

2005-07-11 Thread Uwe Ligges
misc3d is a very new package.
CRAN's windows binary repository for R-2.0.x is no longer updated, hence 
does not contain the package. The corresponding ReadMe tells us: Last 
update: 19.04.2005.

Either upgrade to R-2.1.1 and try again, or compile the package from 
sources yourself.

Uwe Ligges



Mark Edmondson-Jones wrote:

 Hi,
 
 I am trying to install the misc3d package on a Windows (XP) installation of R 
 2.0.1 using install.packages(misc3d) but with no success.  I have used this 
 approach with other packages OK, but for misc3d I get the following output...
 
 trying URL `http://cran.r-project.org/bin/windows/contrib/2.0/PACKAGES'
 Content type `text/plain; charset=iso-8859-1' length 27996 bytes
 opened URL
 downloaded 27Kb
 
 Warning message: 
 No package misc3d on CRAN. in: download.packages(pkgs, destdir = tmpd, 
 available = available,  

 Any help would be much appreciated.
 
 Regards,
 Mark
 
 
 This message has been checked for viruses but the contents of an attachment
 may still contain software viruses, which could damage your computer system:
 you are advised to perform your own checks. Email communications with the
 University of Nottingham may be monitored as permitted by UK legislation.
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] plot(cox.zph()): customize xlab ylab

2005-07-11 Thread Thomas Lumley
On Mon, 11 Jul 2005, Adaikalavan Ramasamy wrote:

 I am not sure if there is an easy way around this. An ugly hack is to
 make a copy the function survival:::plot.cox.zph and make your
 modified function. But there are others in the list who might know
 neater solutions.

If you then send a patch to the package maintainer it stops being an ugly 
hack and turns into an example of collaborative open-source development :)

-thomas

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] class- vs. as()

2005-07-11 Thread Spencer Graves
  What do you want?  Consider the following:

  v - matrix(1:9, 3)
  class(v[1,1])
[1] integer
  class(as.vector(v))
[1] integer
  v2 - v
  dim(v2) - NULL
  class(v2)
[1] integer

  spencer graves

[EMAIL PROTECTED] wrote:

 
 
 
 Dear all,
 
 I would appreciate a lot, if someone could explain to me in a simple
 way, why the assignment class- is not always working and one has to
 take as() like in the example below.
 
 
(v - matrix(1:9, 3))
 
  [,1] [,2] [,3]
 [1,]147
 [2,]258
 [3,]369
 
class(v)
 
 [1] matrix
 
 
class(v) - integer
class(v)
 
 [1] matrix
 
v
 
  [,1] [,2] [,3]
 [1,]147
 [2,]258
 [3,]369
 
 
(vi - as(v, integer))
 
 [1] 1 2 3 4 5 6 7 8 9
 
class(vi)
 
 [1] integer
 
 
 With many thanks,
 
 Stefan
   [[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

-- 
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

[EMAIL PROTECTED]
www.pdf.com http://www.pdf.com
Tel:  408-938-4420
Fax: 408-280-7915

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] randomForest

2005-07-11 Thread Weiwei Shi
Thanks.
Many people pointed that out. (It was due to that I only knew lappy by
that time :).


On 7/11/05, Martin Maechler [EMAIL PROTECTED] wrote:
  Duncan == Duncan Murdoch [EMAIL PROTECTED]
  on Thu, 07 Jul 2005 15:44:38 -0400 writes:
 
 Duncan On 7/7/2005 3:38 PM, Weiwei Shi wrote:
  Hi there:
  I have a question on random foresst:
 
  recently i helped a friend with her random forest and i came with this 
 problem:
  her dataset has 6 classes and since the sample size is pretty small:
  264 and the class distr is like this (Diag is the response variable)
 
  sample.size - lapply(1:6, function(i) sum(Diag==i))
  sample.size
  [[1]]
  [1] 36
 
 
 
 and later you get problems because you didn't know that a *list*
 such as 'sample.size' should be made into a so called
 *atomic vector* {and there's a function  is.atomic(.) ! to test for it}
 and Duncan and others told you about unlist().
 
 Now there are two things I'd want to add:
 
 1) If you had used
 
   s.size - table(Diag)
 
you had used a faster and simpler expression with the same result.
Though in general (when there can be zero counts), to give the
same result, you'd need
 
   s.size - table(factor(Diag, levels = 1:6))
 
Still a bit preferable to the lapply(.) IMO
 
 
 2)  You should get into the habit of using
  sapply(.)   rather than  lapply(.).
 
 sapply() originally was exactly devised for the above task:
 and stands for ``[s]implified lapply''.
 
 It always returns an ``unlisted'' result when appropriate.
 
 Regards,
 Martin Maechler, ETH Zurich
 
 


-- 
Weiwei Shi, Ph.D

Did you always know?
No, I did not. But I believed...
---Matrix III

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] R on kubuntu

2005-07-11 Thread Henrik Andersson
However according to:

http://packages.ubuntu.com/hoary/math/r-base

you will have to live with R version 2.0.1 until the next version of 
Ubuntu is released. (If I understood the Ubuntu policy...)

Or is there some other repository providing more updated binaries for 
Ubuntu 5.04?

Cheers, Henrik Andersson

Lefebure Tristan wrote:
 No problem at all.
 If you allow universe packages, many binary R packages are available (from 
 the GNU/Linux Debian sid). 
 
 
 On Monday 11 July 2005 11:54, Constant Depièreux wrote:
 
Hello all,

I am planning to redeploy my workstation under KUBUNTU.

Does any body has any r experience installing/using r on this platform?

Best regards.
 
 


-- 
-
Henrik Andersson
Netherlands Institute of Ecology -
Centre for Estuarine and Marine Ecology
P.O. Box 140
4400 AC Yerseke
Phone: +31 113 577473
[EMAIL PROTECTED]
http://www.nioo.knaw.nl/ppages/handersson

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] plot(cox.zph()): customize xlab ylab

2005-07-11 Thread Dan Bebber
Dear all,

I've modified the plot.cox.zph function to allow
customized xlab and ylab (see below). Someone might
like to confirm that it works.

Thanks for all the assistance.

Dan
___

plot.cox.zph - function (x, resid = TRUE, se = TRUE,
df = 4, nsmo = 40, var, 
xlab=Time,ylab = paste(Beta(t) for,
dimnames(yy)[[2]]),...) 
{
xx - x$x
yy - x$y
d - nrow(yy)
df - max(df)
nvar - ncol(yy)
pred.x - seq(from = min(xx), to = max(xx), length
= nsmo)
temp - c(pred.x, xx)
lmat - ns(temp, df = df, intercept = TRUE)
pmat - lmat[1:nsmo, ]
xmat - lmat[-(1:nsmo), ]
qmat - qr(xmat)
if (se) {
bk - backsolve(qmat$qr[1:df, 1:df], diag(df))
xtx - bk %*% t(bk)
seval - d * ((pmat %*% xtx) * pmat) %*%
rep(1, df)
}
if (missing(var)) 
var - 1:nvar
else {
if (is.character(var)) 
var - match(var, dimnames(yy)[[2]])
if (any(is.na(var)) || max(var)  nvar ||
min(var)  
1) 
stop(Invalid variable requested)
}
if (x$transform == log) {
xx - exp(xx)
pred.x - exp(pred.x)
}
else if (x$transform != identity) {
xtime - as.numeric(dimnames(yy)[[1]])
apr1 - approx(xx, xtime, seq(min(xx),
max(xx), length = 17)[2 * 
(1:8)])
temp - signif(apr1$y, 2)
apr2 - approx(xtime, xx, temp)
xaxisval - apr2$y
xaxislab - rep(, 8)
for (i in 1:8) xaxislab[i] - format(temp[i])
}
for (i in var) {
y - yy[, i]
yhat - pmat %*% qr.coef(qmat, y)
if (resid) 
yr - range(yhat, y)
else yr - range(yhat)
if (se) {
temp - 2 * sqrt(x$var[i, i] * seval)
yup - yhat + temp
ylow - yhat - temp
yr - range(yr, yup, ylow)
}
if (x$transform == identity) 
plot(range(xx), yr, type = n, xlab =
xlab, ylab = ylab[i], 
...)
else if (x$transform == log) 
plot(range(xx), yr, type = n, xlab =
xlab, ylab = ylab[i], 
log = x, ...)
else {
plot(range(xx), yr, type = n, xlab =
xlab, ylab = ylab[i], 
axes = FALSE, ...)
axis(1, xaxisval, xaxislab)
axis(2)
box()
}
if (resid) 
points(xx, y)
lines(pred.x, yhat)
if (se) {
lines(pred.x, yup, lty = 2)
lines(pred.x, ylow, lty = 2)
}
}
}


--- Thomas Lumley [EMAIL PROTECTED] wrote:

 On Mon, 11 Jul 2005, Adaikalavan Ramasamy wrote:
 
  I am not sure if there is an easy way around this.
 An ugly hack is to
  make a copy the function survival:::plot.cox.zph
 and make your
  modified function. But there are others in the
 list who might know
  neater solutions.
 
 If you then send a patch to the package maintainer
 it stops being an ugly 
 hack and turns into an example of collaborative
 open-source development :)
 
   -thomas
 


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Problems with R on OS X

2005-07-11 Thread Rob J Goedman
Hi Heinz,

Can you send me your version of: '~/Library/Preferences/org.R- 
project.R.plist'?

Most Mac OS questions are posted/answered on R-SIG-Mac.

Thanks,
Rob


On Jul 11, 2005, at 2:42 AM, Heinz Schild wrote:

 I used R on OS X 10.3x quite some time with no serious problems.  
 Sometimes R stopped when I tried to execute a bigger program. After  
 updating to OS X to version 10.4 R worked but I still had the  
 problem with bigger programs. Therefore I re-installed R on top of  
 the existing R version. The installation finished properly but  
 suddenly R did not work. Then I reinstalled OS X 10.4 because I  
 thought some left over registry information from R may caused the  
 problem. R still crashed. I hoped the problem would be overcome  
 with the new R version 1.12. Unfortunately this is not the case.  
 The error report (see appendix) says that Thread 0 crashed.
 What can I do?
 Heinz Schild


 Errror Report.txt

 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting- 
 guide.html

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Off topic -2 Ln Lambda and Chi square

2005-07-11 Thread Liaw, Andy
If you meant the lambda as the likelihood ratio test statistic, the
asymptotic chi-squared distribution comes from the asymptotic normality of
the MLEs.  The proof is in a paper by Abraham Wald in 1943.  See Stuart 
Ord (Kendall's Advanced Statistics) for discussion (e.g., vol. 2, 5th
edition).

Andy

 From: Laura Holt
 
 Dear R :
 
 Sorry for the off topic question, but does anyone know the 
 reference for
 the -2 Ln Lambda following a Chi Square distribution, please?
 
 Possibly one of Bartlett's?
 
 Thanks in advance!
 
 Sincerely,
 Laura Holt
 mailto: [EMAIL PROTECTED]
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 
 


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] R on kubuntu

2005-07-11 Thread Douglas Bates
On 7/11/05, Henrik Andersson [EMAIL PROTECTED] wrote:
 However according to:
 
 http://packages.ubuntu.com/hoary/math/r-base
 
 you will have to live with R version 2.0.1 until the next version of
 Ubuntu is released. (If I understood the Ubuntu policy...)
 
 Or is there some other repository providing more updated binaries for
 Ubuntu 5.04?

It's possible to install the latest packages from the Debian
repository but awkward.  You would need to add a debian/testing or
debian/unstable repository to /etc/apt/sources.list and configure apt
to use the Ubuntu archives unless you explicitly require the testing
or unstable version of a package.  I think that trying to install such
a version will require you to update your entire compiler tool chain,
which may be more effort than you want to make.

I just recently installed Kubuntu on a machine and will, if I manage
to get the hardware working, do another today.  If that works I try
instaling  r-base-dev/unstabl and report back to the r-sig-debian
mailing list.

(Note that there is a special mailing list for Debian-related
questions and I am moving this thread to that list.)
 Cheers, Henrik Andersson
 
 Lefebure Tristan wrote:
  No problem at all.
  If you allow universe packages, many binary R packages are available (from
  the GNU/Linux Debian sid).
 
 
  On Monday 11 July 2005 11:54, Constant Depièreux wrote:
 
 Hello all,
 
 I am planning to redeploy my workstation under KUBUNTU.
 
 Does any body has any r experience installing/using r on this platform?
 
 Best regards.
 
 
 
 
 --
 -
 Henrik Andersson
 Netherlands Institute of Ecology -
 Centre for Estuarine and Marine Ecology
 P.O. Box 140
 4400 AC Yerseke
 Phone: +31 113 577473
 [EMAIL PROTECTED]
 http://www.nioo.knaw.nl/ppages/handersson
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] R on kubuntu

2005-07-11 Thread Lefebure Tristan
On Monday 11 July 2005 16:44, Henrik Andersson wrote:
 Or is there some other repository providing more updated binaries for
 Ubuntu 5.04?

I don't think, but I'm fine with R 2.0.1 ...
The next ubuntu release, Ubuntu 5.10 (The Breezy Badger), is for October 2005 
and will include R 2.1.1 (http://packages.ubuntu.com/breezy/math/r-base).


-- 

Tristan LEFEBURE
Laboratoire d'écologie des hydrosystèmes fluviaux (UMR 5023)
Université Lyon I - Campus de la Doua
Bat. Darwin C 69622 Villeurbanne - France

Phone: (33) (0)4 26 23 44 02
Fax: (33) (0)4 72 43 15 23

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Off topic -2 Ln Lambda and Chi square

2005-07-11 Thread Adelchi Azzalini
On Mon, 11 Jul 2005 11:18:55 -0400, Liaw, Andy wrote:

LA If you meant the lambda as the likelihood ratio test statistic, the
LA asymptotic chi-squared distribution comes from the asymptotic
LA normality of the MLEs.  The proof is in a paper by Abraham Wald in
LA 1943.  See Stuart  Ord (Kendall's Advanced Statistics) for
LA discussion (e.g., vol. 2, 5th edition).
LA 
LA Andy
LA 
LA  From: Laura Holt
LA  
LA  Dear R :
LA  
LA  Sorry for the off topic question, but does anyone know the 
LA  reference for
LA  the -2 Ln Lambda following a Chi Square distribution, please?
LA  

see

S. S Wilks ( 1938)
The large-sample distribution of the likelihood ratio for testing composite 
hypotheses
- Ann. Math. Stat, 9, 60-62

-- 
Adelchi Azzalini  [EMAIL PROTECTED]
Dipart.Scienze Statistiche, Università di Padova, Italia
tel. +39 049 8274147,  http://azzalini.stat.unipd.it/

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Isolating string containing only file name from complete path

2005-07-11 Thread Ken Termiso
Hi all,

What I'd like to do is to is to be able to extract a string corresponding to 
only the file name from a string containing the complete path, i.e. from the 
following path string:

/Users/ken/Desktop/test/runs/file1

I would like to end up with:

file1

This would be most ideally done in a platform-independent way.

Thanks in advance,
-Ken

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] plot(cox.zph()): customize xlab ylab

2005-07-11 Thread Adaikalavan Ramasamy
Dan, I think this works fine now. 

I think the 'ylab' argument in other plot function take only a single
value whereas in the case of plot.cox.zph, it needs a vector. It is
especially confusing when this function does plots everything on a
single page by default and all you see is the last plot.

Thus, can I suggest that you check that the length of the ylab is the
same as the number of terms in the coxph call. To do this you could put
something along the following anywhere after the line yy - x$y

n - length(coefficients(fit))
if( length(ylab) != n){
  warning( paste('ylab' is of length, length(ylab),
 and was expecting a length of, n, elements.) )
}

Regards, Adai



On Mon, 2005-07-11 at 16:04 +0100, Dan Bebber wrote:
 Dear all,
 
 I've modified the plot.cox.zph function to allow
 customized xlab and ylab (see below). Someone might
 like to confirm that it works.
 
 Thanks for all the assistance.
 
 Dan
 ___
 
 plot.cox.zph - function (x, resid = TRUE, se = TRUE,
 df = 4, nsmo = 40, var, 
 xlab=Time,ylab = paste(Beta(t) for,
 dimnames(yy)[[2]]),...) 
 {
 xx - x$x
 yy - x$y
 d - nrow(yy)
 df - max(df)
 nvar - ncol(yy)
 pred.x - seq(from = min(xx), to = max(xx), length
 = nsmo)
 temp - c(pred.x, xx)
 lmat - ns(temp, df = df, intercept = TRUE)
 pmat - lmat[1:nsmo, ]
 xmat - lmat[-(1:nsmo), ]
 qmat - qr(xmat)
 if (se) {
 bk - backsolve(qmat$qr[1:df, 1:df], diag(df))
 xtx - bk %*% t(bk)
 seval - d * ((pmat %*% xtx) * pmat) %*%
 rep(1, df)
 }
 if (missing(var)) 
 var - 1:nvar
 else {
 if (is.character(var)) 
 var - match(var, dimnames(yy)[[2]])
 if (any(is.na(var)) || max(var)  nvar ||
 min(var)  
 1) 
 stop(Invalid variable requested)
 }
 if (x$transform == log) {
 xx - exp(xx)
 pred.x - exp(pred.x)
 }
 else if (x$transform != identity) {
 xtime - as.numeric(dimnames(yy)[[1]])
 apr1 - approx(xx, xtime, seq(min(xx),
 max(xx), length = 17)[2 * 
 (1:8)])
 temp - signif(apr1$y, 2)
 apr2 - approx(xtime, xx, temp)
 xaxisval - apr2$y
 xaxislab - rep(, 8)
 for (i in 1:8) xaxislab[i] - format(temp[i])
 }
 for (i in var) {
 y - yy[, i]
 yhat - pmat %*% qr.coef(qmat, y)
 if (resid) 
 yr - range(yhat, y)
 else yr - range(yhat)
 if (se) {
 temp - 2 * sqrt(x$var[i, i] * seval)
 yup - yhat + temp
 ylow - yhat - temp
 yr - range(yr, yup, ylow)
 }
 if (x$transform == identity) 
 plot(range(xx), yr, type = n, xlab =
 xlab, ylab = ylab[i], 
 ...)
 else if (x$transform == log) 
 plot(range(xx), yr, type = n, xlab =
 xlab, ylab = ylab[i], 
 log = x, ...)
 else {
 plot(range(xx), yr, type = n, xlab =
 xlab, ylab = ylab[i], 
 axes = FALSE, ...)
 axis(1, xaxisval, xaxislab)
 axis(2)
 box()
 }
 if (resid) 
 points(xx, y)
 lines(pred.x, yhat)
 if (se) {
 lines(pred.x, yup, lty = 2)
 lines(pred.x, ylow, lty = 2)
 }
 }
 }
 
 
 --- Thomas Lumley [EMAIL PROTECTED] wrote:
 
  On Mon, 11 Jul 2005, Adaikalavan Ramasamy wrote:
  
   I am not sure if there is an easy way around this.
  An ugly hack is to
   make a copy the function survival:::plot.cox.zph
  and make your
   modified function. But there are others in the
  list who might know
   neater solutions.
  
  If you then send a patch to the package maintainer
  it stops being an ugly 
  hack and turns into an example of collaborative
  open-source development :)
  
  -thomas
  
  
 
 
 
   
   
   
 ___ 
 Yahoo! Messenger - NEW crystal clear PC to PC calling worldwide with 
 voicemail http://uk.messenger.yahoo.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


Re: [R] Isolating string containing only file name from complete path

2005-07-11 Thread Sebastian Luque
Ken Termiso [EMAIL PROTECTED] wrote:
 Hi all,

 What I'd like to do is to is to be able to extract a string
 corresponding to only the file name from a string containing the
 complete path, i.e. from the following path string:

 /Users/ken/Desktop/test/runs/file1

 I would like to end up with:

 file1


?basename

-- 
Sebastian P. Luque

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Isolating string containing only file name from complete path

2005-07-11 Thread Adaikalavan Ramasamy
x - /Users/ken/Desktop/test/runs/file1.txt
basename(x)
[1] file1.txt



On Mon, 2005-07-11 at 16:13 +, Ken Termiso wrote:
 Hi all,
 
 What I'd like to do is to is to be able to extract a string corresponding to 
 only the file name from a string containing the complete path, i.e. from the 
 following path string:
 
 /Users/ken/Desktop/test/runs/file1
 
 I would like to end up with:
 
 file1
 
 This would be most ideally done in a platform-independent way.
 
 Thanks in advance,
 -Ken
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] small first graph of par(3,2), other 5 are correct

2005-07-11 Thread Scot W McNary

Hi,

I'm trying to produce 6 graphs on a single page using code I've borrowed 
from an example by Paul Murrell:

(http://www.stat.auckland.ac.nz/~paul/RGraphics/custombase-xmastree.R).

It involves placing 6 horizontal barplots on one page and adding common 
labels.

The problem is the first graph in my figure (the one in the (1,1) 
position) is smaller than the other 5.  A toy example is included below.

When I compare the par() options set after producing each graph, all of 
the parameters that change look like what I would expect:

 names(after.g1[after.g1%in%after.g2==FALSE])
[1] fig mai mar mfg plt

 what.chgd-names(after.g1[after.g1%in%after.g2==FALSE])

 after.g1[what.chgd]
$fig
[1] 0.000 0.500 0.667 1.000

$mai
[1] 0.0309375 0.0618750 0.3093750 0.3093750

$mar
[1] 0.5 1.0 5.0 5.0

$mfg
[1] 1 1 3 2

$plt
[1] 0.02357143 0.88214285 0.01437097 0.85629032


 after.g2[what.chgd]
$fig
[1] 0.500 1.000 0.667 1.000

$mai
[1] 0.0309375 0.3093750 0.3093750 0.0618750

$mar
[1] 0.5 5.0 5.0 1.0

$mfg
[1] 1 2 3 2

$plt
[1] 0.11785715 0.97642857 0.01437097 0.85629032



One other aspect of this is that the problem does not occur if after I 
create all 6 graphs once, I rerun the code again, but omit only the 
par(mfrow=c(3,2)) statement.  If I don't reset par(), then the first 
graph is then the identical size of the other five.

Up to now I've gotten around this by simply running the code twice, 
omitting the par(mfrow=c(3,2)) the second time and printing the result. 
Now, however, I'm sending the graphs to a pdf file and must set the 
par(mfrow=c(3,2)) after opening the pdf device, so the problem shows up 
with each time.

Any help on fixing this problem so that all 6 graphs appear the same size 
on the figure would be most welcome.

Thanks,

Scot





Here is the toy example that shows the layout and plots in my figure:

par(mfrow=c(3,2))

groups-LETTERS[1:5]

### for graph 1 ###

# data
leftci - c(1:5)
rightci - c(2:6)

# left column graph so:
par(mar=c(0.5, 1, 5, 5))
# right column graph so:
#par(mar=c(0.5, 5, 5, 1))

plot.new()
title(main = Graph 1)
plot.window(xlim=c(0, 8), ylim=c(-1.5, 5.5))

ticks - seq(0, 8, 1)
y - 1:5  # how many spaces on y axis: 1 for each group
h - 0.2  # height?  a function of y?

segments(0, y, 8, y, lty=dotted)  # dotted segments on which bar lies 
(like a grid)
rect(leftci, y-h, rightci, y+h, col=dark grey)
#mtext(groups, at=y, adj=1, side=2, las=2, cex=.75)
par(cex.axis=1.0, mex=0.5)
axis(1, at=ticks, labels=abs(ticks), pos=0, )

box(inner, col=grey)

after.g1-show(par())

### for graph 2 ###

# data:
leftci - c(1:5)
rightci - c(2:6)

# left column graph so:
#par(mar=c(0.5, 1, 5, 5))
# right column graph so:
par(mar=c(0.5, 5, 5, 1))

plot.new()
title(main = Graph 2)

plot.window(xlim=c(0, 8), ylim=c(-1.5, 5.5))

ticks - seq(0, 8, 1)
y - 1:5  # how many spaces on y axis: 1 for each group
h - 0.2  # height?  a function of y?

segments(0, y, 8, y, lty=dotted)  # dotted segments on which bar lies 
(like a grid)
rect(leftci, y-h, rightci, y+h, col=dark grey)
mtext(groups, at=y, adj=.5, side=2, las=2, cex=.75, line = 5) # this line 
only for right column graphs
par(cex.axis=1.0, mex=0.5)
axis(1, at=ticks, labels=abs(ticks), pos=0, )

box(inner, col=grey)

# end graph 2

after.g2-show(par())

### graph 3###

# data
leftci - c(1:5)
rightci - c(2:6)

# left column graph so:
par(mar=c(0.5, 1, 5, 5))
# right column graph so:
#par(mar=c(0.5, 5, 5, 1))

plot.new()
title(main = Graph 3)
plot.window(xlim=c(0, 8), ylim=c(-1.5, 5.5))

ticks - seq(0, 8, 1)
y - 1:5  # how many spaces on y axis: 1 for each group
h - 0.2  # height?  a function of y?

segments(0, y, 8, y, lty=dotted)  # dotted segments on which bar lies 
(like a grid)
rect(leftci, y-h, rightci, y+h, col=dark grey)
#mtext(groups, at=y, adj=1, side=2, las=2, cex=.75)
par(cex.axis=1.0, mex=0.5)
axis(1, at=ticks, labels=abs(ticks), pos=0, )

box(inner, col=grey)

after.g3-show(par())

### for graph 4 ###

# data:
leftci - c(1:5)
rightci - c(2:6)

# left column graph so:
#par(mar=c(0.5, 1, 5, 5))
# right column graph so:
par(mar=c(0.5, 5, 5, 1))

plot.new()
title(main = Graph 4)

plot.window(xlim=c(0, 8), ylim=c(-1.5, 5.5))

ticks - seq(0, 8, 1)
y - 1:5  # how many spaces on y axis: 1 for each group
h - 0.2  # height?  a function of y?

segments(0, y, 8, y, lty=dotted)  # dotted segments on which bar lies 
(like a grid)
rect(leftci, y-h, rightci, y+h, col=dark grey)
mtext(groups, at=y, adj=.5, side=2, las=2, cex=.75, line = 5) # this line 
only for right column graphs
par(cex.axis=1.0, mex=0.5)
axis(1, at=ticks, labels=abs(ticks), pos=0, )

box(inner, col=grey)

after.g4-show(par())

# end graph 4

### graph 5###

# data
leftci - c(1:5)
rightci - c(2:6)

# left column graph so:
par(mar=c(0.5, 1, 5, 5))
# right column graph so:
#par(mar=c(0.5, 5, 5, 1))

plot.new()
title(main = Graph 5)
plot.window(xlim=c(0, 8), ylim=c(-1.5, 5.5))

ticks - seq(0, 8, 1)
y - 1:5  # how many spaces on y axis: 1 for each 

[R] exact values for p-values

2005-07-11 Thread S.O. Nyangoma
Hi there,
If I do an lm, I get p-vlues as

p-value:  2.2e-16

Suppose am interested in exact value such as 

p-value = 1.6e-16 (note = and not )

How do I go about it?

stephen

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Projection Pursuit

2005-07-11 Thread Oliver Lyttelton


Hello,

Just a quick question about ppr in library modreg.

I have looked at Ripley and Venables 2002 and it says that projection
pursuit works by projecting X in M carefully chosen directions

I want to know how it choses the directions? I presume it moves around the
high-dimensional space of unit vectors finding ones that separate the
response variables, but how.

I looked at the code, but can't get deeper in than the Fortran call to
Smart.

I tried to look up Friedman 87, but can't find an electronic reference.

Can anyone help, or point me in the direction of an online reference?

Thanks in advance,

Oliver Lyttelton

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] exact values for p-values - more information.

2005-07-11 Thread S.O. Nyangoma
 Hi there,
 If I do an lm, I get p-vlues as
 
 p-value:  2.2e-16
 
This is obtained from F =39540 with df1 = 1, df2 = 7025.

 Suppose am interested in exact value such as 
 
 p-value = 1.6e-16 (note = and not )
 
 How do I go about it?
 
 stephen

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] R on kubuntu

2005-07-11 Thread ernesto
You can compile the newest version...

EJ


Henrik Andersson wrote:

However according to:

http://packages.ubuntu.com/hoary/math/r-base

you will have to live with R version 2.0.1 until the next version of 
Ubuntu is released. (If I understood the Ubuntu policy...)

Or is there some other repository providing more updated binaries for 
Ubuntu 5.04?

Cheers, Henrik Andersson

Lefebure Tristan wrote:
  

No problem at all.
If you allow universe packages, many binary R packages are available (from 
the GNU/Linux Debian sid). 


On Monday 11 July 2005 11:54, Constant Depièreux wrote:



Hello all,

I am planning to redeploy my workstation under KUBUNTU.

Does any body has any r experience installing/using r on this platform?

Best regards.
  





  


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] exact values for p-values - more information.

2005-07-11 Thread Achim Zeileis
On Mon, 11 Jul 2005, S.O. Nyangoma wrote:

  Hi there,
  If I do an lm, I get p-vlues as

  p-value:  2.2e-16

 This is obtained from F =39540 with df1 = 1, df2 = 7025.

  Suppose am interested in exact value such as

  p-value = 1.6e-16 (note = and not )

  How do I go about it?

You can always extract the `exact' p-value from the summary.lm object or
you can compute it by hand via
  pf(39540, df1 = 1, df2 = 7025, lower.tail = FALSE)
For all practical purposes, the above means that the p-value is 0.
I guess you are on a 32-bit machine, then it also means that the p-value
is smaller than the Machine epsilon
  .Machine$double.eps

So if you want to report the p-value somewhere, I think R's output should
be more than precise enough. If you want to compute some other values that
depend on such a p-value, then it is probably wiser to compute on a log
scale, i.e. instead
  pf(70, df1 = 1, df2 = 7025, lower.tail = FALSE)
use
  pf(70, df1 = 1, df2 = 7025, lower.tail = FALSE, log.p = TRUE)

However, don't expect to be able to evaluate it at such extreme values
such as 39540.
Z

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] demo(scoping)

2005-07-11 Thread Prof Brian Ripley
Why have you sent a message about this, which no indication except the 
non-English word `entercount'?

Note that the message came from within try(), so it was intentional.
If you look at the source it says

try(ross$withdraw(500)) # no way..

It is helpful to learn how the code being demonstrated (here try) works.
If you don't understand it. please explain your difficulty in as much 
detail as you can.


On Mon, 4 Jul 2005, ronggui wrote:

 entercount an error with demo(scoping).

 demo(scoping)


demo(scoping)
 ~~~
 ___snip_


 ross$balance()
 Your balance is 120


 try(ross$withdraw(500))
 Error in ross$withdraw(500) : You don't have that much money!


 version
 _
 platform i486-pc-linux-gnu
 arch i486
 os   linux-gnu
 system   i486, linux-gnu
 status   beta
 major2
 minor1.1
 year 2005
 month06
 day  13
 language R


 -- 
 Department of Sociology
 Fudan University,Shanghai
 Blog:http://sociology.yculblog.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


-- 
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


Re: [R] timezone problems

2005-07-11 Thread Prof Brian Ripley
For the record, this was covered by an answer to your bug report.

The problem is your OS which mishandles a timezone of `GMT', so 
Sys.getenv(TZ) as `GMT' is not actually setting your OS to GMT.
Hence NA is the correct answer.

I know no way to set Windows to GMT as distinct from the timezone of 
London (with summer time).

On Thu, 7 Jul 2005, Martin Keller-Ressel wrote:

 Thank you Don for your hints. I have checked my environment vairable TZ
 again. But everything is set correctly. I think the problem is with
 Sys.timezone(). Maybe it is a conflict between how my system formats the
 time/date and what Sys.timezone() expects.
 This is what I get on my system:

 Sys.getenv(TZ)
TZ
 GMT
 Sys.time()
 [1] 2005-07-07 07:32:39 GMT

 ## everything fine so far

 Sys.timezone()
 [1] NA

 ## This is what Sys.timezone looks like:
 Sys.timezone
 function ()
 {
 z - as.POSIXlt(Sys.time())
 attr(z, tzone)[2 + z$isdst]
 }
 environment: namespace:base

 z - as.POSIXlt(Sys.time())
 attributes(z)
 $names
 [1] sec   min   hour  mday  mon   year  wday  yday  isdst

 $class
 [1] POSIXt  POSIXlt

 $tzone
 [1] GMT

 attr(z,tzone)
 [1] GMT
 z$isdst
 [1] 0
 attr(z,tzone)[2]
 [1] NA

 I dont understand why Sys.timezone doesn't use attr(z,tzone) but tries
 to read its (2+z$isdst)-th element.
 Of course it would be easy to write a workaround, but I wonder why nobody
 else is having this problem.

 best regards,

 Martin Keller-Ressel



 On Wed, 06 Jul 2005 14:45:25 -, Don MacQueen [EMAIL PROTECTED] wrote:

 How did you set the TZ system variable?
 If you did not use Sys.putenv(), try using it instead.
 Otherwise, I think you have to ask the package maintainer.

 You may be misleading yourself by using Sys.time() to test whether TZ is
 set.
 What does Sys.getenv() tell you?

 I get a timezone code from Sys.time() even when TZ is not defined (see
 example below).
 (but I do have a different OS)

  Sys.timezone()
 [1] 
  Sys.time()
 [1] 2005-07-06 07:34:15 PDT
  Sys.getenv('TZ')
 TZ
 
  Sys.putenv(TZ='US/Pacific')
  Sys.timezone()
 [1] US/Pacific
  Sys.getenv('TZ')
TZ
 US/Pacific
  Sys.time()
 [1] 2005-07-06 07:34:38 PDT

  Sys.putenv(TZ='GMT')
  Sys.time()
 [1] 2005-07-06 14:35:45 GMT

  version
   _   platform powerpc-apple-darwin7.9.0
 arch powerpc os   darwin7.9.0
 system   powerpc, darwin7.9.0status
 major2   minor1.1
 year 2005month06
 day  20  language R
 At 9:55 AM + 7/5/05, Martin Keller-Ressel wrote:
 Hi,

 Im using R 2.1.1 and running Code that previously worked (on R 2.1.0 I
 believe) using the 'timeDate' function from the fCalendar package. The
 code now throws an error:

 Error in if (Sys.timezone() != GMT) warning(Set timezone to GMT!)

 However I have read the documentation of the fCalendar package and I
 have set my system variable TZ to GMT.
 I tracked the error down to the function Sys.timezone() which returns
 NA in spite of what Sys.time() returns.

  Sys.timezone()
 [1] NA

  Sys.time()
 [1] 2005-07-05 08:41:53 GMT

 My version:

  version
   _
 platform i386-pc-mingw32
 arch i386
 os   mingw32
 system   i386, mingw32
 status
 major2
 minor1.1
 year 2005
 month06
 day  20
 language R

 Any help is appreciated,

 Martin Keller-Ressel


 ---
 Martin Keller-Ressel
 Research Unit of Financial and Actuarial Mathematics
 TU Vienna

 __
 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





 -- 
 Martin Keller-Ressel
 Research Unit of Financial and Actuarial Mathematics
 TU Vienna

 __
 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


-- 
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


Re: [R] Weighted nls

2005-07-11 Thread Liaw, Andy
Please do your homework, as the posting guide asks you to.  

1.  From the Example section of ?nls:

## weighted nonlinear regression
Treated - Puromycin[Puromycin$state == treated, ]
weighted.MM - function(resp, conc, Vm, K)
{
## Purpose: exactly as white book p.451 -- RHS for nls()
##  Weighted version of Michaelis-Menten model
## 
## Arguments: 'y', 'x' and the two parameters (see book)
## 
## Author: Martin Maechler, Date: 23 Mar 2001, 18:48

pred - (Vm * conc)/(K + conc)
(resp - pred) / sqrt(pred)
}

Pur.wt - nls( ~ weighted.MM(rate, conc, Vm, K), data = Treated,
  start = list(Vm = 200, K = 0.1),
  trace = TRUE)

2.  This has been covered several times on this list.  Search the archive.

Andy

 From: Rundle, Daniel
 
 Dear R Community,
  
 I am attempting to perform a weighted non-linear least 
 squares fit.  It has already been noted that the weights 
 option is not yet implemented for the nls function, but no 
 one seems to offer any suggestions for getting around this 
 problem.  I am still curious if a) anyone has code they have 
 written which includes a weight options for nls, or b) if 
 there is another model which would allow least squares 
 fitting to a non-linear function.
  
 Thank you for reading this,
  
 Dan
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 
 


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] exact values for p-values - more information.

2005-07-11 Thread Spencer Graves
  I just checked:

  pf(39540, 1, 7025, lower.tail=FALSE, log.p=TRUE)
[1] -Inf

  This is not correct.  With 7025 denominator degrees of freedom, we 
might use the chi-square approximation to the F distribution:

  pchisq(39540, 1, lower.tail=FALSE, log.p=TRUE)
[1] -19775.52

  In sum, my best approximation to  pf(39540, 1, 7025, 
lower.tail=FALSE, log.p=TRUE), given only a minute to work on this, is 
exp(pchisq(39540, 1, lower.tail=FALSE, log.p=TRUE)) = exp(-19775.52).

  I'm confident that many violations of assumptions would likely be 
more important than the differences between p-value:  2.2e-16 and 
these other two answers.  However, I have also used numbers like 
exp(-19775.52) to guestimate relative degrees of plausibility for 
different alternatives.  That doesn't mean they are right, only the best 
I can get with the available resources.

  spencer graves

Achim Zeileis wrote:

 On Mon, 11 Jul 2005, S.O. Nyangoma wrote:
 
 
 Hi there,
 If I do an lm, I get p-vlues as

 p-value:  2.2e-16

This is obtained from F =39540 with df1 = 1, df2 = 7025.

 Suppose am interested in exact value such as

 p-value = 1.6e-16 (note = and not )

 How do I go about it?
 
 
 You can always extract the `exact' p-value from the summary.lm object or
 you can compute it by hand via
   pf(39540, df1 = 1, df2 = 7025, lower.tail = FALSE)
 For all practical purposes, the above means that the p-value is 0.
 I guess you are on a 32-bit machine, then it also means that the p-value
 is smaller than the Machine epsilon
   .Machine$double.eps
 
 So if you want to report the p-value somewhere, I think R's output should
 be more than precise enough. If you want to compute some other values that
 depend on such a p-value, then it is probably wiser to compute on a log
 scale, i.e. instead
   pf(70, df1 = 1, df2 = 7025, lower.tail = FALSE)
 use
   pf(70, df1 = 1, df2 = 7025, lower.tail = FALSE, log.p = TRUE)
 
 However, don't expect to be able to evaluate it at such extreme values
 such as 39540.
 Z
 
 __
 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

-- 
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

[EMAIL PROTECTED]
www.pdf.com http://www.pdf.com
Tel:  408-938-4420
Fax: 408-280-7915

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Projection Pursuit

2005-07-11 Thread Liaw, Andy
Google can be helpful.  The 10th hit I got from projection pursuit
regression is http://www.scs.gmu.edu/~jgentle/csi991/03f/ppreg1024.rtf,
which gives some rough outline of the algorithm.  If you want more detail,
Ripley (1996) PRNN would suffice, I believe.

Andy

 From: Oliver Lyttelton
 
 Hello,
 
 Just a quick question about ppr in library modreg.
 
 I have looked at Ripley and Venables 2002 and it says that projection
 pursuit works by projecting X in M carefully chosen directions
 
 I want to know how it choses the directions? I presume it 
 moves around the
 high-dimensional space of unit vectors finding ones that separate the
 response variables, but how.
 
 I looked at the code, but can't get deeper in than the Fortran call to
 Smart.
 
 I tried to look up Friedman 87, but can't find an electronic 
 reference.
 
 Can anyone help, or point me in the direction of an online reference?
 
 Thanks in advance,
 
 Oliver Lyttelton
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 
 


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Weighted nls

2005-07-11 Thread Peter Dalgaard
Rundle, Daniel [EMAIL PROTECTED] writes:

 Dear R Community,
  
 I am attempting to perform a weighted non-linear least squares fit.
 It has already been noted that the weights option is not yet
 implemented for the nls function, but no one seems to offer any
 suggestions for getting around this problem. I am still curious if
 a) anyone has code they have written which includes a weight options
 for nls, or b) if there is another model which would allow least
 squares fitting to a non-linear function.
  

gnls in the nlme package. Or follow the example on help(nls) immediately
after the  

## weighted nonlinear regression

(!)

-- 
   O__   Peter Dalgaard Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark  Ph:  (+45) 35327918
~~ - ([EMAIL PROTECTED])  FAX: (+45) 35327907

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] small first graph of par(3,2), other 5 are correct

2005-07-11 Thread Brahm, David
Scot,

Here is your toy example in more condensed form:

x11()
par(mar=c(0.5, 1, 5, 5))
par(mfrow=c(3,2))
plot(1:10)
par(mex)   # mex=1.0 here
par(cex.axis=1.0, mex=0.5)   # Now you change it
for (i in 1:5) plot(1:10)

When you build your first plot (effectively at the plot.new() command),
your mex parameter is 1, and that leaves lots of space in the margins.
Then you set mex=0.5, which holds for all subsequent plots, and they get
less margin space.

-- David Brahm ([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


Re: [R] exact values for p-values - more information.

2005-07-11 Thread Adaikalavan Ramasamy
Compare the following

   t.test( 1:100, 101:200 )$p.value 
   t.test( 1:100, 101:200 )

In the latter, the print method truncates to 2.2e-16. 
You can go as far as (depending on your machine)

   .Machine$double.xmin
   [1] 2.225074e-308

before it becomes indistinguishable from zero.

But there are good reasons to truncate it at 2.2e-16 such as the
difficulty in trying to accurately estimate the extreme tail
probabilities.

Regards, Adai



On Mon, 2005-07-11 at 18:52 +0200, S.O. Nyangoma wrote:
  Hi there,
  If I do an lm, I get p-vlues as
  
  p-value:  2.2e-16
  
 This is obtained from F =39540 with df1 = 1, df2 = 7025.
 
  Suppose am interested in exact value such as 
  
  p-value = 1.6e-16 (note = and not )
  
  How do I go about it?
  
  stephen
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] indexing into and modifying dendrograms

2005-07-11 Thread Jenny Bryan
I would like to be able to exert certain types of control over the
plotting of dendrograms (representing hierarchical clusterings) that I
think is best achieved by modifying the dendrogram object
prior to plotting.  I am using the dendrogram class and associated
methods.

Define the cluster number of each cluster formed as the corresponding
row of the merge object.  So, if you are clustering m objects, the
cluster numbers range from 1 to m-1, with cluster m-1 containing all m
objects, by definition.  I basically want a way to index into an
object of class dendrogram using the above definition of cluster
number and/or to act on a dendrogram, where I specify the target node
using cluster number.

The first application would be to 'flip' the two elements in target
node of the dendrogram (made clear in the small example below). (The
setting is genomics and I have applications where I want to man-handle
my dendrograms to make certain features of the clustering more obvious
to the naked eye.)  I could imagine other, related actions that would
be useful in decorating dendrograms.

I think I need a function that takes a dendrogram and cluster
number(s) as input and returns the relevant part(s) of the dendrogram
object -- but in a form that makes it easy to then, say, set certain
attributes (perhaps recursively) for the target nodes (and perhaps
those contained in it).  I'm including a small example below that
hopefully illustrates this (it looks long, but it's mostly comments!).

Any help would be appreciated.

Jenny Bryan

## get ready for step-by-step figures
par(mfrow = c(2,2))

## get 5 objects, with 2-dimensional features
pts - rbind(c(2,1.6),
  c(1.8,2.4),
  c(2.1, 2.7),
  c(5,2.6),
  c(4.7,3.1))
plot(pts, xlim = c(0,6), ylim = c(0,4),type = n,
  xlab = Feature 1, ylab = Feature 2)
points(pts,pch = as.character(1:5))

## build a hierarhical tree, store as a dendrogram
aggTree - hclust(dist(pts), method = single)
(dend1 - JB.as.dendrogram.hclust(aggTree))
## NOTE: only thing I added to official version of
## as.dendrogram.hclust:
## each node has an attribute cNum, which gives
## the merge step at which it was formed,
## i.e. gives the row of the merge object which
## describes the formation of that node
## one new line near end of nMerge loop:
## ***
## *** 51,56 
## --- 51,60 
##   attr(z[[x[2]]], midpoint))/2
##  }
##   attr(zk, height) - oHgt[k]
## +
## + ## JB added July 6 2005
## + attr(zk, cNum) - k
## +
##  z[[k - as.character(k)]] - zk
##   }
##   z - z[[k]]
attributes(dend1)
attributes(dend1[[1]])
## here's a table relating dend1 and the cNum attribute
## dend1   cNum
## -
## dend14
## dend1[[1]]   2
## dend1[[2]]   3
## dend1[[2]][[1]]  not set
## dend1[[2]][[1]]  1

## use cNum attribute in edgetext
## following example in dendrogram documentation
## would really rather associate with the node than the edge
## but current plotting function has no notion of nodetext
addE - function(n) {
   if(!is.leaf(n)) {
 attr(n, edgePar) - list(p.col=plum)
 attr(n, edgetext) - attr(n,cNum)
   }
   n
}
dend2 - dendrapply(dend1, addE)
## overlays the cNum (cluster number) attribute on dendrogram
plot(dend2, main = dend2)
## why does no plum polygon appear around the '4' for the root
## edge?

## swap order of clusters 2 and 3,
## i.e. 'flip' cluster 4
dend3 - dend2
dend3[[1]] - dend2[[2]]
dend3[[2]] - dend2[[1]]
plot(dend3, main = dend3)
## wish I could achieve with 'dend3 - flip(dend2, cNum = 4)

## swap order of cluster 1 and object 1,
## i.e. 'flip' cluster 3
dend4 - dend2
dend4[[2]][[1]] - dend2[[2]][[2]]
dend4[[2]][[2]] - dend2[[2]][[1]]
plot(dend4, main = dend4)
## wish I could achieve with 'dend4 - flip(dend2, cNum = 3)

## finally, it's clear that the midpoint attribute would also
## need to be modified by 'flip'

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] building packages on Windows

2005-07-11 Thread Sundar Dorai-Raj
Hi, all,

I just recently upgraded my computer though I'm using the same OS (XP).
But now I'm having difficulty building packages and I cannot seem to
solve the problem. I'm using R-2.1.1pat on Windows XP.

Here is what I tried:

D:\Users\sundard\slib\sundar\RR CMD CHECK sundar
* checking for working latex ... OK
* using log directory 'D:/Users/sundard/slib/sundar/R/sundar.Rcheck'
* using R version 2.1.1, 2005-06-21
* checking for file 'sundar/DESCRIPTION' ... OK
* this is package 'sundar' version '1.1'
* checking if this is a source package ... OK

installing R.css in D:/Users/sundard/slib/sundar/R/sundar.Rcheck


-- Making package sundar 
   adding build stamp to DESCRIPTION
   installing NAMESPACE file and metadata
Error in file(file, r) : unable to open connection
In addition: Warning message:
cannot open file
'D:/Users/sundard/slib/sundar/R/sundar.Rcheck/sundar/NAMESPACE'

Execution halted
make[2]: *** [nmspace] Error 1
make[1]: *** [all] Error 2
make: *** [pkg-sundar] Error 2
*** Installation of sundar failed ***

Removing 'D:/Users/sundard/slib/sundar/R/sundar.Rcheck/sundar'
  ERROR
Installation failed.

I've also tried to remove the NAMESPACE which then passes `check' but
fails on `build' with the following:

* checking for file 'sundar/DESCRIPTION' ... OK
* preparing 'sundar':
* checking DESCRIPTION meta-information ... OK
* removing junk files
Error: cannot open file 'sundar/DESCRIPTION' for reading

If this is relevant, here is my PATH:
C:\WINDOWS\system32;
C:\WINDOWS;
.;
D:\R\rw2011pat\bin;
D:\R\tools\bin;
D:\Perl\bin;
D:\Tcl\bin;
D:\TeXLive\bin\win32;
D:\Mingw\bin;
C:\Program Files\Insightful\splus62\cmd;
C:\Program Files\HTML Help Workshop

Here, the D:\R\tools\bin directory contains the utilities for building R
from source, downloaded today. Am I missing something obvious?

Thanks,

--sundar

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] building packages on Windows

2005-07-11 Thread Duncan Murdoch
On 7/11/2005 3:21 PM, Sundar Dorai-Raj wrote:
 Hi, all,
 
 I just recently upgraded my computer though I'm using the same OS (XP).
 But now I'm having difficulty building packages and I cannot seem to
 solve the problem. I'm using R-2.1.1pat on Windows XP.
 
 Here is what I tried:
 
 D:\Users\sundard\slib\sundar\RR CMD CHECK sundar
 * checking for working latex ... OK
 * using log directory 'D:/Users/sundard/slib/sundar/R/sundar.Rcheck'
 * using R version 2.1.1, 2005-06-21
 * checking for file 'sundar/DESCRIPTION' ... OK
 * this is package 'sundar' version '1.1'
 * checking if this is a source package ... OK
 
 installing R.css in D:/Users/sundard/slib/sundar/R/sundar.Rcheck
 
 
 -- Making package sundar 
adding build stamp to DESCRIPTION
installing NAMESPACE file and metadata
 Error in file(file, r) : unable to open connection
 In addition: Warning message:
 cannot open file
 'D:/Users/sundard/slib/sundar/R/sundar.Rcheck/sundar/NAMESPACE'
 
 Execution halted
 make[2]: *** [nmspace] Error 1
 make[1]: *** [all] Error 2
 make: *** [pkg-sundar] Error 2
 *** Installation of sundar failed ***
 
 Removing 'D:/Users/sundard/slib/sundar/R/sundar.Rcheck/sundar'
   ERROR
 Installation failed.
 
 I've also tried to remove the NAMESPACE which then passes `check' but
 fails on `build' with the following:
 
 * checking for file 'sundar/DESCRIPTION' ... OK
 * preparing 'sundar':
 * checking DESCRIPTION meta-information ... OK
 * removing junk files
 Error: cannot open file 'sundar/DESCRIPTION' for reading
 
 If this is relevant, here is my PATH:
 C:\WINDOWS\system32;
 C:\WINDOWS;
 .;
 D:\R\rw2011pat\bin;
 D:\R\tools\bin;
 D:\Perl\bin;
 D:\Tcl\bin;
 D:\TeXLive\bin\win32;
 D:\Mingw\bin;
 C:\Program Files\Insightful\splus62\cmd;
 C:\Program Files\HTML Help Workshop
 
 Here, the D:\R\tools\bin directory contains the utilities for building R
 from source, downloaded today. Am I missing something obvious?

Yes, you need to put the Windows directories later in your path.  There 
are some programs in the toolset which have like-named commands in 
Windows; you need to use ours, not theirs.

There are other differences between your path and the one recommended in 
the Admin and Installation manual; they may also be causing problems.

Duncan Murdoch

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] building packages on Windows

2005-07-11 Thread Sundar Dorai-Raj


Duncan Murdoch wrote:
 On 7/11/2005 3:21 PM, Sundar Dorai-Raj wrote:
 
 Hi, all,

 I just recently upgraded my computer though I'm using the same OS (XP).
 But now I'm having difficulty building packages and I cannot seem to
 solve the problem. I'm using R-2.1.1pat on Windows XP.

 Here is what I tried:

 D:\Users\sundard\slib\sundar\RR CMD CHECK sundar
 * checking for working latex ... OK
 * using log directory 'D:/Users/sundard/slib/sundar/R/sundar.Rcheck'
 * using R version 2.1.1, 2005-06-21
 * checking for file 'sundar/DESCRIPTION' ... OK
 * this is package 'sundar' version '1.1'
 * checking if this is a source package ... OK

 installing R.css in D:/Users/sundard/slib/sundar/R/sundar.Rcheck


 -- Making package sundar 
adding build stamp to DESCRIPTION
installing NAMESPACE file and metadata
 Error in file(file, r) : unable to open connection
 In addition: Warning message:
 cannot open file
 'D:/Users/sundard/slib/sundar/R/sundar.Rcheck/sundar/NAMESPACE'

 Execution halted
 make[2]: *** [nmspace] Error 1
 make[1]: *** [all] Error 2
 make: *** [pkg-sundar] Error 2
 *** Installation of sundar failed ***

 Removing 'D:/Users/sundard/slib/sundar/R/sundar.Rcheck/sundar'
   ERROR
 Installation failed.

 I've also tried to remove the NAMESPACE which then passes `check' but
 fails on `build' with the following:

 * checking for file 'sundar/DESCRIPTION' ... OK
 * preparing 'sundar':
 * checking DESCRIPTION meta-information ... OK
 * removing junk files
 Error: cannot open file 'sundar/DESCRIPTION' for reading

 If this is relevant, here is my PATH:
 C:\WINDOWS\system32;
 C:\WINDOWS;
 .;
 D:\R\rw2011pat\bin;
 D:\R\tools\bin;
 D:\Perl\bin;
 D:\Tcl\bin;
 D:\TeXLive\bin\win32;
 D:\Mingw\bin;
 C:\Program Files\Insightful\splus62\cmd;
 C:\Program Files\HTML Help Workshop

 Here, the D:\R\tools\bin directory contains the utilities for building R
 from source, downloaded today. Am I missing something obvious?
 
 
 Yes, you need to put the Windows directories later in your path.  There 
 are some programs in the toolset which have like-named commands in 
 Windows; you need to use ours, not theirs.
 
 There are other differences between your path and the one recommended in 
 the Admin and Installation manual; they may also be causing problems.
 
 Duncan Murdoch


Hi, Duncan,

Thanks for the reply. I moved things around so that my PATH now looks like:

.;
D:\R\tools\bin;
D:\Perl\bin;
D:\Mingw\bin;
D:\TeXLive\bin\win32;
C:\Program Files\HTML Help Workshop;
D:\R\rw2011pat\bin;
D:\Tcl\bin;
C:\Program Files\Insightful\splus62\cmd;
C:\WINDOWS;
C:\WINDOWS\system32

However, the problem was actually to do with permissions. Apparently 
when I when I do the following:

chmod 666 DESCRIPTION NAMESPACE

The permissions must have been changed when upgrading my computer since 
I didn't have any problem on my previous machine. However, now another 
problem has arisen:

D:\Users\sundard\slib\sundar\RR CMD install sundar


-- Making package sundar 
   adding build stamp to DESCRIPTION
   installing NAMESPACE file and metadata

output snipped

Compiling d:\Users\sundard\slib\sundar\R\sundar\chm\sundar.chm

HHC5003: Error: Compilation failed while compiling logo.jpg.
HHC5003: Error: Compilation failed while compiling Rchm.css.


The following files were not compiled:
logo.jpg
Rchm.css
   adding MD5 sums

* DONE (sundar)

This still installs fine except the chmhelp file is not created as it 
should be. I've deleted the chm subdirectory and tried to reinstall and 
I get the same error. This seems not to be an R problem, but possibly 
some configuration I'm not aware of, but I'd be grateful for any insight.

Thanks,

--sundar

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] building packages on Windows

2005-07-11 Thread Gabor Grothendieck
On 7/11/05, Duncan Murdoch [EMAIL PROTECTED] wrote:
 On 7/11/2005 3:21 PM, Sundar Dorai-Raj wrote:
  Hi, all,
 
  I just recently upgraded my computer though I'm using the same OS (XP).
  But now I'm having difficulty building packages and I cannot seem to
  solve the problem. I'm using R-2.1.1pat on Windows XP.
 
  Here is what I tried:
 
  D:\Users\sundard\slib\sundar\RR CMD CHECK sundar
  * checking for working latex ... OK
  * using log directory 'D:/Users/sundard/slib/sundar/R/sundar.Rcheck'
  * using R version 2.1.1, 2005-06-21
  * checking for file 'sundar/DESCRIPTION' ... OK
  * this is package 'sundar' version '1.1'
  * checking if this is a source package ... OK
 
  installing R.css in D:/Users/sundard/slib/sundar/R/sundar.Rcheck
 
 
  -- Making package sundar 
 adding build stamp to DESCRIPTION
 installing NAMESPACE file and metadata
  Error in file(file, r) : unable to open connection
  In addition: Warning message:
  cannot open file
  'D:/Users/sundard/slib/sundar/R/sundar.Rcheck/sundar/NAMESPACE'
 
  Execution halted
  make[2]: *** [nmspace] Error 1
  make[1]: *** [all] Error 2
  make: *** [pkg-sundar] Error 2
  *** Installation of sundar failed ***
 
  Removing 'D:/Users/sundard/slib/sundar/R/sundar.Rcheck/sundar'
ERROR
  Installation failed.
 
  I've also tried to remove the NAMESPACE which then passes `check' but
  fails on `build' with the following:
 
  * checking for file 'sundar/DESCRIPTION' ... OK
  * preparing 'sundar':
  * checking DESCRIPTION meta-information ... OK
  * removing junk files
  Error: cannot open file 'sundar/DESCRIPTION' for reading
 
  If this is relevant, here is my PATH:
  C:\WINDOWS\system32;
  C:\WINDOWS;
  .;
  D:\R\rw2011pat\bin;
  D:\R\tools\bin;
  D:\Perl\bin;
  D:\Tcl\bin;
  D:\TeXLive\bin\win32;
  D:\Mingw\bin;
  C:\Program Files\Insightful\splus62\cmd;
  C:\Program Files\HTML Help Workshop
 
  Here, the D:\R\tools\bin directory contains the utilities for building R
  from source, downloaded today. Am I missing something obvious?
 
 Yes, you need to put the Windows directories later in your path.  There
 are some programs in the toolset which have like-named commands in
 Windows; you need to use ours, not theirs.
 
 There are other differences between your path and the one recommended in
 the Admin and Installation manual; they may also be causing problems.

In 
http://cran.r-project.org/contrib/extra/batchfiles/ 

there is a batch file called Rfind.bat that searches your system and lists 
the paths of the various tools used by R.  It does not actually change 
any environment variables nor does it change any other aspect of your 
system -- its display only so its generally safe to run.

This might help you in configuration.

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] indexing into and modifying dendrograms

2005-07-11 Thread Adaikalavan Ramasamy
Probably not answering your questions here, but have you considered the
functions prune.tree and snip.tree from package tree or prune.rpart and
snip.rpart from the package rpart ?



On Mon, 2005-07-11 at 11:48 -0700, Jenny Bryan wrote:
 I would like to be able to exert certain types of control over the
 plotting of dendrograms (representing hierarchical clusterings) that I
 think is best achieved by modifying the dendrogram object
 prior to plotting.  I am using the dendrogram class and associated
 methods.
 
 Define the cluster number of each cluster formed as the corresponding
 row of the merge object.  So, if you are clustering m objects, the
 cluster numbers range from 1 to m-1, with cluster m-1 containing all m
 objects, by definition.  I basically want a way to index into an
 object of class dendrogram using the above definition of cluster
 number and/or to act on a dendrogram, where I specify the target node
 using cluster number.
 
 The first application would be to 'flip' the two elements in target
 node of the dendrogram (made clear in the small example below). (The
 setting is genomics and I have applications where I want to man-handle
 my dendrograms to make certain features of the clustering more obvious
 to the naked eye.)  I could imagine other, related actions that would
 be useful in decorating dendrograms.
 
 I think I need a function that takes a dendrogram and cluster
 number(s) as input and returns the relevant part(s) of the dendrogram
 object -- but in a form that makes it easy to then, say, set certain
 attributes (perhaps recursively) for the target nodes (and perhaps
 those contained in it).  I'm including a small example below that
 hopefully illustrates this (it looks long, but it's mostly comments!).
 
 Any help would be appreciated.
 
 Jenny Bryan
 
 ## get ready for step-by-step figures
 par(mfrow = c(2,2))
 
 ## get 5 objects, with 2-dimensional features
 pts - rbind(c(2,1.6),
   c(1.8,2.4),
   c(2.1, 2.7),
   c(5,2.6),
   c(4.7,3.1))
 plot(pts, xlim = c(0,6), ylim = c(0,4),type = n,
   xlab = Feature 1, ylab = Feature 2)
 points(pts,pch = as.character(1:5))
 
 ## build a hierarhical tree, store as a dendrogram
 aggTree - hclust(dist(pts), method = single)
 (dend1 - JB.as.dendrogram.hclust(aggTree))
 ## NOTE: only thing I added to official version of
 ## as.dendrogram.hclust:
 ## each node has an attribute cNum, which gives
 ## the merge step at which it was formed,
 ## i.e. gives the row of the merge object which
 ## describes the formation of that node
 ## one new line near end of nMerge loop:
 ## ***
 ## *** 51,56 
 ## --- 51,60 
 ## attr(z[[x[2]]], midpoint))/2
 ##}
 ## attr(zk, height) - oHgt[k]
 ## +
 ## + ## JB added July 6 2005
 ## + attr(zk, cNum) - k
 ## +
 ##z[[k - as.character(k)]] - zk
 ##   }
 ##   z - z[[k]]
 attributes(dend1)
 attributes(dend1[[1]])
 ## here's a table relating dend1 and the cNum attribute
 ## dend1   cNum
 ## -
 ## dend14
 ## dend1[[1]]   2
 ## dend1[[2]]   3
 ## dend1[[2]][[1]]  not set
 ## dend1[[2]][[1]]  1
 
 ## use cNum attribute in edgetext
 ## following example in dendrogram documentation
 ## would really rather associate with the node than the edge
 ## but current plotting function has no notion of nodetext
 addE - function(n) {
if(!is.leaf(n)) {
  attr(n, edgePar) - list(p.col=plum)
  attr(n, edgetext) - attr(n,cNum)
}
n
 }
 dend2 - dendrapply(dend1, addE)
 ## overlays the cNum (cluster number) attribute on dendrogram
 plot(dend2, main = dend2)
 ## why does no plum polygon appear around the '4' for the root
 ## edge?
 
 ## swap order of clusters 2 and 3,
 ## i.e. 'flip' cluster 4
 dend3 - dend2
 dend3[[1]] - dend2[[2]]
 dend3[[2]] - dend2[[1]]
 plot(dend3, main = dend3)
 ## wish I could achieve with 'dend3 - flip(dend2, cNum = 4)
 
 ## swap order of cluster 1 and object 1,
 ## i.e. 'flip' cluster 3
 dend4 - dend2
 dend4[[2]][[1]] - dend2[[2]][[2]]
 dend4[[2]][[2]] - dend2[[2]][[1]]
 plot(dend4, main = dend4)
 ## wish I could achieve with 'dend4 - flip(dend2, cNum = 3)
 
 ## finally, it's clear that the midpoint attribute would also
 ## need to be modified by 'flip'
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] exact values for p-values - more information.

2005-07-11 Thread S.O. Nyangoma
Hi there,
Actually my aim was to compare anumber of extreme values (e.g. 39540) 
with df1=1, df2=7025 via p-values.

Spencer mentions that 

However, I have also used numbers like 
exp(-19775.52) to guestimate relative degrees of plausibility for 
different alternatives.

Can someone point to me an article using this method?

Regards. Stephen.

 

- Original Message -
From: Spencer Graves [EMAIL PROTECTED]
Date: Monday, July 11, 2005 7:39 pm
Subject: Re: [R] exact values for p-values - more information.

  I just checked:
 
  pf(39540, 1, 7025, lower.tail=FALSE, log.p=TRUE)
 [1] -Inf
 
  This is not correct.  With 7025 denominator degrees of 
 freedom, we 
 might use the chi-square approximation to the F distribution:
   
  pchisq(39540, 1, lower.tail=FALSE, log.p=TRUE)
 [1] -19775.52
 
  In sum, my best approximation to  pf(39540, 1, 7025, 
 lower.tail=FALSE, log.p=TRUE), given only a minute to work on 
 this, is 
 exp(pchisq(39540, 1, lower.tail=FALSE, log.p=TRUE)) = exp(-19775.52).
 
  I'm confident that many violations of assumptions would 
 likely be 
 more important than the differences between p-value:  2.2e-16 
 and 
 That doesn't mean they are right, only 
 the best 
 I can get with the available resources.
 
  spencer graves
 
 Achim Zeileis wrote:
 
  On Mon, 11 Jul 2005, S.O. Nyangoma wrote:
  
  
  Hi there,
  If I do an lm, I get p-vlues as
 
  p-value:  2.2e-16
 
 This is obtained from F =39540 with df1 = 1, df2 = 7025.
 
  Suppose am interested in exact value such as
 
  p-value = 1.6e-16 (note = and not )
 
  How do I go about it?
  
  
  You can always extract the `exact' p-value from the summary.lm 
 object or
  you can compute it by hand via
pf(39540, df1 = 1, df2 = 7025, lower.tail = FALSE)
  For all practical purposes, the above means that the p-value is 0.
  I guess you are on a 32-bit machine, then it also means that the 
 p-value
  is smaller than the Machine epsilon
.Machine$double.eps
  
  So if you want to report the p-value somewhere, I think R's 
 output should
  be more than precise enough. If you want to compute some other 
 values that
  depend on such a p-value, then it is probably wiser to compute 
 on a log
  scale, i.e. instead
pf(70, df1 = 1, df2 = 7025, lower.tail = FALSE)
  use
pf(70, df1 = 1, df2 = 7025, lower.tail = FALSE, log.p = TRUE)
  
  However, don't expect to be able to evaluate it at such extreme 
 values such as 39540.
  Z
  
  __
  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
 
 -- 
 Spencer Graves, PhD
 Senior Development Engineer
 PDF Solutions, Inc.
 333 West San Carlos Street Suite 700
 San Jose, CA 95110, USA
 
 [EMAIL PROTECTED]
 www.pdf.com http://www.pdf.com
 Tel:  408-938-4420
 Fax: 408-280-7915
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-
 guide.html

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Dependence of bundle dse on setRNG (was Misbehaviour of DSE)

2005-07-11 Thread Prof Brian Ripley
As the posting guide says, there is no R 2.1.  The first message suggests 
this is R 2.1.0, and the posting guide does ask you to use the latest 
version (and to quote versions accurately).

The dse bundle depends on package setRNG, which you have not installed,
so you need to do that.  Look at e.g.

 library(help=dse2)

which the posting guide actually asks you to do and include in your 
posting.  (And the message you got is pretty clearcut.)

(Suggestion to Paul Gilbert: it would work better to have the setRNG 
dependence in the top-level DESCRIPTION and not in the dependence for each 
package.  That way install.packages() can work out it is required from the 
information on CRAN, which is only at bundle level.)

On Mon, 11 Jul 2005, Ajay Shah wrote:

 I am finding problems with using dse:

 library(dse1)
 Loading required package: tframe
 Error: c(package '%s' required by '%s' could not be found, setRNG, dse1)
 library(dse2)
 Loading required package: setRNG
 Error: package 'setRNG' could not be loaded
 In addition: Warning message:
 there is no package called 'setRNG' in: library(pkg, character.only = TRUE, 
 logical = TRUE, lib.loc = lib.loc)

 This is on R 2.1 on an Apple ibook (OS X) panther.

-- 
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


Re: [R] Boxplot philosophy {was Boxplot in R}

2005-07-11 Thread Ted Harding
On 11-Jul-05 Martin Maechler wrote:
 AdaiR == Adaikalavan Ramasamy [EMAIL PROTECTED]
 on Mon, 11 Jul 2005 03:04:44 +0100 writes:
 
 AdaiR Just an addendum on the philosophical aspect of doing
 AdaiR this.  By selecting the 5% and 95% quantiles, you are
 AdaiR always going to get 10% of the data as extreme and
 AdaiR these points may not necessarily outliers.  So when
 AdaiR you are comparing information from multiple columns
 AdaiR (i.e.  boxplots), it is harder to say which column
 AdaiR contains more extreme value compared to others etc.
 
 Yes, indeed!
 
 People {and software implementations} have several times provided
 differing definitions of how the boxplot whiskers should be defined.
 
 I strongly believe that this is very often a very bad idea!!
 
 A boxplot should be a universal mean communication and so one
 should be *VERY* reluctant redefining the outliers.
 
 I just find that Matlab (in their statistics toolbox)
 does *NOT* use such a silly 5% / 95% definition of the whiskers,
 at least not according to their documentation.
 That's very good (and I wonder where you, Larry, got the idea of
 the 5 / 95 %).
 Using such a fixed percentage is really a very inferior idea to
 John Tukey's definition {the one in use in all implementations
 of S (including R) probably for close to 20 years now}.
 
 I see one flaw in Tukey's definition {which is shared of course
 by any silly percentage based ``outlier'' definition}:
 
The non-dependency on the sample size.
 
 If you have a 1000 (or even many more) points,
 you'll get more and more `outliers' even for perfectly normal data.
 
 But then, I assume John Tukey would have told us to do more
 sophisticated things {maybe things like the violin plots} than
 boxplot  if you have really very many data points, you may want
 to see more features -- or he would have agreed to use 
boxplot(*,  range = monotone_slowly_growing(n) )
 for largish sample sizes n.
 
 Martin Maechler, ETH Zurich

I happily agree with Martin's essay on Boxplot philiosophy.

It would cerainly confuse boxplot watchers if the interpretation
of what they saw had to vary from case to case. The fact that
careful (and necessarily detailed) explanations of what was
different this time would be necessary in the text would not
help much, and would defeat the primary objective of the boxplot
which is to present a summary of features of the data in a form
which can be grasped visually very quickly indeed.

I'm sure many of us have at times felt some frustration at the
rigidly precise numerical interpretations which Tukey imposed
on the elements of his many EDA techniques; but this did ensure
that the viewer really knew, at a glance, what he was looking at.

EDA brilliantly combined several aspects of looking at data:
selection of features of the data; highly efficient encoding of
these, and of their inter-relationships, into a medium directly
adapted to visual perception; robustness (so that the perceptions
were not unstable with respect to wondering just what the underlying
distribution might be); accessibility (in the sense of being truly
understood) to non-theoreticians; and capacity to be implemented on
primitive information technology.

Indeed, one might say that the core team of EDA consists of the
techniques for which you need only pencil and paper.

Nevertheless, Tukey was no rigid dogmatist. His objective was
always to give a good representation of the data, and he would
happily shift his ground, or adapt a technique (albeit probably
giving it a different name), or devise a new one, if that would
be useful for the case in hand.

Best wishes to all,
Ted.



E-Mail: (Ted Harding) [EMAIL PROTECTED]
Fax-to-email: +44 (0)870 094 0861
Date: 11-Jul-05   Time: 22:19:47
-- 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


Re: [R] R on kubuntu

2005-07-11 Thread Hans Gardfjell
Dear R users,

It's possible to install newer versions of R from the Ubuntu depositories,
but the new version of R can only be found in 'breezy', the 'unstable'
version of Ubuntu. You can however create a mixed system with most of your 
packages
from the old 'hoary' and only some (like R) from the unstable 'breezy' by using 
the
pinning mechanism in apt-get. Check out the instructions on pinning in the
Ubuntu wiki: https://wiki.ubuntu.com/PinningHowto?highlight=%28apt-get+-t%29

In short it works like this:

1. Put in links to other depositories in /etc/apt/sources.list
2. Set the pinning priorities in /etc/apt/preferences
3. Install newer R-base with apt-get -t breezy install r-base

Cheers,

Hans Gardfjell
Dept. of Ecology and Environmental Science
Umeå University, Sweden


 Hello all,
 
 I am planning to redeploy my workstation under KUBUNTU.
 
 Does any body has any r experience installing/using r on this platform?
 
 Best regards.
 
 
 --
 Constant Depièreux
 Managing Director
 Applied QUality Technologies Europe sprl
 Rue des Déportés 123, B-4800 Verviers
 (Tel) +32 87 292175 - (Fax) +32 87 292171 - (Mobile) +32 475 555 818
 (Web) http://www.aqte.be - (Courriel) constant.depiereux at aqte.be 
https://stat.ethz.ch/mailman/listinfo/r-help
 (Skype) cdepiereux

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Boxplot philosophy {was Boxplot in R}

2005-07-11 Thread Berton Gunter
FWIW:

I have been an enthusiastic user of boxplots for decades. Of course, the
issue of how to handle the whiskers (outliers] is a valid one, and indeed
sample size related. Dogma is always dangerous. I got to know John Tukey
somewhat (I used to chauffer him to and from meetings with a group of Merck
statisticians), and I,too,think he would have been the first to agree that
some flexibility here is wise. 

HOWEVER, the chief advantage of boxplots is their simplicity at displaying
simultaneously and easily **several** important aspects of the data, of
which outliers are probably the most problematic (as they often result in
severe distortion of the plots without careful scaling). Even with dozens of
boxplots, center, scale, and skewness are easy to discern and compare. I
think this would NOT be true of violin plots and other more complex
versions -- simplicity can be a virtue.

Finally, a tidbit for boxplot afficianados: how does one detect bimodality
from a boxplot?

-- Bert Gunter
Genentech Non-Clinical Statistics
South San Francisco, CA
 
The business of the statistician is to catalyze the scientific learning
process.  - George E. P. Box
 
 

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of Ted Harding
 Sent: Monday, July 11, 2005 2:52 PM
 To: r-help@stat.math.ethz.ch
 Subject: Re: [R] Boxplot philosophy {was Boxplot in R}
 
 On 11-Jul-05 Martin Maechler wrote:
  AdaiR == Adaikalavan Ramasamy [EMAIL PROTECTED]
  on Mon, 11 Jul 2005 03:04:44 +0100 writes:
  
  AdaiR Just an addendum on the philosophical aspect of doing
  AdaiR this.  By selecting the 5% and 95% quantiles, you are
  AdaiR always going to get 10% of the data as extreme and
  AdaiR these points may not necessarily outliers.  So when
  AdaiR you are comparing information from multiple columns
  AdaiR (i.e.  boxplots), it is harder to say which column
  AdaiR contains more extreme value compared to others etc.
  
  Yes, indeed!
  
  People {and software implementations} have several times provided
  differing definitions of how the boxplot whiskers should be defined.
  
  I strongly believe that this is very often a very bad idea!!
  
  A boxplot should be a universal mean communication and so one
  should be *VERY* reluctant redefining the outliers.
  
  I just find that Matlab (in their statistics toolbox)
  does *NOT* use such a silly 5% / 95% definition of the whiskers,
  at least not according to their documentation.
  That's very good (and I wonder where you, Larry, got the idea of
  the 5 / 95 %).
  Using such a fixed percentage is really a very inferior idea to
  John Tukey's definition {the one in use in all implementations
  of S (including R) probably for close to 20 years now}.
  
  I see one flaw in Tukey's definition {which is shared of course
  by any silly percentage based ``outlier'' definition}:
  
 The non-dependency on the sample size.
  
  If you have a 1000 (or even many more) points,
  you'll get more and more `outliers' even for perfectly normal data.
  
  But then, I assume John Tukey would have told us to do more
  sophisticated things {maybe things like the violin plots} than
  boxplot  if you have really very many data points, you may want
  to see more features -- or he would have agreed to use 
 boxplot(*,  range = monotone_slowly_growing(n) )
  for largish sample sizes n.
  
  Martin Maechler, ETH Zurich
 
 I happily agree with Martin's essay on Boxplot philiosophy.
 
 It would cerainly confuse boxplot watchers if the interpretation
 of what they saw had to vary from case to case. The fact that
 careful (and necessarily detailed) explanations of what was
 different this time would be necessary in the text would not
 help much, and would defeat the primary objective of the boxplot
 which is to present a summary of features of the data in a form
 which can be grasped visually very quickly indeed.
 
 I'm sure many of us have at times felt some frustration at the
 rigidly precise numerical interpretations which Tukey imposed
 on the elements of his many EDA techniques; but this did ensure
 that the viewer really knew, at a glance, what he was looking at.
 
 EDA brilliantly combined several aspects of looking at data:
 selection of features of the data; highly efficient encoding of
 these, and of their inter-relationships, into a medium directly
 adapted to visual perception; robustness (so that the perceptions
 were not unstable with respect to wondering just what the underlying
 distribution might be); accessibility (in the sense of being truly
 understood) to non-theoreticians; and capacity to be implemented on
 primitive information technology.
 
 Indeed, one might say that the core team of EDA consists of the
 techniques for which you need only pencil and paper.
 
 Nevertheless, Tukey was no rigid dogmatist. His objective was
 always to give a good representation of the data, and he would
 happily shift his ground, or adapt a 

[R] CIs in predict?

2005-07-11 Thread Guy Forrester
Dear All,


I am trying to put some Confidence intervals on some regressions from a linear 
model with no luck.  I can extract the fitted values using 'predict', but am 
having difficulty in getting at the confidence intervals, or the standard 
errors.

Any suggestions would be welcome

Cheers

Guy

Using Version 2.1.0  (2005-04-18) on a PC




vol.mod3 - lm(log.volume~log.area*lake,data=vol)
summary(vol.mod3)

plot(c(1.3,2.5),c(-0.7,0.45),type=n,xlab=Log area,ylab=Log volume)

areapred.a - seq(min(vol$log.area[vol$lake==a]), 
max(vol$log.area[vol$lake==a]), length=100)
areapred.b - seq(min(vol$log.area[vol$lake==b]), 
max(vol$log.area[vol$lake==b]), length=100)


preda - predict(vol.mod3, data.frame(log.area=areapred.a,interval=confidence 
,lake=rep(a,100)))

#This gives the fitted values as predicted, but no CIs
 preda
   123456   
 789 
-0.562577529 -0.553263576 -0.543949624 -0.534635671 -0.525321718 -0.516007765 
-0.506693813 -0.497379860 -0.488065907 
  10   11   12   13   14   15   
16   17   18 
-0.478751955 -0.469438002 -0.460124049 -0.450810097 -0.441496144 -0.432182191 
-0.422868239 -0.413554286 -0.404240333 
  19   20   21   22   23   24   
25   26   27 
-0.394926380 -0.385612428 -0.376298475 ETC ETC

#As does this, but with no SEs
 preda - predict(vol.mod3, data.frame(log.area=areapred.a,se.fit=T 
 ,lake=rep(a,100)))
 preda
   123456   
 789   10 
-0.562577529 -0.553263576 -0.543949624 -0.534635671 -0.525321718 -0.516007765 
-0.506693813 -0.497379860 -0.488065907 -0.478751955 
  11   12   13   14   15   16   
17   18   19   20 
-0.469438002 -0.460124049 -0.450810097 ETC ETC





Guy J Forrester
Biometrician
Manaaki Whenua - Landcare Research
PO Box 69, Lincoln, New Zealand.
Tel. +64 3 325 6701 x3738
Fax +64 3 325 2418
E-mail [EMAIL PROTECTED] 
www.LandcareResearch.co.nz 



WARNING: This email and any attachments may be confidential ...{{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


Re: [R] CIs in predict?

2005-07-11 Thread Peter Dalgaard
Guy Forrester [EMAIL PROTECTED] writes:

 Dear All,
 
 
 I am trying to put some Confidence intervals on some regressions from a 
 linear model with no luck.  I can extract the fitted values using 'predict', 
 but am having difficulty in getting at the confidence intervals, or the 
 standard errors.
 
 Any suggestions would be welcome

help(predict.lm) should get you there soon enough.

-- 
   O__   Peter Dalgaard Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark  Ph:  (+45) 35327918
~~ - ([EMAIL PROTECTED])  FAX: (+45) 35327907

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] CIs in predict?

2005-07-11 Thread Simon Blomberg

At 08:40 AM 12/07/2005, Guy Forrester wrote:
Dear All,


I am trying to put some Confidence intervals on some regressions from a 
linear model with no luck.  I can extract the fitted values using 
'predict', but am having difficulty in getting at the confidence 
intervals, or the standard errors.

Any suggestions would be welcome

Cheers

Guy

Using Version 2.1.0  (2005-04-18) on a PC




vol.mod3 - lm(log.volume~log.area*lake,data=vol)
summary(vol.mod3)

plot(c(1.3,2.5),c(-0.7,0.45),type=n,xlab=Log area,ylab=Log volume)

areapred.a - seq(min(vol$log.area[vol$lake==a]), 
max(vol$log.area[vol$lake==a]), length=100)
areapred.b - seq(min(vol$log.area[vol$lake==b]), 
max(vol$log.area[vol$lake==b]), length=100)


preda - predict(vol.mod3, 
data.frame(log.area=areapred.a,interval=confidence ,lake=rep(a,100)))

You have interval=confidence inside your call to data.frame, not inside 
your call to predict. Hence you are creating a data frame with a variable 
called interval, with one level called confidence, and predict does not see 
interval=confidence at all! See ?predict.lm.

HTH,

Simon.


#This gives the fitted values as predicted, but no CIs
  preda
12345 
6789
-0.562577529 -0.553263576 -0.543949624 -0.534635671 -0.525321718 
-0.516007765 -0.506693813 -0.497379860 -0.488065907
   10   11   12   13   14 
   15   16   17   18
-0.478751955 -0.469438002 -0.460124049 -0.450810097 -0.441496144 
-0.432182191 -0.422868239 -0.413554286 -0.404240333
   19   20   21   22   23 
   24   25   26   27
-0.394926380 -0.385612428 -0.376298475 ETC ETC

#As does this, but with no SEs
  preda - predict(vol.mod3, data.frame(log.area=areapred.a,se.fit=T 
 ,lake=rep(a,100)))
  preda
12345 
6789   10
-0.562577529 -0.553263576 -0.543949624 -0.534635671 -0.525321718 
-0.516007765 -0.506693813 -0.497379860 -0.488065907 -0.478751955
   11   12   13   14   15 
   16   17   18   19   20
-0.469438002 -0.460124049 -0.450810097 ETC ETC





Guy J Forrester
Biometrician
Manaaki Whenua - Landcare Research
PO Box 69, Lincoln, New Zealand.
Tel. +64 3 325 6701 x3738
Fax +64 3 325 2418
E-mail [EMAIL PROTECTED]
www.LandcareResearch.co.nz



WARNING: This email and any attachments may be confidential ...{{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

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Problems with R on OS X

2005-07-11 Thread Bill Northcott
On 11/07/2005, at 8:00 PM, Heinz Schild wrote:
 I used R on OS X 10.3x quite some time with no serious problems.
 Sometimes R stopped when I tried to execute a bigger program. After
 updating to OS X to version 10.4 R worked but I still had the problem
 with bigger programs. Therefore I re-installed R on top of the
 existing R version. The installation finished properly but suddenly R
 did not work. Then I reinstalled OS X 10.4 because I thought some
 left over registry information from R may caused the problem. R still
 crashed. I hoped the problem would be overcome with the new R version
 1.12. Unfortunately this is not the case. The error report (see
 appendix) says that Thread 0 crashed.
 What can I do?
 Heinz Schild

The best place for MacOS issues is the r-sig-mac list to which I have  
copied this mesage.

MacOS X does not have a registry!  Configuration files for MacOS  
applications live in ~/Library/Preferences.  The file for the R GUI  
app is called org.R-project.R.plist.  The rest of R on a Mac is just  
like any other UNIX platform with settings and history in ~/.Rdata  
and ~/.Rhistory.  (Files starting with a period '.' are invisible in  
the MacOS desktop.)  If you want to trash the configuration just  
delete these three files.

The version 1.12 you quote is only for the GUI, which is just a front  
end for the R framework.  The current version of the framework binary  
is 2.1.1, which you will see on the 'R Console' window when you start  
the application.  If you try the 1.12 GUI with a version of the  
framework other than 2.1.1 it is very likely to crash.  If you  
install the complete package from CRAN it will install both the GUI  
application and the framework.

BIll Northcott

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] nlme plot

2005-07-11 Thread R V
Hello,

I am running this script from Pinheiro  Bates book in R Version 2.1.1 (WinXP).
But, I can't plot Figure 2.3.
What's wrong?


TIA.
Rod.

-
library(nlme)
 names( Orthodont )
[1] distance age  Subject  Sex 
 levels( Orthodont$Sex )
[1] Male   Female
 OrthoFem - Orthodont[ Orthodont$Sex == Female, ]
 
 fm1OrthF - lme( distance ~ age, data = OrthoFem, random = ~ 1 | Subject )
 fm2OrthF - update( fm1OrthF, random = ~ age | Subject )
 orthLRTsim - simulate.lme( fm1OrthF, fm2OrthF, nsim = 1000 )
 plot( orthLRTsim, df = c(1, 2) )# produces Figure 2.3
Error in if ((dfType - as.double(names(x)[1])) == 1) { : 
argument is of length zero
Execution halted

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] simulate.lme plot

2005-07-11 Thread R V
Hello,
I am running this script from Pinheiro  Bates book in R (Version 2.1.1,WinXP).
But, I can't plot Figure 2.3.
What's wrong?
Thanks, Rod.

-
library(nlme)
 names( Orthodont )
[1] distance age  Subject  Sex
 levels( Orthodont$Sex )
[1] Male   Female
 OrthoFem - Orthodont[ Orthodont$Sex == Female, ]

 fm1OrthF - lme( distance ~ age, data = OrthoFem, random = ~ 1 | Subject )
 fm2OrthF - update( fm1OrthF, random = ~ age | Subject )
 orthLRTsim - simulate.lme( fm1OrthF, fm2OrthF, nsim = 1000 )
 plot( orthLRTsim, df = c(1, 2) )# produces Figure 2.3
Error in if ((dfType - as.double(names(x)[1])) == 1) { :
   argument is of length zero
Execution halted

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] transition matrix and discretized data

2005-07-11 Thread Dimitri Joe
Hi there,

I have data on earnings of 12000 individuals at two points in time. I intend to 
construct a transition matrix, where the typical element, p_ij, gives the 
probability that an individual ends at the j-th decile of the earnings 
distribution given that he was was initially at the i-th decile. Thus, this is 
a bi-stochastic matrix. The problem is that the income data is nearly discrete 
in the sense that many individuals hold the same income level at each point. 
For instance, there are 1400 individuals who earned in period-one the minimum 
positive income (say, $100). Therefore, in the first decile there will be more 
than 10% of the individuals. This happens for both periods, and for a few 
income levels. As a result, the transition matrix won't have both rows and 
columns summing to one.

The solution I've found for this problem was to generate a uniform random 
vector, with entries ranging from, say, -.0001 to .0001, and ad it to both 
earnings vectors and compute the transition matrix. Repeat the procedure 1000 
times and get the mean of the resulting matrices. The thing is I'm totally new 
to simulations. Here's part of what I'm trying to do:

# X is a two-comun data frame. column 1 is the period-one individuals' earnings 
and column two is the period-two.

n - nrow(X)#12000
sim - runif(n, -.0001, .0001)
X - X + sim
q - 10# in order to compute deciles. it could be 
quintiles, quartiles, whatever
p - seq(0,1,1/q)
f.x - quantile(X[,1],  p, names=F)
f.y - quantile(X[,2],  p, names=F)
f.x[1] - 0; f.y[1] - 0
a -  cut(X[,1], f.x, right=T)
b -  cut(X[,2], f.y, right=T)
P - table(a,b)
P - P/rowSums(P)
P

The point is that I don't know how to store the matrix P efficiently so that it 
can be averaged with the remianing 999.
Also, suggestions on how to solve the problem of discretized data are welcome.

Thanks a lot,
Dimitri


[[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


Re: [R] nlme plot

2005-07-11 Thread Douglas Bates
On 7/11/05, R V [EMAIL PROTECTED] wrote:
 Hello,
 
 I am running this script from Pinheiro  Bates book in R Version 2.1.1 
 (WinXP).
 But, I can't plot Figure 2.3.
 What's wrong?

There was a change in the way that R handles assignments of names of
components and that affected the construction of this plot. I've
rewritten the plot construction code (the logic in the current version
was bizarre) and will upload a new version of nlme after I test this
out.

By the way, there is a simpler way of reproducing the examples in our book.  Use

source(system.file(scripts/ch02.R, package = nlme))

but don't try that with the current version of the nlme package. 
There is another infelicity (to use Bill Venables' term) that I will
correct.


 
 
 TIA.
 Rod.
 
 -
 library(nlme)
  names( Orthodont )
 [1] distance age  Subject  Sex
  levels( Orthodont$Sex )
 [1] Male   Female
  OrthoFem - Orthodont[ Orthodont$Sex == Female, ]
 
  fm1OrthF - lme( distance ~ age, data = OrthoFem, random = ~ 1 | Subject )
  fm2OrthF - update( fm1OrthF, random = ~ age | Subject )
  orthLRTsim - simulate.lme( fm1OrthF, fm2OrthF, nsim = 1000 )
  plot( orthLRTsim, df = c(1, 2) )# produces Figure 2.3
 Error in if ((dfType - as.double(names(x)[1])) == 1) { :
 argument is of length zero
 Execution halted
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] exact values for p-values - more information.

2005-07-11 Thread Spencer Graves
  I don't have a reference, but you could look under data mining and p 
values / Bonferroni.

  spencer graves

S.O. Nyangoma wrote:

 Hi there,
 Actually my aim was to compare anumber of extreme values (e.g. 39540) 
 with df1=1, df2=7025 via p-values.
 
 Spencer mentions that 
 
 However, I have also used numbers like 
 exp(-19775.52) to guestimate relative degrees of plausibility for 
 different alternatives.
 
 Can someone point to me an article using this method?
 
 Regards. Stephen.
 
  
 
 - Original Message -
 From: Spencer Graves [EMAIL PROTECTED]
 Date: Monday, July 11, 2005 7:39 pm
 Subject: Re: [R] exact values for p-values - more information.
 
 
 I just checked:


pf(39540, 1, 7025, lower.tail=FALSE, log.p=TRUE)

[1] -Inf

 This is not correct.  With 7025 denominator degrees of 
freedom, we 
might use the chi-square approximation to the F distribution:
  

pchisq(39540, 1, lower.tail=FALSE, log.p=TRUE)

[1] -19775.52

 In sum, my best approximation to  pf(39540, 1, 7025, 
lower.tail=FALSE, log.p=TRUE), given only a minute to work on 
this, is 
exp(pchisq(39540, 1, lower.tail=FALSE, log.p=TRUE)) = exp(-19775.52).

 I'm confident that many violations of assumptions would 
likely be 
more important than the differences between p-value:  2.2e-16 
and 
 
  That doesn't mean they are right, only 
 
the best 
I can get with the available resources.

 spencer graves

Achim Zeileis wrote:


On Mon, 11 Jul 2005, S.O. Nyangoma wrote:



Hi there,
If I do an lm, I get p-vlues as

p-value:  2.2e-16

This is obtained from F =39540 with df1 = 1, df2 = 7025.

Suppose am interested in exact value such as

p-value = 1.6e-16 (note = and not )

How do I go about it?


You can always extract the `exact' p-value from the summary.lm 

object or

you can compute it by hand via
  pf(39540, df1 = 1, df2 = 7025, lower.tail = FALSE)
For all practical purposes, the above means that the p-value is 0.
I guess you are on a 32-bit machine, then it also means that the 

p-value

is smaller than the Machine epsilon
  .Machine$double.eps

So if you want to report the p-value somewhere, I think R's 

output should

be more than precise enough. If you want to compute some other 

values that

depend on such a p-value, then it is probably wiser to compute 

on a log

scale, i.e. instead
  pf(70, df1 = 1, df2 = 7025, lower.tail = FALSE)
use
  pf(70, df1 = 1, df2 = 7025, lower.tail = FALSE, log.p = TRUE)

However, don't expect to be able to evaluate it at such extreme 

values such as 39540.

Z

__
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

-- 
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

[EMAIL PROTECTED]
www.pdf.com http://www.pdf.com
Tel:  408-938-4420
Fax: 408-280-7915

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-
guide.html
 
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

-- 
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

[EMAIL PROTECTED]
www.pdf.com http://www.pdf.com
Tel:  408-938-4420
Fax: 408-280-7915

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Boxplot philosophy {was Boxplot in R}

2005-07-11 Thread Spencer Graves
  I'll bite:  How does one detect bimodalidty from a boxplot?

  spencer graves

Berton Gunter wrote:

 FWIW:
 
 I have been an enthusiastic user of boxplots for decades. Of course, the
 issue of how to handle the whiskers (outliers] is a valid one, and indeed
 sample size related. Dogma is always dangerous. I got to know John Tukey
 somewhat (I used to chauffer him to and from meetings with a group of Merck
 statisticians), and I,too,think he would have been the first to agree that
 some flexibility here is wise. 
 
 HOWEVER, the chief advantage of boxplots is their simplicity at displaying
 simultaneously and easily **several** important aspects of the data, of
 which outliers are probably the most problematic (as they often result in
 severe distortion of the plots without careful scaling). Even with dozens of
 boxplots, center, scale, and skewness are easy to discern and compare. I
 think this would NOT be true of violin plots and other more complex
 versions -- simplicity can be a virtue.
 
 Finally, a tidbit for boxplot afficianados: how does one detect bimodality
 from a boxplot?
 
 -- Bert Gunter
 Genentech Non-Clinical Statistics
 South San Francisco, CA
  
 The business of the statistician is to catalyze the scientific learning
 process.  - George E. P. Box
  
  
 
 
-Original Message-
From: [EMAIL PROTECTED] 
[mailto:[EMAIL PROTECTED] On Behalf Of Ted Harding
Sent: Monday, July 11, 2005 2:52 PM
To: r-help@stat.math.ethz.ch
Subject: Re: [R] Boxplot philosophy {was Boxplot in R}

On 11-Jul-05 Martin Maechler wrote:

AdaiR == Adaikalavan Ramasamy [EMAIL PROTECTED]
on Mon, 11 Jul 2005 03:04:44 +0100 writes:

AdaiR Just an addendum on the philosophical aspect of doing
AdaiR this.  By selecting the 5% and 95% quantiles, you are
AdaiR always going to get 10% of the data as extreme and
AdaiR these points may not necessarily outliers.  So when
AdaiR you are comparing information from multiple columns
AdaiR (i.e.  boxplots), it is harder to say which column
AdaiR contains more extreme value compared to others etc.

Yes, indeed!

People {and software implementations} have several times provided
differing definitions of how the boxplot whiskers should be defined.

I strongly believe that this is very often a very bad idea!!

A boxplot should be a universal mean communication and so one
should be *VERY* reluctant redefining the outliers.

I just find that Matlab (in their statistics toolbox)
does *NOT* use such a silly 5% / 95% definition of the whiskers,
at least not according to their documentation.
That's very good (and I wonder where you, Larry, got the idea of
the 5 / 95 %).
Using such a fixed percentage is really a very inferior idea to
John Tukey's definition {the one in use in all implementations
of S (including R) probably for close to 20 years now}.

I see one flaw in Tukey's definition {which is shared of course
by any silly percentage based ``outlier'' definition}:

   The non-dependency on the sample size.

If you have a 1000 (or even many more) points,
you'll get more and more `outliers' even for perfectly normal data.

But then, I assume John Tukey would have told us to do more
sophisticated things {maybe things like the violin plots} than
boxplot  if you have really very many data points, you may want
to see more features -- or he would have agreed to use 
   boxplot(*,  range = monotone_slowly_growing(n) )
for largish sample sizes n.

Martin Maechler, ETH Zurich

I happily agree with Martin's essay on Boxplot philiosophy.

It would cerainly confuse boxplot watchers if the interpretation
of what they saw had to vary from case to case. The fact that
careful (and necessarily detailed) explanations of what was
different this time would be necessary in the text would not
help much, and would defeat the primary objective of the boxplot
which is to present a summary of features of the data in a form
which can be grasped visually very quickly indeed.

I'm sure many of us have at times felt some frustration at the
rigidly precise numerical interpretations which Tukey imposed
on the elements of his many EDA techniques; but this did ensure
that the viewer really knew, at a glance, what he was looking at.

EDA brilliantly combined several aspects of looking at data:
selection of features of the data; highly efficient encoding of
these, and of their inter-relationships, into a medium directly
adapted to visual perception; robustness (so that the perceptions
were not unstable with respect to wondering just what the underlying
distribution might be); accessibility (in the sense of being truly
understood) to non-theoreticians; and capacity to be implemented on
primitive information technology.

Indeed, one might say that the core team of EDA consists of the
techniques for which you need only pencil and paper.

Nevertheless, Tukey was no rigid dogmatist. His objective was
always to give a good representation of the data, and he would
happily shift his 

[R] R CMD INSTALL use differenct c++ compiler

2005-07-11 Thread Luke
Dear R Users,

When I installed e1071 use R CMD INSTALL, I got

configure: WARNING: g++ 2.96 cannot reliably be used with this package.
configure: error: Please use a different C++ compiler.

But how to let R CMD INSTALL use a different C++ compiler? and which
C++ compiler is good?

-Luke

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] ---- Message Quarantined! ----

2005-07-11 Thread MessageRelease
MailMarshal has quarantined the following message you sent:

   Message: B42d3485f0001.0001.0002.mml
   From:[EMAIL PROTECTED]
   To:  [EMAIL PROTECTED]
   Subject: Good day


With the proliferation of mass mailer type viruses, Massey University has 
had to implement a policy that quarantines emails with attachments that could 
be potential viruses, e.g. VBS, SCR and EXE.

However, this message is eligible for Self-Service Message Release.  If this 
attachment is genuine and you are certain that it does not contain any viruses, 
please simply reply to this message without editing it.

The quarantined email will be automatically released to [EMAIL PROTECTED] 
upon receipt of your reply, else it will be deleted after 14 days.

Please note that your reply will be processed by the automated release system. 
If you have any queries or comments, please send them to [EMAIL PROTECTED]


Message release code: MMRelCode
ZYoJxwg/WQoP+jKD7y9MEPkkJa7kC4lbnw7+GjUyOmMbCbdfbt
MkS6FoVQwZpNHhLosQRcgvZoUocQ+3XTWkjwMxfFQw0Ea0zz/z
VtE7ZIuQXVj8ZKBR0Q4YvUEN2Qa/6Q5POCRmm/qTc6d+Wg==
MMEndRelCode

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[R] How to obtain Frequency Tables for Histogram outputs and Frequency Polygons

2005-07-11 Thread Tan Hui Hui Jenny
Couple of questions:
 
1. How can I obtain the frequency tables for a histogram chart?
2. Is there a short cut to obtain the frequency polygons directly without 
having to generate the frequency table and doing a line plot?
 
Thanks ina dvance for any reply.
 
j.
 
Example data:
0.3,0.229,0.218,0.211,0.207,0.21,0.212,0.232,0.312,0.289,0.486,0.475,0.274,0.478,0.284,0.477,0.281,0.29,0.29,0.32,0.28,0.404,0.489,0.406,0.267,0.456,0.443,0.457,0.411,0.444,0.401,0.468,0.309,0.425,0.312,0.3,0.464,0.457,0.442,0.469,0.423,0.403,0.432,0.379,0.361,0.37,0.47,0.374,0.387,0.381,0.311,0.257,0.381,0.162,0.229,0.221,0.25,0.164,0.368,0.367,0.332,0.165,0.152,0.478,0.489,0.482,0.167,0.164,0.164,0.15,0.149,0.368,0.366,0.361,0.359,0.359,0.364,0.364

__
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] How to obtain Frequency Tables for Histogram outputs and Frequency Polygons

2005-07-11 Thread Gabor Grothendieck
On 7/12/05, Tan Hui Hui Jenny [EMAIL PROTECTED] wrote:
 Couple of questions:
 
 1. How can I obtain the frequency tables for a histogram chart?
 2. Is there a short cut to obtain the frequency polygons directly without 
 having to generate the frequency table and doing a line plot?

res - hist(x)  # res contains various data about histogram

# create red polygons using rect
plot(range(res$breaks), range(res$counts), type=n)
rect(res$breaks[-length(res$breaks)],0,res$breaks[-1],res$counts, col=red)

__
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