[R] Re : View code of function

2006-11-24 Thread justin bem
Hi,
About code read this document by Mr Uwe Ligges
http://www.statistik.uni-dortmund.de/~ligges/R_Help_Desk_preview.pdf
 
Justin BEM
Elève Ingénieur Statisticien Economiste
BP 294 Yaoundé.
Tél (00237)9597295.



- Message d'origine 
De : Wee-Jin Goh [EMAIL PROTECTED]
À : R-Help r-help@stat.math.ethz.ch
Envoyé le : Vendredi, 24 Novembre 2006, 7h57mn 31s
Objet : Re: [R] View code of function


Hi,

Just type mean.default without the () and the code of the function  
will be displayed.

Wee-Jin

On 24 Nov 2006, at 06:31, Lize van der Merwe wrote:

 Dear list

 I need to see the code behind a function.   I used to be able to  
 see the
 code behind a function by typing e.g. mean().  Now I get the error  
 message:
 Error in mean.default() : argument x is missing, with no  
 default.   Please
 advise.

 Regards

 Lize van der Merwe








 -- 
 This e-mail and its contents are subject to the
 South African Medical Research Council
 e-mail legal notice available at http://www.mrc.ac.za/about/ 
 EmailLegalNotice.htm


 [[alternative HTML version deleted]]

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

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






___ 
Découvrez une nouvelle façon d'obtenir des réponses à toutes vos questions ! 
Profitez des connaissances, des opinions et des expériences des internaut

[[alternative HTML version deleted]]

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


[R] Code of plot.TukeyHSD

2006-11-24 Thread Lize van der Merwe


Dear list
I need to see the code of plot.TukeyHSD.  There is no .default for it.
Please help.
Lize





-- 
This e-mail and its contents are subject to the 
South African Medical Research Council
e-mail legal notice available at http://www.mrc.ac.za/about/EmailLegalNotice.htm

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


[R] barplot help needed

2006-11-24 Thread Antje
hello,

I would like to create the following barplot:

I have 4 different data sets (same length + stddev for each data point)

data1
sd1
data2
sd2
data3
sd3
data4
sd4

now, I'd like to plot in the following way:

data1[1],data2[1],data3[1],data4[1] with it's sd-values side-by-side at 
one x-axis label (named position 1) and each bar in different colors.

data1[2],data2[2],data3[2],data4[2] at the next x-axis label (named 
position 2) with the same color scheme

and so on over the whole length.

I managed to plot one set in the following way:

par(mai=c(1.5,1,1,0.6))
plotInfo - barplot(data1, las=2, ylim = c(0,plotMax+1), ylab = 
Percentage)
arrows(plotInfo,data1,plotInfo,  data1 + sd1, length=0.1, angle=90)
arrows(plotInfo,data1,plotInfo,  data1 - sd1, length=0.1, angle=90)

could anybody give me a help on this?

Antje

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


Re: [R] Re : View code of function

2006-11-24 Thread Uwe Ligges


justin bem wrote:
 Hi,
 About code read this document by Mr Uwe Ligges
 http://www.statistik.uni-dortmund.de/~ligges/R_Help_Desk_preview.pdf

I should remove it, because the updated version has been published in 
the most recent R Newsletter. Please check out that one.

Thanks,
Uwe Ligges


 Justin BEM
 Elève Ingénieur Statisticien Economiste
 BP 294 Yaoundé.
 Tél (00237)9597295.
 
 
 
 - Message d'origine 
 De : Wee-Jin Goh [EMAIL PROTECTED]
 À : R-Help r-help@stat.math.ethz.ch
 Envoyé le : Vendredi, 24 Novembre 2006, 7h57mn 31s
 Objet : Re: [R] View code of function
 
 
 Hi,
 
 Just type mean.default without the () and the code of the function  
 will be displayed.
 
 Wee-Jin
 
 On 24 Nov 2006, at 06:31, Lize van der Merwe wrote:
 
 Dear list

 I need to see the code behind a function.   I used to be able to  
 see the
 code behind a function by typing e.g. mean().  Now I get the error  
 message:
 Error in mean.default() : argument x is missing, with no  
 default.   Please
 advise.

 Regards

 Lize van der Merwe








 -- 
 This e-mail and its contents are subject to the
 South African Medical Research Council
 e-mail legal notice available at http://www.mrc.ac.za/about/ 
 EmailLegalNotice.htm


 [[alternative HTML version deleted]]

 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting- 
 guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 
   
 
   
   
 ___ 
 Découvrez une nouvelle façon d'obtenir des réponses à toutes vos questions ! 
 Profitez des connaissances, des opinions et des expériences des internaut
 
   [[alternative HTML version deleted]]
 
 
 
 
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] Code of plot.TukeyHSD

2006-11-24 Thread Prof Brian Ripley
getAnywhere(plot.TukeyHSD)

On Fri, 24 Nov 2006, Lize van der Merwe wrote:

 Dear list
 I need to see the code of plot.TukeyHSD.  There is no .default for it.
 Please help.
 Lize

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

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


Re: [R] barplot help needed

2006-11-24 Thread Jacques VESLOT
tab - do.call(rbind, list(data1, data2, data3, data4))
etype - rep(c(sd1, sd2, sd3, sd4), length(data1))
b - barplot(tab, beside=T)
arrows(unlist(b), unlist(tab) - etype, unlist(b),  unlist(tab) + etype, code=3)
---
Jacques VESLOT

CNRS UMR 8090
I.B.L (2ème étage)
1 rue du Professeur Calmette
B.P. 245
59019 Lille Cedex

Tel : 33 (0)3.20.87.10.44
Fax : 33 (0)3.20.87.10.31

http://www-good.ibl.fr
---

Antje a écrit :
 hello,
 
 I would like to create the following barplot:
 
 I have 4 different data sets (same length + stddev for each data point)
 
 data1
 sd1
 data2
 sd2
 data3
 sd3
 data4
 sd4
 
 now, I'd like to plot in the following way:
 
 data1[1],data2[1],data3[1],data4[1] with it's sd-values side-by-side at 
 one x-axis label (named position 1) and each bar in different colors.
 
 data1[2],data2[2],data3[2],data4[2] at the next x-axis label (named 
 position 2) with the same color scheme
 
 and so on over the whole length.
 
 I managed to plot one set in the following way:
 
 par(mai=c(1.5,1,1,0.6))
 plotInfo - barplot(data1, las=2, ylim = c(0,plotMax+1), ylab = 
 Percentage)
 arrows(plotInfo,data1,plotInfo,  data1 + sd1, length=0.1, angle=90)
 arrows(plotInfo,data1,plotInfo,  data1 - sd1, length=0.1, angle=90)
 
 could anybody give me a help on this?
 
 Antje
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


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


Re: [R] t.test()

2006-11-24 Thread Robin Hankin

On 23 Nov 2006, at 13:46, Peter Dalgaard wrote:

 Robin Hankin [EMAIL PROTECTED] writes:

 Hi

 I have a vector x of length n.   I am interested in x[1]
 being different from the other observations (ie x[-1]).

[snip]


 What arguments do I need to send to t.test() to test my null?



[snip]

 Alternatively, just write up the formula for the t statistic:

 x - c(23,25,29,27,30,30)
 (x[1]-mean(x[-1]))/sqrt(var(x[-1])*(1+1/(length(x)-1)))



The gotcha in Peter Dalgaard's formula is that the maximum  
likelihood estimate for the
variance, with a sample size of one, is zero.  This is why var[x[1]]  
doesn't appear.

[R reports var(1) as NA because it uses the unbiased formula with  
(n-1) on the denominator,
as documented]

Last night I derived the likelihood test for testing my null.  Consider

H1: x~N(mu_x,s^2);  y~N(mu_y,s^2)
H2: x,y~N(mu,s^2)


The support gained by allowing the two means to differ [ie compare H2  
to H1] is:

\[
E=
\frac{n}{2}\ln\left(
\frac{\sum(z_i-\overline{z})^2}{
\sum(x_i-\overline{x})^2+
\sum(y_i-\overline{y})^2
}\right)
\]

where z=c(x,y) is both sets of observations taken together.  This  
formula supposes one
uses the appropriate maximum likelihood estimate for the (common)  
variance.  Note that
the MLE  for the variance is different on H1 and H2.

Thus if E  2 we can reject H2 and if  !(E 2) we can accept (sic) H2.





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

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


Re: [R] statistics help

2006-11-24 Thread Ingmar Visser
Dear Rohan, 
Why would you want to simulate these probabilities? As far as I can tell by
your description these are all solvable analytically, see eg Kemeney,
Mirkil, Snell  Thompson, 1958, Finite Mathematical Structures. There are
undoubtedly more recent publications that cover first passage times in
Markov models, which is what you are looking for.
Hth, Ingmar

 From: Rohan Saldanha [EMAIL PROTECTED]
 Date: Thu, 23 Nov 2006 19:07:54 +
 To: r-help@stat.math.ethz.ch
 Subject: [R] statistics help
 
 hi
 
 im a bioinformatics student as i have never had any previous programming
 experience
 i need help
 
 this is the question i need to answer:
 
 Random walk model
 
 we want to model a random walk where you take a step to the left with
 probability p
 and one to the right with probability 1-p. Now assume that there is a line
 of 11 squares.
 once you are in square 0 or in square 10 the walk ends. The aim of this
 problem set is
 to write a simulation for the random walk and analyse its dynamics
 
 1. write a simulation for a random walk which allows you to calculate the
 probability of ending
  up in square 0 starting from any other square.
 
 2. analyse the probability of ending up in square 0 starting from any other
 square. Also calculate the mean time until you have reached square 1 for p=
 0.1,0.2,0.3,0.4 and 0.5. What is the probability of reaching square 10 for
 these parameters
 hint simulate each scenario 1000 times and plot on histogram
 
 
 this is the code that i have come up with but its not working very well.
 
 rw-function(sw,p,nrep){
 Z=0
 T=0
 count=0
for (i in 1:nrep)
{
   n-0
   s=sw
 
   while (s0  s10)
  {
  x-runif(1, min=0, max=1)
  if (xp) {s-s-1}
  else{s-s+1}
  print (s)
  n-n+1
  }
  count-count+n
  print (count)
 
  if (s==0) {Z-Z+1}
  else {T-T+1}
}
  PrZ=Z/nrep
  print(c(PrZ,PrZ),quote=FALSE)
 }
 
 if you could shed some light on this it woul be really helpful. please let
 me know how much
 you would like as payment aswell.
 
 thanks
 rohan
 
 _
 Eat well and eat right. Get tips on nutrition from Naini Setalvad
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] how to loop this?[Solved]

2006-11-24 Thread antonio rodriguez
Petr Pikal escribió:
 Hi

 when working with lists, considerable option is to look on lapply, 
 sapply, mapply or other *apply functions. If I am not mistaken you 
 can do

 ttt-lapply(lasker, function(x) 
 data.frame(table(substr(names(subset(x, x =4)), 1, 7

 to get a list with desired output.

 HTH
 Petr
   

Dear Petr,

It works great!!

Many Thanks,

Antonio


 On 23 Nov 2006 at 21:06, antonio rodriguez wrote:

 Date sent:Thu, 23 Nov 2006 21:06:18 +0100
 From: antonio rodriguez [EMAIL PROTECTED]
 To:   R-Help r-help@stat.math.ethz.ch
 Subject:  [R] how to loop this?

   
 Hi,

 I have the next procedure:

 t1-data.frame(table(substr(names(subset(lasker[[1]], lasker[[1]] =
 4)), 1, 7)))

 t1[1:5,]

  Var1 Freq
 1 1988-023
 2 1988-031
 3 1988-041
 4 1988-052
 5 1988-063

 How to make a new list?, dataframe? having 189 elements in the
 'lasker' list:

   str(lasker[[1]])
  'table' int [, 1:1274] 1 1 3 2 1 5 4 1 1 4 ...
  - attr(*, dimnames)=List of 1
   ..$ : chr [1:1274] 1988-01-13 1988-01-16 1988-01-20
   1988-01-25 ...
 .
 .
   str(lasker[[189]])
  'table' int [, 1:464] 1 1 4 1 4 4 6 2 3 3 ...
  - attr(*, dimnames)=List of 1
   ..$ : chr [1:464] 1999-07-21 1999-07-23 1999-07-25
   1999-07-31 ...


 I've tried:

 windows-vector(list,189)

 for (i in 1:189){
 windows[[i]]-data.frame(table(substr(names(subset(lasker[[i]], 
 lasker[[i]] = 4)), 1, 7)))
 }

 Erro en rep.int(rep.int(seq_len(nx), rep.int(rep.fac, nx)), orep) :
 invalid number of copies in rep.int()

 I think this is due because each new element (like t1) has 2 columns.
 But what kind of structure is needed to accomodate the new elements?

 Thanks

 Antonio


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


Re: [R] command option for R CMD BATCH

2006-11-24 Thread Ramon Diaz-Uriarte
Thanks.

R.

On Thursday 23 November 2006 16:32, Prof Brian Ripley wrote:
 On Thu, 23 Nov 2006, Ramon Diaz-Uriarte wrote:
  On Thursday 23 November 2006 15:44, Prof Brian Ripley wrote:
  Try this:
 
  gannet% cat month.R
  x - commandArgs()
  print(x[length(x)])
 
  gannet% R --slave --args January  month.R
  [1] January
 
  Is the above
  R --slave --args January  month.R
  the preferred way of using it?

 Yes it is.  That's exactly what --args was added to allow.

  I tend to use
 
  R --slave  month.R January
 
  instead (as a consequence of reconverting former scripts that used R CMD
  BATCH). The second call produces a ARGUMENT 'January' __ignored__ but
  otherwise seems to do the same thing.

-- 
Ramón Díaz-Uriarte
Bioinformatics 
Centro Nacional de Investigaciones Oncológicas (CNIO)
(Spanish National Cancer Center)
Melchor Fernández Almagro, 3
28029 Madrid (Spain)
Fax: +-34-91-224-6972
Phone: +-34-91-224-6900

http://ligarto.org/rdiaz
PGP KeyID: 0xE89B3462
(http://ligarto.org/rdiaz/0xE89B3462.asc)



**NOTA DE CONFIDENCIALIDAD** Este correo electrónico, y en s...{{dropped}}

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


[R] VB linked to R

2006-11-24 Thread Khaled Radhouane
help
I'd like to develop an VB 6 application and use R as
the backgroud calculator.
How do I call R with its functions and link it to the VB interface?

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


[R] Advice on Grid Data

2006-11-24 Thread Lorenzo Isella
Dear All,
I would like to automate the analysis and plotting of data taken from a grid.
Typically I deal with 2 spatial coordinates and a scalar f(x,y), but
the spatial grid is not evenly spaced at all and usually given in this
form:
x y f(x,y)
0.0  0.048979383  2.7659438106975056
0.0  0.044986665  2.603891585041688
0.0023807306 0.04787451   2.715949356768243
0.0  0.040993948  2.469223979694342
0.0023807306 0.043881793  2.5625191444824265
0.004761461  0.046769638  2.6629703119429022
0.0  0.03700123   2.361940994655468
0.0023807306 0.039889075  2.436480700580665
0.004761461  0.04277692   2.517958884562618
0.0071421918 0.045664765  2.606844303834078
0.0  0.060880877  3.470808435449538
0.0  0.05691371   3.1907723461238686
0.0020467786 0.059650626  3.3672237912200016
0.0  0.05294655   2.9558174712065237
0.0020467786 0.055683464  3.1075221272152054
0.004093557  0.05842038   3.268965886866726
0.0020467786 0.051716298  2.8929170062920653
0.004093557  0.054453213  3.029154848759894
0.006140336  0.057190128  3.1754081073235088
0.0  0.02473231   2.138648866573983
0.0  0.020556964  2.092324627395541

I tried the image plot and lattice but unsuccessfully. Now I am
reading about the sp package (
http://cran.r-project.org/src/contrib/Descriptions/sp.html ), but I
mainly would like a piece of advice about what tools to use and how to
read and plot these data (I suppose it must be common e.g. in
geography to deal with this kind of problems).
Kind Regards

Lorenzo

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


Re: [R] VB linked to R

2006-11-24 Thread Prof Brian Ripley
See the rw-FAQ Q2.18, since this must surely be under Windows (unstated).

On Thu, 23 Nov 2006, Khaled Radhouane wrote:

 help
 I'd like to develop an VB 6 application and use R as
 the backgroud calculator.
 How do I call R with its functions and link it to the VB interface?

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


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

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


Re: [R] show unique combination of two factor

2006-11-24 Thread Chuck Cleland
Aimin Yan wrote:
 p factor have 5 levels
 aa factor have 19 levels.
 totally it should have 95 combinations.
 but I just find there are 92 combinations.
 Does anyone know how to code to find what combinations are missed?

  Here is an example with fewer factor levels of one way you might do this:

df - data.frame(p = rep(c(A,B,C,D), each=10),
 aa = rep(c(Yes,No), 20))

df$aa - replace(df$aa, df$p == D, No)

table(df)
   aa
p   No Yes
  A  5   5
  B  5   5
  C  5   5
  D 10   0

names(which(with(df, table(interaction(p, aa))) == 0))
[1] D.Yes

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

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

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


Re: [R] show unique combination of two factor

2006-11-24 Thread Michael Dewey
At 22:02 23/11/2006, Aimin Yan wrote:
p factor have 5 levels
aa factor have 19 levels.
totally it should have 95 combinations.
but I just find there are 92 combinations.
Does anyone know how to code to find what combinations are missed?

Does
which(table(p,aa) == 0, arr.ind=TRUE)
do what you wanted.



Thanks,

Aimin



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

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


Re: [R] random effect question and glm

2006-11-24 Thread Michael Dewey
At 21:52 23/11/2006, Aimin Yan wrote:
consider p as random effect with 5 levels, what is difference between these
two models?

   p5.random.p - lmer(Y
~p+(1|p),data=p5,family=binomial,control=list(usePQL=FALSE,msV=1))
   p5.random.p1 - lmer(Y
~1+(1|p),data=p5,family=binomial,control=list(usePQL=FALSE,msV=1))

Well, try them and see. Then if you cannot understand the output tell us
a) what you found
b) how it differed from what you expected

in addtion, I try these two models, it seems they are same.
what is the difference between these two model. Is this a cell means model?
m00 - glm(sc ~aa-1,data = p5)
m000 - glm(sc ~1+aa-1,data = p5)

See above


thanks,

Aimin Yan



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

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


Re: [R] Advice on Grid Data

2006-11-24 Thread Roger Bivand
On Fri, 24 Nov 2006, Lorenzo Isella wrote:

 Dear All,
 I would like to automate the analysis and plotting of data taken from a grid.
 Typically I deal with 2 spatial coordinates and a scalar f(x,y), but
 the spatial grid is not evenly spaced at all and usually given in this
 form:
 x y f(x,y)
 0.0  0.048979383  2.7659438106975056
 0.0  0.044986665  2.603891585041688
 0.0023807306 0.04787451   2.715949356768243
 0.0  0.040993948  2.469223979694342
 0.0023807306 0.043881793  2.5625191444824265
 0.004761461  0.046769638  2.6629703119429022
 0.0  0.03700123   2.361940994655468
 0.0023807306 0.039889075  2.436480700580665
 0.004761461  0.04277692   2.517958884562618
 0.0071421918 0.045664765  2.606844303834078
 0.0  0.060880877  3.470808435449538
 0.0  0.05691371   3.1907723461238686
 0.0020467786 0.059650626  3.3672237912200016
 0.0  0.05294655   2.9558174712065237
 0.0020467786 0.055683464  3.1075221272152054
 0.004093557  0.05842038   3.268965886866726
 0.0020467786 0.051716298  2.8929170062920653
 0.004093557  0.054453213  3.029154848759894
 0.006140336  0.057190128  3.1754081073235088
 0.0  0.02473231   2.138648866573983
 0.0  0.020556964  2.092324627395541
 
 I tried the image plot and lattice but unsuccessfully. Now I am
 reading about the sp package (
 http://cran.r-project.org/src/contrib/Descriptions/sp.html ), but I
 mainly would like a piece of advice about what tools to use and how to
 read and plot these data (I suppose it must be common e.g. in
 geography to deal with this kind of problems).

Your data are not on a 2D grid, there are (here) 7 unique x values, but 21
unique y values of 21. You can treat the data as a SpatialPointsDataFrame
(see note in R News in 2005), but if you want to display them on an actual
grid, you will have to interpolate. For more ideas, perhaps try the 
R-sig-geo mailing list.

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

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

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


[R] Diebold Mariano Test

2006-11-24 Thread Graham Leask
Dear List

Has anyone used R to distnguish between alternative forecasting models? In 
particular 
is the Diebold Mariano test available for use within R.

Any assistance would be greatly appreciated.



Graham Leask
Lecturer in Strategy
Economics  Strategy Group
Aston University
Aston Triangle
Birmingham
B4 7ET

Tel: 0121 204 3150
E Mail: [EMAIL PROTECTED]

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


[R] axis in an image() plot

2006-11-24 Thread Ricardo Rodríguez - Your XEN ICT Team
Hi all,

I've found quite usefull colored-grid created by image() but I'm facing a doubt 
I am not able to solve.

Given the following data rectangle...

 1  2  3  4  5  6  7  8  9 10 11 12 13 14
  1 12 22  0  7  2  1  0  2  0  2  6 -3  0  3
  2  0 -1  0  9  3 -4  0  0  0  0  3  0  0  0
  3 29 45  6 12 16 85 -2  0 -3 -4 89 -1 -1  1
  4  2  9  3  6 17  3 -2 -9 -2  8 -1  0  0  0
  5 44 16 -3 21 23  3  2  1  0 -2 13 18 -5  2

I am not able to draw x and y axis labeled 1 to 14 and 1 to 5 by 1. I've tried 
a number of options by using axis() to no avail.

It will be also very helpfull to be able to draw the value of each x,y 
combination within its position in the grid. Text() seems to be the answer, but 
I am not able yet to get the correct position for each label.

Please, could you point me in the right direction or offer some example?

Thank you in advance,

--
Ricardo Rodríguez
Your XEN ICT Team

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


Re: [R] soil.texture() function with bug?

2006-11-24 Thread Sander Oom
Hi Patrick,

Not sure what the problem is from your email. Can you send a code
example which reproduces the error? Make sure you mention the version of
R you are using!

Also send your question/reply to/cc R-help as well. Then the rest of the
world is there to help you too.

Greetings,

Sander.

Patrick Kuss wrote:
 Dear Sander Oom and Jim Lemon,
 
 thanks for putting the soils.texture() function into R. However, for
 whatever reason I am not able to display the triangle correctly. Each of
 the 27 tick labels shows as c(10,20,30,40,50,60,70,80,90) and thus
 basically cover the whole triangle. Plotting points and adding graphical
 paramters like 'col.symbols' works fine. I am also able to plot my soil
 data using triax.plot() without any problems, but of course, the nice
 feature of soil type areas is not available.
 
 Do you have any insights?
 
 Thanks a lot and cheers from Alaska
 
 Patrick


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


Re: [R] barplot help needed

2006-11-24 Thread Antje
Thank you very much for your help.
I just don't understand the following line (which also gives me a 
dimension error later in the arrows command)

etype - rep(c(sd1, sd2, sd3, sd4), length(data1))

Antje

(I don't see my emails to the mailinglist anymore... just the answers 
from other people... I don't understand???)


Jacques VESLOT schrieb:
 tab - do.call(rbind, list(data1, data2, data3, data4))
 etype - rep(c(sd1, sd2, sd3, sd4), length(data1))
 b - barplot(tab, beside=T)
 arrows(unlist(b), unlist(tab) - etype, unlist(b),  unlist(tab) + etype, 
 code=3)
 ---
 Jacques VESLOT
 
 CNRS UMR 8090
 I.B.L (2ème étage)
 1 rue du Professeur Calmette
 B.P. 245
 59019 Lille Cedex
 
 Tel : 33 (0)3.20.87.10.44
 Fax : 33 (0)3.20.87.10.31
 
 http://www-good.ibl.fr
 ---
 
 Antje a écrit :
 hello,

 I would like to create the following barplot:

 I have 4 different data sets (same length + stddev for each data point)

 data1
 sd1
 data2
 sd2
 data3
 sd3
 data4
 sd4

 now, I'd like to plot in the following way:

 data1[1],data2[1],data3[1],data4[1] with it's sd-values side-by-side 
 at one x-axis label (named position 1) and each bar in different 
 colors.

 data1[2],data2[2],data3[2],data4[2] at the next x-axis label (named 
 position 2) with the same color scheme

 and so on over the whole length.

 I managed to plot one set in the following way:

 par(mai=c(1.5,1,1,0.6))
 plotInfo - barplot(data1, las=2, ylim = c(0,plotMax+1), ylab = 
 Percentage)
 arrows(plotInfo,data1,plotInfo,  data1 + sd1, length=0.1, angle=90)
 arrows(plotInfo,data1,plotInfo,  data1 - sd1, length=0.1, angle=90)

 could anybody give me a help on this?

 Antje

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



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


Re: [R] Diebold Mariano Test

2006-11-24 Thread Hannu Kahra
Not using R, but you can have a look at the RATS source code
DMARIANO.SRCavailable on page
http://www.estima.com/Forecasting.shtml and translate it into R.

Hannu

On 11/24/06, Graham Leask [EMAIL PROTECTED] wrote:

 Dear List

 Has anyone used R to distnguish between alternative forecasting models? In
 particular
 is the Diebold Mariano test available for use within R.

 Any assistance would be greatly appreciated.



 Graham Leask
 Lecturer in Strategy
 Economics  Strategy Group
 Aston University
 Aston Triangle
 Birmingham
 B4 7ET

 Tel: 0121 204 3150
 E Mail: [EMAIL PROTECTED]

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


[[alternative HTML version deleted]]

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


Re: [R] Using roman and italic fonts in an xlab expression for a plot

2006-11-24 Thread Gabor Grothendieck
or even shorter:

plot(1, xlab = Level is ~ italic(M))


On 11/23/06, Gabor Grothendieck [EMAIL PROTECTED] wrote:
 Try this:

  plot(1, xlab = quote(Level is ~ italic(M)))

 or

   plot(1, xlab = quote(Level ~ is ~ italic(M))


 On 11/23/06, Philip Boland [EMAIL PROTECTED] wrote:
  Just wondering if it is possible to put something like
 
  The excess level is M
 
  in the xlab position of a plot but where the M is in italic font and
  the rest in ordinary (or Roman) font.
 
  Thanks - Phil Boland
 
  __
  R-help@stat.math.ethz.ch mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
 


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


[R] Plot with two seperate y axis

2006-11-24 Thread Thorsten Muehge

Hello,
I would like to plot the following matrix:
Wk=x achsis.
Para 1 = left y-axis as a barplot
para 2 right y-axis as a normal scatter plat.

I could not find such a solution in any of my documentation.

Can someone help me?

Thanks a lot
Thorsten

WkPara 1  Para 2
312000  99.8
322005  99.0
332002  98.0
341090  98.5
352001  99.1
362010  97.0
372010  98.8

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


Re: [R] barplot help needed

2006-11-24 Thread Jacques VESLOT
thought sd1, sd2... were scalars but if not just do:
etype - c(sd1, sd2, sd3, sd4)
---
Jacques VESLOT

CNRS UMR 8090
I.B.L (2ème étage)
1 rue du Professeur Calmette
B.P. 245
59019 Lille Cedex

Tel : 33 (0)3.20.87.10.44
Fax : 33 (0)3.20.87.10.31

http://www-good.ibl.fr
---


Antje a écrit :
 Thank you very much for your help.
 I just don't understand the following line (which also gives me a 
 dimension error later in the arrows command)
 
 etype - rep(c(sd1, sd2, sd3, sd4), length(data1))
 
 Antje
 
 (I don't see my emails to the mailinglist anymore... just the answers 
 from other people... I don't understand???)
 
 
 Jacques VESLOT schrieb:
 
tab - do.call(rbind, list(data1, data2, data3, data4))
etype - rep(c(sd1, sd2, sd3, sd4), length(data1))
b - barplot(tab, beside=T)
arrows(unlist(b), unlist(tab) - etype, unlist(b),  unlist(tab) + etype, 
code=3)
---
Jacques VESLOT

CNRS UMR 8090
I.B.L (2ème étage)
1 rue du Professeur Calmette
B.P. 245
59019 Lille Cedex

Tel : 33 (0)3.20.87.10.44
Fax : 33 (0)3.20.87.10.31

http://www-good.ibl.fr
---

Antje a écrit :

hello,

I would like to create the following barplot:

I have 4 different data sets (same length + stddev for each data point)

data1
sd1
data2
sd2
data3
sd3
data4
sd4

now, I'd like to plot in the following way:

data1[1],data2[1],data3[1],data4[1] with it's sd-values side-by-side 
at one x-axis label (named position 1) and each bar in different 
colors.

data1[2],data2[2],data3[2],data4[2] at the next x-axis label (named 
position 2) with the same color scheme

and so on over the whole length.

I managed to plot one set in the following way:

par(mai=c(1.5,1,1,0.6))
plotInfo - barplot(data1, las=2, ylim = c(0,plotMax+1), ylab = 
Percentage)
arrows(plotInfo,data1,plotInfo,  data1 + sd1, length=0.1, angle=90)
arrows(plotInfo,data1,plotInfo,  data1 - sd1, length=0.1, angle=90)

could anybody give me a help on this?

Antje

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


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


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


Re: [R] Advice on Grid Data

2006-11-24 Thread Lorenzo Isella
On 24/11/06, Roger Bivand [EMAIL PROTECTED] wrote:
 On Fri, 24 Nov 2006, Lorenzo Isella wrote:

  Dear All,
  I would like to automate the analysis and plotting of data taken from a 
  grid.
  Typically I deal with 2 spatial coordinates and a scalar f(x,y), but
  the spatial grid is not evenly spaced at all and usually given in this
  form:
  x y 
  f(x,y)
  0.0  0.048979383  2.7659438106975056
  0.0  0.044986665  2.603891585041688
  0.0023807306 0.04787451   2.715949356768243
  0.0  0.040993948  2.469223979694342
  0.0023807306 0.043881793  2.5625191444824265
  0.004761461  0.046769638  2.6629703119429022
  0.0  0.03700123   2.361940994655468
  0.0023807306 0.039889075  2.436480700580665
  0.004761461  0.04277692   2.517958884562618
  0.0071421918 0.045664765  2.606844303834078
  0.0  0.060880877  3.470808435449538
  0.0  0.05691371   3.1907723461238686
  0.0020467786 0.059650626  3.3672237912200016
  0.0  0.05294655   2.9558174712065237
  0.0020467786 0.055683464  3.1075221272152054
  0.004093557  0.05842038   3.268965886866726
  0.0020467786 0.051716298  2.8929170062920653
  0.004093557  0.054453213  3.029154848759894
  0.006140336  0.057190128  3.1754081073235088
  0.0  0.02473231   2.138648866573983
  0.0  0.020556964  2.092324627395541
 
  I tried the image plot and lattice but unsuccessfully. Now I am
  reading about the sp package (
  http://cran.r-project.org/src/contrib/Descriptions/sp.html ), but I
  mainly would like a piece of advice about what tools to use and how to
  read and plot these data (I suppose it must be common e.g. in
  geography to deal with this kind of problems).

 Your data are not on a 2D grid, there are (here) 7 unique x values, but 21
 unique y values of 21. You can treat the data as a SpatialPointsDataFrame
 (see note in R News in 2005), but if you want to display them on an actual
 grid, you will have to interpolate. For more ideas, perhaps try the
 R-sig-geo mailing list.

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

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



Thanks for the advice, but I hope I have not been misleading: only
part of the grid coordinates are reported here (the whole list is very
long).
Does your conclusion hold anyway?
Cheers

Lorenzo

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


Re: [R] barplot help needed

2006-11-24 Thread ONKELINX, Thierry
Is the length of all your datasets equal? If not try 
etype - rep(c(sd1, sd2, sd3, sd4), c(length(data1), length(data2), 
length(data3), length(data4))

Cheers,

Thierry


ir. Thierry Onkelinx

Instituut voor natuur- en bosonderzoek / Reseach Institute for Nature and Forest

Cel biometrie, methodologie en kwaliteitszorg / Section biometrics, methodology 
and quality assurance

Gaverstraat 4

9500 Geraardsbergen

Belgium

tel. + 32 54/436 185

[EMAIL PROTECTED]

www.inbo.be 

 

Do not put your faith in what statistics say until you have carefully 
considered what they do not say.  ~William W. Watt

A statistical analysis, properly conducted, is a delicate dissection of 
uncertainties, a surgery of suppositions. ~M.J.Moroney


-Oorspronkelijk bericht-
Van: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Namens Antje
Verzonden: vrijdag 24 november 2006 13:42
Aan: r-help@stat.math.ethz.ch
Onderwerp: Re: [R] barplot help needed

Thank you very much for your help.
I just don't understand the following line (which also gives me a 
dimension error later in the arrows command)

etype - rep(c(sd1, sd2, sd3, sd4), length(data1))

Antje

(I don't see my emails to the mailinglist anymore... just the answers 
from other people... I don't understand???)


Jacques VESLOT schrieb:
 tab - do.call(rbind, list(data1, data2, data3, data4))
 etype - rep(c(sd1, sd2, sd3, sd4), length(data1))
 b - barplot(tab, beside=T)
 arrows(unlist(b), unlist(tab) - etype, unlist(b),  unlist(tab) + etype, 
 code=3)
 ---
 Jacques VESLOT
 
 CNRS UMR 8090
 I.B.L (2ème étage)
 1 rue du Professeur Calmette
 B.P. 245
 59019 Lille Cedex
 
 Tel : 33 (0)3.20.87.10.44
 Fax : 33 (0)3.20.87.10.31
 
 http://www-good.ibl.fr
 ---
 
 Antje a écrit :
 hello,

 I would like to create the following barplot:

 I have 4 different data sets (same length + stddev for each data point)

 data1
 sd1
 data2
 sd2
 data3
 sd3
 data4
 sd4

 now, I'd like to plot in the following way:

 data1[1],data2[1],data3[1],data4[1] with it's sd-values side-by-side 
 at one x-axis label (named position 1) and each bar in different 
 colors.

 data1[2],data2[2],data3[2],data4[2] at the next x-axis label (named 
 position 2) with the same color scheme

 and so on over the whole length.

 I managed to plot one set in the following way:

 par(mai=c(1.5,1,1,0.6))
 plotInfo - barplot(data1, las=2, ylim = c(0,plotMax+1), ylab = 
 Percentage)
 arrows(plotInfo,data1,plotInfo,  data1 + sd1, length=0.1, angle=90)
 arrows(plotInfo,data1,plotInfo,  data1 - sd1, length=0.1, angle=90)

 could anybody give me a help on this?

 Antje

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



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

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


Re: [R] Plot with two seperate y axis

2006-11-24 Thread David Barron
Assuming the data are in a data frame called dt, this should work:

 plot(dt$Wk,dt$Para1,type=h)
 par(new=TRUE)
 plot(dt$Wk,dt$Para2,yaxt=n)
 axis(4,at=97:100)


On 24/11/06, Thorsten Muehge [EMAIL PROTECTED] wrote:

 Hello,
 I would like to plot the following matrix:
 Wk=x achsis.
 Para 1 = left y-axis as a barplot
 para 2 right y-axis as a normal scatter plat.

 I could not find such a solution in any of my documentation.

 Can someone help me?

 Thanks a lot
 Thorsten

 WkPara 1  Para 2
 312000  99.8
 322005  99.0
 332002  98.0
 341090  98.5
 352001  99.1
 362010  97.0
 372010  98.8

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



-- 
=
David Barron
Said Business School
University of Oxford
Park End Street
Oxford OX1 1HP

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


[R] How to find AUC in SVM (kernlab package)

2006-11-24 Thread Muhammad Subianto
Dear all,
I was wondering if someone can help me. I am learning SVM for 
classification in my research with kernlab package. I want to know about 
classification performance using Area Under Curve (AUC). I know ROCR 
package can do this job but I found all example in ROCR package have 
include prediction, for example, ROCR.hiv {ROCR}. My problem is how to 
produce prediction in SVM and to find AUC.

Here is a simple example:

library(MASS)
library(kernlab)
library(ROCR)

pimamodel - ksvm(type ~ .,data=Pima.tr,type=C-svc,C=10,prob.model=TRUE)
pimamodel
fitted(pimamodel)

pima.pred - predict(pimamodel, Pima.te[,-8], type=probabilities)
pima.pred

# try to find AUC
#predid.no  - prediction(pima.pred[,1], Pima.te[,8])
#predid.yes - prediction(pima.pred[,2], Pima.te[,8])
predid - prediction(pima.pred, Pima.te[,8])
perfid - performance(predid,tpr,fpr)
perfid.auc - performance(predid,auc)
perfid.auc

Thank you very much for your help.

Best wishes, Muhammad Subianto

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


[R] Fwd: Dates Conversion/write.foreign

2006-11-24 Thread Shubha Karanth
-- Forwarded message --
From: Shubha Vishwanath Karanth [EMAIL PROTECTED]
Date: Nov 24, 2006 7:54 PM
Subject: Dates Conversion/write.foreign
To: Shubha Karanth [EMAIL PROTECTED], Shubha Vishwanath
Karanth [EMAIL PROTECTED]




Hi R experts,



I need an urgent help...



I have an a dataframe caled idat. Below i give a snapshot of it.



   Datetime Volume_a
Volume_b



(11/21/06 12:50:00) (11/21/06 12:50:00)   00

(11/21/06 13:00:00) (11/21/06 13:00:00)   00

(11/21/06 13:10:00) (11/21/06 13:10:00)   00

(11/21/06 13:20:00) (11/21/06 13:20:00)   00

(11/21/06 13:30:00) (11/21/06 13:30:00)  342113 0

(11/21/06 13:40:00) (11/21/06 13:40:00)  695071 0

(11/21/06 13:50:00) (11/21/06 13:50:00)  470943  4690

(11/21/06 14:00:00) (11/21/06 14:00:00)  870072 0

(11/21/06 14:10:00) (11/21/06 14:10:00) 1010101  2000

(11/21/06 14:20:00) (11/21/06 14:20:00)  71428750

(11/21/06 14:30:00) (11/21/06 14:30:00)  388716  1780

(11/21/06 14:40:00) (11/21/06 14:40:00)  380038  1245



The type of each variable is given below:





idat : 'data.frame':144 obs. of  3 variables:

 $ Datetime   :Classes 'chron', 'dates', 'times'  atomic
[1:144] 13473 13473 13473 13473 13473 ...

 $ Volume_a   : num  0 0 0 0 0 0 0 0 0 0 ...

 $ Volume_b   : num  0 0 0 0 0 0 0 0 0 0 ...





The variable i am interested is Datetime.



I exported this dataframe into SAS by the write.foreign command,



write.foreign(idat,Z:\\New\\idat_d.sas7bdat,Z:New\\idat_c.sas,package=SAS,dataname=f.idat)



When I export, i see only the date value (Note: no time) in my SAS
data set. Why is this?



The write.foreign package says:



Numeric variables and factors are supported for all packages, dates
and times (Date, dates, date, and POSIXt classes) are also supported
for SAS and characters are supported for SPSS.



So, what should i do to get my time value also apart from the
datevalue? Do i need to convert the chron object, Datetime into
POSIXt...or anything alse...so that SAS can read it as datetime.
format?





Thank you,

Shubha.

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


Re: [R] Fwd: Dates Conversion/write.foreign

2006-11-24 Thread ONKELINX, Thierry
Correct me if I'm wrong, but a few weeks ago my professor from the
statistical computing class told us that SAS sometimes stores dates
including the time but only displays the date. So it looks like the time
isn't stored. This was the case with data imported from Excel. In our
exercise we had to explicitly define the date as a day, month and year.

So maybe you don't have a problem at all ;-)

Cheers,

Thierry




ir. Thierry Onkelinx

Instituut voor natuur- en bosonderzoek / Reseach Institute for Nature
and Forest

Cel biometrie, methodologie en kwaliteitszorg / Section biometrics,
methodology and quality assurance

Gaverstraat 4

9500 Geraardsbergen

Belgium

tel. + 32 54/436 185

[EMAIL PROTECTED]

www.inbo.be 

 

Do not put your faith in what statistics say until you have carefully
considered what they do not say.  ~William W. Watt

A statistical analysis, properly conducted, is a delicate dissection of
uncertainties, a surgery of suppositions. ~M.J.Moroney


-Oorspronkelijk bericht-
Van: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] Namens Shubha Karanth
Verzonden: vrijdag 24 november 2006 15:32
Aan: r-help@stat.math.ethz.ch; [EMAIL PROTECTED]; [EMAIL PROTECTED]
Onderwerp: [R] Fwd: Dates Conversion/write.foreign

-- Forwarded message --
From: Shubha Vishwanath Karanth [EMAIL PROTECTED]
Date: Nov 24, 2006 7:54 PM
Subject: Dates Conversion/write.foreign
To: Shubha Karanth [EMAIL PROTECTED], Shubha Vishwanath
Karanth [EMAIL PROTECTED]




Hi R experts,



I need an urgent help...



I have an a dataframe caled idat. Below i give a snapshot of it.



   Datetime Volume_a
Volume_b



(11/21/06 12:50:00) (11/21/06 12:50:00)   00

(11/21/06 13:00:00) (11/21/06 13:00:00)   00

(11/21/06 13:10:00) (11/21/06 13:10:00)   00

(11/21/06 13:20:00) (11/21/06 13:20:00)   00

(11/21/06 13:30:00) (11/21/06 13:30:00)  342113 0

(11/21/06 13:40:00) (11/21/06 13:40:00)  695071 0

(11/21/06 13:50:00) (11/21/06 13:50:00)  470943  4690

(11/21/06 14:00:00) (11/21/06 14:00:00)  870072 0

(11/21/06 14:10:00) (11/21/06 14:10:00) 1010101  2000

(11/21/06 14:20:00) (11/21/06 14:20:00)  71428750

(11/21/06 14:30:00) (11/21/06 14:30:00)  388716  1780

(11/21/06 14:40:00) (11/21/06 14:40:00)  380038  1245



The type of each variable is given below:





idat : 'data.frame':144 obs. of  3 variables:

 $ Datetime   :Classes 'chron', 'dates', 'times'  atomic
[1:144] 13473 13473 13473 13473 13473 ...

 $ Volume_a   : num  0 0 0 0 0 0 0 0 0 0 ...

 $ Volume_b   : num  0 0 0 0 0 0 0 0 0 0 ...





The variable i am interested is Datetime.



I exported this dataframe into SAS by the write.foreign command,



write.foreign(idat,Z:\\New\\idat_d.sas7bdat,Z:New\\idat_c.sas,packag
e=SAS,dataname=f.idat)



When I export, i see only the date value (Note: no time) in my SAS
data set. Why is this?



The write.foreign package says:



Numeric variables and factors are supported for all packages, dates
and times (Date, dates, date, and POSIXt classes) are also supported
for SAS and characters are supported for SPSS.



So, what should i do to get my time value also apart from the
datevalue? Do i need to convert the chron object, Datetime into
POSIXt...or anything alse...so that SAS can read it as datetime.
format?





Thank you,

Shubha.

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

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


Re: [R] Q. about the solve function in _package_ Matrix

2006-11-24 Thread Martin Maechler
 yyan == yyan liu [EMAIL PROTECTED]
 on Wed, 22 Nov 2006 15:04:33 -0800 (PST) writes:

yyan Hi:
yyan I have some problems when I use the function solve function in a 
loop. In the following code, I have a diagonal martix ttt whose  elements 
change in every iteration in a loop. I defined a dpoMatrixclass  before the 
loop so I do not need to define this class every time in the loop. The reason 
is to save some computing time. The code is below. The inverse of the matrix 
ttt
yyan should change according to the change of ttt in the loop. However, 
the values in sigma.dpo.solve, which is the inverse of ttt does not change. 
 Can anybody figure out what is wrong with my code? 
 
Well, you should not assign to the 'x' slot in this case
since if you look at  
  str(sigma.dpo)
you see that the matrix also has cached its Cholesky factor,
and the direct slot assignment does not ``see the need'' for
recomputing the Cholesky factor.

Of course, one could argue this to be a bit unfortunate,
but then, you really should only directly assign to the slots of
an S4 object if you know what you are doing ... which you did
not :-) ;-)
A more extreme view would say this to be a design flaw in the Matrix 
*package* (not library) that the authors have to consider.
Ideally, almost any slot reassignment of such a Matrix should
automatically annihilate the 'factors' slot.

  yyan  Thanks a lot in advance!

you're welcome.
Martin Maechler, ETH Zurich

BTW, I hope that your real application does not work with
diagonal matrices, because these are much more efficiently
handled by Diagonal() and the diagonalMatrix class objects it
produces.

Here's your code amended  to do what you want.

## ---
library(Matrix)

ttt - diag(1,2)
str( sigma.dpo - as(ttt, dpoMatrix) )

## use one of these:
see.more - TRUE
see.more - FALSE

for(i in 1:5) {
cat(\n-,i,---\n)
ttt - diag(i,2)
##  assigning to 'x' slot is ``dangerous''
## If you really want to do this, you also need to
## eliminate the cashed Cholesky factor:
[EMAIL PROTECTED] - as.vector(ttt)
[EMAIL PROTECTED] - list()
##
Isigma.dpo - solve(sigma.dpo)
print(sigma.dpo)
if(see.more) ## see what's going on:
str(sigma.dpo)
print(Isigma.dpo)
}

1
## ---

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


[R] odfTable: how to escape

2006-11-24 Thread Francesco Cernuto
Dear All,

I'm appreciating odfWeave as a nice reporting tool, but I had some pain 
in producing tables with odfTable command
where the first column began with  or  such as in age class heading, 
for example:
35
35-39
40-49
50-50
 60
In this case, to avoid a content.xml error, I had to change 35 in less 
than 35 and  60 in over 60.

Anyone knows how to escape those characters while producing a table in R 
or in a document chunk?


Thanks in advance,

Francesco Cernuto

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


[R] don't put line in plot

2006-11-24 Thread fernando espindola

Hi R-user,

I have a problem when try to run the next code:

plot(prueba$IC, type=l)

but plot with type=p, there not problem, I don't know what is the 
problem?? Anybody can help me


Regards,

Fernando

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


Re: [R] barplot help needed

2006-11-24 Thread Antje
Still, there is one problem. The SD-Values don't fit to the bar they 
belong to. I made the following experiment:

  data1 - c(2,4,6,2,5)
  data2 - data1
  sd1 - c(0.5,1,1.5,1,2)
  sd2 - sd1
  tab - do.call(rbind, list(data1, data2))
  etype - c(sd1,sd2)
  b - barplot(tab, beside=T)
  arrows(unlist(b), unlist(tab) - etype, unlist(b),  unlist(tab) + 
etype, code=3)

I expect the bars with the same height and the same stddev. The height 
is okay, but the stddev is messed up...

if I do it like this:

etype - matrix(c(sd1,sd2), nrow=2, byrow=TRUE)

it works (but maybe there is an easier way...)

Antje




Jacques VESLOT schrieb:
 thought sd1, sd2... were scalars but if not just do:
 etype - c(sd1, sd2, sd3, sd4)
 ---
 Jacques VESLOT
 
 CNRS UMR 8090
 I.B.L (2ème étage)
 1 rue du Professeur Calmette
 B.P. 245
 59019 Lille Cedex
 
 Tel : 33 (0)3.20.87.10.44
 Fax : 33 (0)3.20.87.10.31
 
 http://www-good.ibl.fr
 ---
 
 
 Antje a écrit :
 Thank you very much for your help.
 I just don't understand the following line (which also gives me a 
 dimension error later in the arrows command)

 etype - rep(c(sd1, sd2, sd3, sd4), length(data1))

 Antje

 (I don't see my emails to the mailinglist anymore... just the answers 
 from other people... I don't understand???)


 Jacques VESLOT schrieb:

 tab - do.call(rbind, list(data1, data2, data3, data4))
 etype - rep(c(sd1, sd2, sd3, sd4), length(data1))
 b - barplot(tab, beside=T)
 arrows(unlist(b), unlist(tab) - etype, unlist(b),  unlist(tab) + 
 etype, code=3)
 ---
 Jacques VESLOT

 CNRS UMR 8090
 I.B.L (2ème étage)
 1 rue du Professeur Calmette
 B.P. 245
 59019 Lille Cedex

 Tel : 33 (0)3.20.87.10.44
 Fax : 33 (0)3.20.87.10.31

 http://www-good.ibl.fr
 ---

 Antje a écrit :

 hello,

 I would like to create the following barplot:

 I have 4 different data sets (same length + stddev for each data point)

 data1
 sd1
 data2
 sd2
 data3
 sd3
 data4
 sd4

 now, I'd like to plot in the following way:

 data1[1],data2[1],data3[1],data4[1] with it's sd-values side-by-side 
 at one x-axis label (named position 1) and each bar in different 
 colors.

 data1[2],data2[2],data3[2],data4[2] at the next x-axis label (named 
 position 2) with the same color scheme

 and so on over the whole length.

 I managed to plot one set in the following way:

 par(mai=c(1.5,1,1,0.6))
 plotInfo - barplot(data1, las=2, ylim = c(0,plotMax+1), ylab = 
 Percentage)
 arrows(plotInfo,data1,plotInfo,  data1 + sd1, length=0.1, angle=90)
 arrows(plotInfo,data1,plotInfo,  data1 - sd1, length=0.1, angle=90)

 could anybody give me a help on this?

 Antje

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



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



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


[R] Error Message saying .Call(R_lazyLoadDBfetch, etc.

2006-11-24 Thread Jacques VESLOT
Hi,

I got the following error message when running a function of mine doing 
intensive computations:
Erreur dans .Call(R_lazyLoadDBfetch, key, file, compressed, hook, PACKAGE = 
base) :
 référence d'argument par défaut récursive

I haven't found neither where the problem lies nor what could be recursive.

Besides, I wonder whether it might or not be a problem of memory size.

Could you please give me any suggestion on how to interpret this message?

Thanks in advance,

jacques

---
Jacques VESLOT

CNRS UMR 8090
I.B.L (2ème étage)
1 rue du Professeur Calmette
B.P. 245
59019 Lille Cedex

Tel : 33 (0)3.20.87.10.44
Fax : 33 (0)3.20.87.10.31

http://www-good.ibl.fr

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


Re: [R] barplot help needed

2006-11-24 Thread ONKELINX, Thierry
The arrows are messed up because they are partially outside the borders of the 
barplot. Try adding a good ylim to the barplot. Something like:

b - barplot(tab, beside=T, ylim = c(0, 8))

Cheers,

Thierry


ir. Thierry Onkelinx

Instituut voor natuur- en bosonderzoek / Reseach Institute for Nature and Forest

Cel biometrie, methodologie en kwaliteitszorg / Section biometrics, methodology 
and quality assurance

Gaverstraat 4

9500 Geraardsbergen

Belgium

tel. + 32 54/436 185

[EMAIL PROTECTED]

www.inbo.be 

 

Do not put your faith in what statistics say until you have carefully 
considered what they do not say.  ~William W. Watt

A statistical analysis, properly conducted, is a delicate dissection of 
uncertainties, a surgery of suppositions. ~M.J.Moroney


-Oorspronkelijk bericht-
Van: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] Namens Antje
Verzonden: vrijdag 24 november 2006 16:17
Aan: r-help@stat.math.ethz.ch
Onderwerp: Re: [R] barplot help needed

Still, there is one problem. The SD-Values don't fit to the bar they 
belong to. I made the following experiment:

  data1 - c(2,4,6,2,5)
  data2 - data1
  sd1 - c(0.5,1,1.5,1,2)
  sd2 - sd1
  tab - do.call(rbind, list(data1, data2))
  etype - c(sd1,sd2)
  b - barplot(tab, beside=T)
  arrows(unlist(b), unlist(tab) - etype, unlist(b),  unlist(tab) + 
etype, code=3)

I expect the bars with the same height and the same stddev. The height 
is okay, but the stddev is messed up...

if I do it like this:

etype - matrix(c(sd1,sd2), nrow=2, byrow=TRUE)

it works (but maybe there is an easier way...)

Antje




Jacques VESLOT schrieb:
 thought sd1, sd2... were scalars but if not just do:
 etype - c(sd1, sd2, sd3, sd4)
 ---
 Jacques VESLOT
 
 CNRS UMR 8090
 I.B.L (2ème étage)
 1 rue du Professeur Calmette
 B.P. 245
 59019 Lille Cedex
 
 Tel : 33 (0)3.20.87.10.44
 Fax : 33 (0)3.20.87.10.31
 
 http://www-good.ibl.fr
 ---
 
 
 Antje a écrit :
 Thank you very much for your help.
 I just don't understand the following line (which also gives me a 
 dimension error later in the arrows command)

 etype - rep(c(sd1, sd2, sd3, sd4), length(data1))

 Antje

 (I don't see my emails to the mailinglist anymore... just the answers 
 from other people... I don't understand???)


 Jacques VESLOT schrieb:

 tab - do.call(rbind, list(data1, data2, data3, data4))
 etype - rep(c(sd1, sd2, sd3, sd4), length(data1))
 b - barplot(tab, beside=T)
 arrows(unlist(b), unlist(tab) - etype, unlist(b),  unlist(tab) + 
 etype, code=3)
 ---
 Jacques VESLOT

 CNRS UMR 8090
 I.B.L (2ème étage)
 1 rue du Professeur Calmette
 B.P. 245
 59019 Lille Cedex

 Tel : 33 (0)3.20.87.10.44
 Fax : 33 (0)3.20.87.10.31

 http://www-good.ibl.fr
 ---

 Antje a écrit :

 hello,

 I would like to create the following barplot:

 I have 4 different data sets (same length + stddev for each data point)

 data1
 sd1
 data2
 sd2
 data3
 sd3
 data4
 sd4

 now, I'd like to plot in the following way:

 data1[1],data2[1],data3[1],data4[1] with it's sd-values side-by-side 
 at one x-axis label (named position 1) and each bar in different 
 colors.

 data1[2],data2[2],data3[2],data4[2] at the next x-axis label (named 
 position 2) with the same color scheme

 and so on over the whole length.

 I managed to plot one set in the following way:

 par(mai=c(1.5,1,1,0.6))
 plotInfo - barplot(data1, las=2, ylim = c(0,plotMax+1), ylab = 
 Percentage)
 arrows(plotInfo,data1,plotInfo,  data1 + sd1, length=0.1, angle=90)
 arrows(plotInfo,data1,plotInfo,  data1 - sd1, length=0.1, angle=90)

 could anybody give me a help on this?

 Antje

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



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



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

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

[R] Using multiple objects in a for loop

2006-11-24 Thread Daniel Yanosky
I have a situation where I want to perform the same manipulations to multiple R 
objects in series. I have constructed a vector () to serve as a list of the 
objects. However, my assign() assigns the object name (on[j]) as a character 
string to temp and not the object itself.

===

fn-c(FT1.1.01.RData,FT1.1.02.RData)

on-c(FT1101,FT1102)

b-1
e-2

for (i in 1:2){

for (j in b:e){

load(fn[j])
assign(temp,on[j])

.
.
.
rm(paste(stderr[i],.t,sep=))
rm(temp)
rm(fn[i])
}

}

===

Daniel

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


Re: [R] Rgraphviz -404 Page not found

2006-11-24 Thread j.joshua thomas
Again i have problem in locating the package for clique-graphs
I tried with BioConductor under Browse for packages, it doesn't work atall.

Kindly guid me

Thanks
JJ



On 8/23/06, Seth Falcon [EMAIL PROTECTED] wrote:

 j.joshua thomas [EMAIL PROTECTED] writes:

  Dear Robert,
 
  Thanks for your time.
  I have downloaded Rgraphviz (windows binary) from www.bioconductor.org
  and put inside R2.3.0 library then  i installed from the local zip
  its says package 'graph' couldnot be loaded.
 
  Am i doing the installation correctly? Still the new user.
 
  Can you guide me sir?

 Questions about BioC packages are best directed to the bioconductor
 mailing list.

 I would recommend trying:

 source(http://bioconductor.org/biocLite.R;)
 biocLite(Rgraphviz)


 Best,

 + seth

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




-- 
Lecturer J. Joshua Thomas
KDU College Penang Campus
Research Student,
University Sains Malaysia

[[alternative HTML version deleted]]

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


Re: [R] barplot help needed

2006-11-24 Thread Jacques VESLOT
tab - do.call(rbind, list(data1, data2, data3, data4))
etype - do.call(rbind, list(sd1, sd2, sd3, sd4))
b - barplot(tab, beside=T, ylim=c(0,max(tab+etype)))
arrows(as.vector(b), as.vector(tab) - as.vector(etype), as.vector(b),  
as.vector(tab) + 
as.vector(etype), code=3)

unlist() is not correct - sorry - since all are matrices - not data frames !
so use as.vector()
---
Jacques VESLOT

CNRS UMR 8090
I.B.L (2ème étage)
1 rue du Professeur Calmette
B.P. 245
59019 Lille Cedex

Tel : 33 (0)3.20.87.10.44
Fax : 33 (0)3.20.87.10.31

http://www-good.ibl.fr
---


Antje a écrit :
 Still, there is one problem. The SD-Values don't fit to the bar they 
 belong to. I made the following experiment:
 
   data1 - c(2,4,6,2,5)
   data2 - data1
   sd1 - c(0.5,1,1.5,1,2)
   sd2 - sd1
   tab - do.call(rbind, list(data1, data2))
   etype - c(sd1,sd2)
   b - barplot(tab, beside=T)
   arrows(unlist(b), unlist(tab) - etype, unlist(b),  unlist(tab) + 
 etype, code=3)
 
 I expect the bars with the same height and the same stddev. The height 
 is okay, but the stddev is messed up...
 
 if I do it like this:
 
 etype - matrix(c(sd1,sd2), nrow=2, byrow=TRUE)
 
 it works (but maybe there is an easier way...)
 
 Antje
 
 
 
 
 Jacques VESLOT schrieb:
 
thought sd1, sd2... were scalars but if not just do:
etype - c(sd1, sd2, sd3, sd4)
---
Jacques VESLOT

CNRS UMR 8090
I.B.L (2ème étage)
1 rue du Professeur Calmette
B.P. 245
59019 Lille Cedex

Tel : 33 (0)3.20.87.10.44
Fax : 33 (0)3.20.87.10.31

http://www-good.ibl.fr
---


Antje a écrit :

Thank you very much for your help.
I just don't understand the following line (which also gives me a 
dimension error later in the arrows command)

etype - rep(c(sd1, sd2, sd3, sd4), length(data1))

Antje

(I don't see my emails to the mailinglist anymore... just the answers 
from other people... I don't understand???)


Jacques VESLOT schrieb:


tab - do.call(rbind, list(data1, data2, data3, data4))
etype - rep(c(sd1, sd2, sd3, sd4), length(data1))
b - barplot(tab, beside=T)
arrows(unlist(b), unlist(tab) - etype, unlist(b),  unlist(tab) + 
etype, code=3)
---
Jacques VESLOT

CNRS UMR 8090
I.B.L (2ème étage)
1 rue du Professeur Calmette
B.P. 245
59019 Lille Cedex

Tel : 33 (0)3.20.87.10.44
Fax : 33 (0)3.20.87.10.31

http://www-good.ibl.fr
---

Antje a écrit :


hello,

I would like to create the following barplot:

I have 4 different data sets (same length + stddev for each data point)

data1
sd1
data2
sd2
data3
sd3
data4
sd4

now, I'd like to plot in the following way:

data1[1],data2[1],data3[1],data4[1] with it's sd-values side-by-side 
at one x-axis label (named position 1) and each bar in different 
colors.

data1[2],data2[2],data3[2],data4[2] at the next x-axis label (named 
position 2) with the same color scheme

and so on over the whole length.

I managed to plot one set in the following way:

par(mai=c(1.5,1,1,0.6))
plotInfo - barplot(data1, las=2, ylim = c(0,plotMax+1), ylab = 
Percentage)
arrows(plotInfo,data1,plotInfo,  data1 + sd1, length=0.1, angle=90)
arrows(plotInfo,data1,plotInfo,  data1 - sd1, length=0.1, angle=90)

could anybody give me a help on this?

Antje

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


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


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


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


Re: [R] How to find AUC in SVM (kernlab package)

2006-11-24 Thread Amir Safari
Hi
  you need predict.ksvm() function.
  
  for more information see The kernlab package here:
   http://lib.stat.cmu.edu/R/CRAN/doc/packages/kernlab.pdf
  
  cheers,
  Amir
  

Muhammad Subianto [EMAIL PROTECTED] wrote:  Dear all,
I was wondering if someone can help me. I am learning SVM for 
classification in my research with kernlab package. I want to know about 
classification performance using Area Under Curve (AUC). I know ROCR 
package can do this job but I found all example in ROCR package have 
include prediction, for example, ROCR.hiv {ROCR}. My problem is how to 
produce prediction in SVM and to find AUC.

Here is a simple example:

library(MASS)
library(kernlab)
library(ROCR)

pimamodel - ksvm(type ~ .,data=Pima.tr,type=C-svc,C=10,prob.model=TRUE)
pimamodel
fitted(pimamodel)

pima.pred - predict(pimamodel, Pima.te[,-8], type=probabilities)
pima.pred

# try to find AUC
#predid.no  - prediction(pima.pred[,1], Pima.te[,8])
#predid.yes - prediction(pima.pred[,2], Pima.te[,8])
predid - prediction(pima.pred, Pima.te[,8])
perfid - performance(predid,tpr,fpr)
perfid.auc - performance(predid,auc)
perfid.auc

Thank you very much for your help.

Best wishes, Muhammad Subianto

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


 
-

[[alternative HTML version deleted]]

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


Re: [R] How to find AUC in SVM (kernlab package)

2006-11-24 Thread Muhammad Subianto
On this day 11/24/2006 05:03 PM, Amir Safari wrote:
 Hi
 you need predict.ksvm() function.

Yes, like this is below I do to predict (kernlab package, 
http://www.jstatsoft.org/v11/i09/v11i09.pdf):
pima.pred - predict(pimamodel, Pima.te[,-8], type=probabilities)
pima.pred

My problems is how to find AUC (with ROCR package, or other ROC 
functions) from predict above.

Regards, Muhammad Subianto

 for more information see The kernlab package here:
  http://lib.stat.cmu.edu/R/CRAN/doc/packages/kernlab.pdf

 cheers,
 Amir


 */Muhammad Subianto [EMAIL PROTECTED]/* wrote:

 Dear all,
 I was wondering if someone can help me. I am learning SVM for
 classification in my research with kernlab package. I want to know
 about
 classification performance using Area Under Curve (AUC). I know ROCR
 package can do this job but I found all example in ROCR package have
 include prediction, for example, ROCR.hiv {ROCR}. My problem is
 how to
 produce prediction in SVM and to find AUC.

 Here is a simple example:

 library(MASS)
 library(kernlab)
 library(ROCR)

 pimamodel - ksvm(type ~
 .,data=Pima.tr,type=C-svc,C=10,prob.model=TRUE)
 pimamodel
 fitted(pimamodel)

 pima.pred - predict(pimamodel, Pima.te[,-8], type=probabilities)
 pima.pred

 # try to find AUC
 #predid.no - prediction(pima.pred[,1], Pima.te[,8])
 #predid.yes - prediction(pima.pred[,2], Pima.te[,8])
 predid - prediction(pima.pred, Pima.te[,8])
 perfid - performance(predid,tpr,fpr)
 perfid.auc - performance(predid,auc)
 perfid.auc

 Thank you very much for your help.

 Best wishes, Muhammad Subianto

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


 
 Cheap Talk? Check out 
 http://us.rd.yahoo.com/mail_us/taglines/postman8/*http://us.rd.yahoo.com/evt=39663/*http://voice.yahoo.com
  
 Yahoo! Messenger's low PC-to-Phone call rates.

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


Re: [R] Error Message saying .Call(R_lazyLoadDBfetch, etc.

2006-11-24 Thread Peter Dalgaard
Jacques VESLOT [EMAIL PROTECTED] writes:

 Hi,
 
 I got the following error message when running a function of mine doing 
 intensive computations:
 Erreur dans .Call(R_lazyLoadDBfetch, key, file, compressed, hook, PACKAGE = 
 base) :
  référence d'argument par défaut récursive
 
 I haven't found neither where the problem lies nor what could be recursive.
 
 Besides, I wonder whether it might or not be a problem of memory size.
 
 Could you please give me any suggestion on how to interpret this message?

The prototypical situation is this

 f-function(a=a)a
 f(2)
[1] 2
 f()
Error in f() : recursive default argument reference

The default for argument a ends up needing the value of a. This can
also happens more indirectly when multiple arguments refer to
eachother. 

BTW: Is that error message proper French? Should it not be

référence récursive d'argument par défaut 

or so.

-- 
   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
and provide commented, minimal, self-contained, reproducible code.


[R] Splitting criterion in tree and rpart

2006-11-24 Thread Paolo Radaelli
Dear all,
I'm interested in fitting a classification tree by using a user-specified 
splitting criterion. If I am right,  Tree package allows to use only deviance 
or gini while rpart package offers  anova, poisson, class or exp. and 
tries to make an intelligent guess if method is missing. 
I also found that, in rpart, method can be a list of functions named init, 
split and eval but I really don't know how can I do it.  Does someone can 
suggest me some documentation I can have a look do define my splitting 
criterion and pass it to rpart ?
Thanks
Paolo

Paolo Radaelli
Dipartimento di Metodi Quantitativi per le Scienze Economiche ed Aziendali
Facoltà di Economia
Università degli Studi di Milano-Bicocca
e-mail [EMAIL PROTECTED]

[[alternative HTML version deleted]]

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


Re: [R] Rgraphviz -404 Page not found

2006-11-24 Thread Seth Falcon
j.joshua thomas [EMAIL PROTECTED] writes:

 Again i have problem in locating the package for clique-graphs
 I tried with BioConductor under Browse for packages, it doesn't work atall.

 Kindly guid me

I think you already got an answer on the Bioconductor list: RBGL is
likely the package you are looking for.

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


Re: [R] SNA packages/network package

2006-11-24 Thread Seth Falcon
 Bagatti Davide wrote:
 Hello everyone,

 I am an italian student who is working with packages SNA (social network
 analysis) and network.
 I ask  you if there is a simple way to write a R-script using these packages
 to extract from an adjacency matrix the following things:

 -number of cliques which are in the network;
 -number of k-cores (e.g. 3-cores, 4-cores);
 -cycles (e.g. 3-cycles, 4-cycles)
 -hub  authorities.

You might take a look at the 'graph' and 'RBGL' packages (both are
Bioconductor packages).  

There is a graphAM class that can be initialized with an adjacency
matrix.  Most functions in RBGL take either any graph instance or a
graphNEL.  So you might have (approx. and untested):

library(graph)
library(RBGL)
myAdjMat
g1 = new(graphAM, myAdjMat)
g2 = as(g1, graphNEL)
## call some RBGL function on g2

+ seth

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


[R] Error in Calling C++ function from R!!!

2006-11-24 Thread gsmatos1
Hello, 

I tried to call an external function of R from the following code in C++: 

void prodgdot(double *x, double *y, int *n, double *output) 

{ 

int i; 

*output=0; 

for (i=0;i*n;i++) 

{ 

*output+=x[i]*y[i]; 

} 

} 

I compiled it using from my working directory in linux terminal and I think 
it's ok: 

[EMAIL PROTECTED]:~/mysrc/meus_testes_iniciais$ R CMD SHLIB codigoprova.cc 
# The output 
g++ -I/usr/share/R/include -I/usr/share/R/include  -fpic  -g -O2 -c 
codigoprova.cc -o codigoprova.o 
g++ -shared  -o codigoprova.so codigoprova.o   -L/usr/lib/R/lib -lR 

After, I tried to call R from my working directory an Error occurs: 

 x-c(1,4,6,2) 
 y-c(3,2.4,1,9) 
 dyn.load(codigoprova.so) 
 is.loaded(codigoprova.so) 
[1] FALSE 
 product-.C(prodgdot,myx=x,muy=y,myn=NROW(x),myoutput=as.double(0)) 
Error in .C(prodgdot, myx = x, muy = y, myn = NROW(x), myoutput = 
as.double(0)) : 
C symbol name prodgdot not in load table 
 

Does anyone know what is the problem? 

Thank's in advance! 
Gilberto. 


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


[R] McQuitty method

2006-11-24 Thread Jeff Miller
Is anyone able to explain the McQuitty method in hclust?

 

Thanks in advance,

Jeff Miller

 

 

 


[[alternative HTML version deleted]]

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


Re: [R] Error in Calling C++ function from R!!!

2006-11-24 Thread Barry Rowlingson

product-.C(prodgdot,myx=x,muy=y,myn=NROW(x),myoutput=as.double(0)) 
 
 Error in .C(prodgdot, myx = x, muy = y, myn = NROW(x), myoutput = 
 as.double(0)) : 
 C symbol name prodgdot not in load table 
 
 
 Does anyone know what is the problem? 

  C++ name mangling?

http://cran.r-project.org/doc/manuals/R-exts.html#Interfacing-C_002b_002b-code

  Solution: wrap in extern C { ... } as discussed there.

Barry

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


[R] Random effects ANOVA

2006-11-24 Thread Luis Salasar
Hi,

I used the command lmer, in package lme4, to fit a random effects ANOVA. But i 
didn't get the p-values of significance tests of variance components. Does 
anyone know how  to do it? Thanks,

Luis Ernesto


-

[[alternative HTML version deleted]]

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


Re: [R] using nested ifelse and rowSums to create new variable?

2006-11-24 Thread Tony N. Brown
Thanks Dimitris and Marc. The solution to my problem was found in a 
a combination of your suggestions. I created a count vector 
separately and referenced it as part of the nested ifelse 
statements (see below).

Cheers,
Tony

 count - rowSums(DF[, 3:6]  4)
 count
 1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20
 3  1  2  2  3  2  2  0  1  2  3  2  0  3  3  2  3  2  2  4
 DF2-data.frame(DF, count)
 DF2
   x1 x2 x3 x4 x5 x6 count
1   5  1  1  1  2  5 3
2   5  5  5  5  4  1 1
3   5  1  1  5  4  1 2
4   5  5  1  4  5  1 2
5   1  1  1  2  4  1 3
6   1  1  1  4  1  5 2
7   5  5  1  4  1  5 2
8   1  5  5  5  5  5 0
9   5  1  5  4  2  5 1
10  5  1  1  1  4  5 2
11  1  1  1  1  4  1 3
12  5  5  5  1  2  5 2
13  5  5  5  4  4  5 0
14  5  1  1  2  1  5 3
15  1  5  1  1  2  5 3
16  1  1  5  5  2  1 2
17  1  1  5  2  2  1 3
18  1  5  5  2  4  1 2
19  5  1  1  2  5  5 2
20  1  1  1  2  2  1 4
 dep-with(DF2,
+   ifelse((x1==5)  (x2==5), -2,
+   ifelse((x1==1  x2==1), -1,
+   ifelse((x1==1  x2==5) | (x1==5  x2==1) 
+   (count==0), 0,
+   ifelse((x1==1  x2==5) | (x1==5  x2==1) 
+   (count==1), 1,
+   ifelse((x1==1  x2==5) | (x1==5  x2==1) 
+   (count==2), 2,
+   ifelse((x1==1  x2==5) | (x1==5  x2==1) 
+   (count==3), 3,
+   ifelse((x1==1  x2==5) | (x1==5  x2==1) 
+   (count==4), 4,
+ 99
 table(dep)
dep
-2 -1  0  1  2  3
 5  6  3  1  3  2


--On Tuesday, November 21, 2006 3:42 PM -0600 Marc Schwartz 
[EMAIL PROTECTED] wrote:

 On Tue, 2006-11-21 at 15:24 -0600, Marc Schwartz wrote:
 On Tue, 2006-11-21 at 14:26 -0600, Tony N. Brown wrote:
  Dear R-help community,
 
  If I have a data.frame df as follows:
 
   df
 x1 x2 x3 x4 x5 x6
  1   5  5  1  1  2  1
  2   5  5  5  5  1  5
  3   1  5  5  5  5  5
  4   5  5  1  4  5  5
  5   5  1  5  2  4  1
  6   5  1  5  4  5  1
  7   5  1  5  4  4  5
  8   5  1  1  1  1  5
  9   1  5  1  1  2  5
  10  5  1  5  4  5  5
  11  1  5  5  2  1  1
  12  5  5  5  4  4  1
  13  1  5  1  4  4  1
  14  1  1  5  4  5  5
  15  1  5  5  4  5  1
  16  1  1  5  5  5  1
  17  5  5  5  2  2  5
  18  1  5  1  5  5  5
  19  5  5  5  2  4  5
  20  1  1  5  2  4  5
 
  How can I create a variable that captures the pattern of
  responses  and counts across rows?
 
  I used the ifelse function and that works fine for the first
  two  conditions (see R code below). But I need help figuring
  out how to  count the number of scores in each row for columns
  x3, x4, x5, and  x6 that are less than 4, conditional upon an
  ifelse. I then want to  assign a value to the new variable
  based upon the count.
 
  The new variable I want to create is called dep. Here's my R
  code:
 
  dep-with(df,
ifelse((x1==5)  (x2==5), 0,
ifelse((x1==1  x2==1), 1,
 
ifelse((x1==1  x2==5) | (x1==5  x2==1) 
(rowSums(df[ ,c(x3, x4, x5, x6)]4) ==1), 2,
ifelse((x1==1  x2==5) | (x1==5  x2==1) 
(rowSums(df[ ,c(x3, x4, x5, x6)]4) ==2), 3,
ifelse((x1==1  x2==5) | (x1==5  x2==1) 
(rowSums(df[ ,c(x3, x4, x5, x6)]4) ==3), 4,
 
  99))
 
  dep
   0  1  2 99
   6  3  6  5
 
  I expected dep to range from 0 to 4 and its length to be equal
  to  20.
 
  Thanks in advance for your help and suggestions.

 I may be a bit off in what your desired end result is, but
 perhaps this will provide some insight. I renamed your data
 frame to DF, since df is a function:

  apply(DF[, -(1:2)], 1, function(x) sum(x  4))
  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20
  4  1  0  1  2  1  0  3  3  0  3  1  2  0  1  1  2  1  1  1


 Thanks to Dimitris' post tweaking my neurons, the above can of
 course be simplified to:

 rowSums(DF[, -(1:2)]  4)
  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20
  4  1  0  1  2  1  0  3  3  0  3  1  2  0  1  1  2  1  1  1

 Note that rowSums() will accept an all numeric data frame as an
 argument, so the initial coercion is not required.

 Hence:

 TAB - table(rowSums(DF[, -(1:2)]  4),
  factor(rowSums(DF[, 1:2]),
 labels = c((1,1), (5,1)|(1,5), (5,5

 TAB - addmargins(TAB)


 TAB

   (1,1) (5,1)|(1,5) (5,5) Sum
   0   1   3 0   4
   1   2   3 4   9
   2   0   2 1   3
   3   0   3 0   3
   4   0   0 1   1
   Sum 3  11 6  20


 HTH,

 Marc
 Time for more coffee...





---
Tony N. Brown, Ph.D.
Assistant Professor of Sociology
Vanderbilt University
Phone:  (615) 322-7518
Fax:(615) 322-7505
Email:  [EMAIL PROTECTED]

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


Re: [R] don't put line in plot

2006-11-24 Thread Carmen Meier
fernando espindola schrieb:
 Hi R-user,

 I have a problem when try to run the next code:

 plot(prueba$IC, type=l)

 but plot with type=p, there not problem, I don't know what is the 
 problem?? Anybody can help me

normally there is no problem with type=l maybe you should send a 
reproducible code

With regards Carmen



[[alternative HTML version deleted]]

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


[R] low-variance warning in lmer

2006-11-24 Thread Ben Bolker

  For block effects with small variance, lmer will sometimes
estimate the variance as being very close to zero and issue
a warning.  I don't have a problem with this -- I've explored
things a bit with some simulations (see below) and conclude that
this is probably inevitable when trying to incorporate
random effects with not very much data (the means and medians
of estimates are plausibly close to the nominal values:
the fraction of runs with warnings/near-zero estimates drops
from about 50% when the between-block variance is 1% of
the total (with 2 treatments, 12 blocks nested within treatment,
3 replicates per block), to 15% when between=30% of total,
to near zero when between=50% of total.

  My question is: what should I suggest to students when this
situation comes up?  Can anyone point me to appropriate references?
(I haven't found anything relevant in Pinheiro and Bates, but
I may not have looked in the right place ...)

  thanks
Ben Bolker


self-contained but unnecessarily complicated simulation
code/demonstration:
---
  library(lme4)
library(lattice)

simfun - function(reefeff,ntreat=2,nreef=12,
   nreefpertreat=3,
   t.eff=10,
   totvar=25,seed=NA) {
  if (!is.na(seed)) set.seed(seed)
  ntot = nreef*nreefpertreat
  npertreat=ntot/ntreat
  reef = gl(nreef,nreefpertreat)
  treat = gl(ntreat,npertreat)
  r.sd = sqrt(totvar*reefeff)
  e.sd = sqrt(totvar*(1-reefeff))
  y.det = ifelse(treat==1,0,t.eff)
  r.vals = rnorm(nreef,sd=r.sd)
  e.vals = rnorm(ntot,sd=e.sd)
  y - y.det+r.vals[as.numeric(reef)]+e.vals
  data.frame(y,treat,reef)
}

getranvar - function(x) {
  vc - VarCorr(x)
  c(diag(vc[[1]]),attr(vc,sc)^2)
}

estfun - function(reefeff,...) {
  x - simfun(reefeff=reefeff,...)
  ow - options(warn=2)
  f1 - try(lmer(y~treat+(1|reef),data=x))
  w - (class(f1)==try-error  length(grep(effectively zero,f1))0)
  options(ow)
  f2 - lmer(y~treat+(1|reef),data=x)
  c(getranvar(f2),as.numeric(w))
}


rvec - rep(c(0.01,0.05,0.1,0.15,0.2,0.3,0.5),each=100)
X - t(sapply(rvec,estfun))
colnames(X) - c(reefvar,resvar,warn)
rfrac - X[,reefvar]/(X[,reefvar]+X[,resvar])
fracwarn - tapply(X[,warn],rvec,mean)
est.mean - tapply(rfrac,rvec,mean)


op - par(mfrow=c(1,2))
plot(rvec,rfrac,type=n,xlim=c(-0.02,0.55),axes=FALSE)
axis(side=2)
box()
boxplot(rfrac~rvec,at=unique(rvec),add=TRUE,pars=list(boxwex=0.03),
col=gray)
points(jitter(rvec),rfrac,col=X[,warn]+1)
lines(unique(rvec),fracwarn,col=blue,type=b,lwd=2)
text(0.4,0.1,fraction\nzero,col=blue)
abline(a=0,b=1,lwd=2)
points(unique(rvec),est.mean,cex=1.5,col=5)
##
plot(rvec,rfrac,type=n,xlim=c(-0.02,0.55),axes=FALSE,log=y)
axis(side=2)
axis(side=1,at=unique(rvec))
box()
points(jitter(rvec),rfrac,col=X[,warn]+1)
curve(1*x,add=TRUE,lwd=2)

print(densityplot(~rfrac,groups=rvec,from=0,to=0.9,
panel=function(...) {
  panel.densityplot(...)
  cols - trellis.par.get(superpose.line)$col
  lpoints(unique(rvec),rep(8,length(rvec)),type=h,
  col=cols[1:length(rvec)])
}))




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


[R] [R-pkgs] New package `np' - nonparametric kernel smoothing methods for mixed datatypes

2006-11-24 Thread Jeffrey Racine
Dear R users,

A new package titled `np' is now available from CRAN.

The package implements recently developed kernel methods that seamlessly
handle the mix of continuous, unordered, and ordered factor datatypes
often found in applied settings. 

The package also allows users to create their own
nonparametric/semiparametric routines using high-level function calls
(via the function npksum()) rather than writing their own C or Fortran
code. Much of the code underlying the package is written in C including
the function npksum().

Currently, a range of methods can be found in the package including

- multivariate nonparametric unconditional and conditional density
estimation

- multivariate nonparametric conditional mean and gradient estimation
(local constant and local linear)

- multivariate nonparametric conditional distribution, quantile and
gradient estimation 

- nonparametric model specification tests for testing correctness of
parametric regression models

- nonparametric significance tests for nonparametric regression

- semiparametric multivariate regression methods (partially linear,
index models, average derivative estimation, varying/smooth coefficient
models)

A function npplot() is implemented that allows users to visualize the
resulting objects.

A variety of methods for computing standard errors and error bounds are
implemented including asymptotic and resampling-based methods (the
latter employing the `boot' package which is required).

A range of examples can be found in the examples section of the
accompanying help files.

We caution potential users that multivariate data-driven bandwidth
selection methods can be numerically intensive. For this reason we are
now using some of the functionality contained in the Rmpi package to
develop an MPI-aware version of the np package that we have tentatively
titled the `npRmpi' package.

The np package implements a number of methods found in the newly
released publication by Li and Racine (2007) titled Nonparametric
Econometrics: Theory and Practice, Princeton University Press.

See press.princeton.edu/titles/8355.html for further details.

Best regards, Jeff.

-- 
Professor J. S. Racine Phone:  (905) 525 9140 x 23825
Department of EconomicsFAX:(905) 521-8232
McMaster Universitye-mail: [EMAIL PROTECTED]
1280 Main St. W.,Hamilton, URL:
http://www.economics.mcmaster.ca/racine/
Ontario, Canada. L8S 4M4

`The generation of random numbers is too important to be left to chance'

___
R-packages mailing list
R-packages@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-packages

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


[R] Nonlinear statistical modeling -- a comparison of R and AD Model Builder

2006-11-24 Thread dave fournier

There has recently been some discussion on the list about
AD Model builder and the suitability of R for constructing the
types of models used in fisheries management.

   https://stat.ethz.ch/pipermail/r-help/2006-January/086841.html

   https://stat.ethz.ch/pipermail/r-help/2006-January/086858.html

I  think that many R users understimate the numerical challenges
that some of the typical nonlinear statistical model used in different
fields present. R may not be a suitable platform for development for
such models.

Around 10 years ago John Schnute, Laura Richards, and Norm Olsen
with Canadian federal fisheries undertook an investigation
comparing various statistical modeling packages for a simple
age-structured statistical model of the type commonly used in
fisheries. They compared AD Mdel Builder, Gauss, Matlab, and
Splus. Unfortunately a working model could not be produced with Splus
so its times could not be included in the comparison. It is possible
to produce a working model with the present day version of R so that
R can now be directly compared with AD Model Builder for this type of model.

I have put the results of the test together with the original
Schnute and Richards paper and the working R and AD Model Builder
codes on Otter's web site

 http://otter-rsch.ca/tresults.htm

The results are that AD Model builder is roughly 1000 times faster than
R for this problem. ADMB takes about 2 seconds to converge while
R takes over 90 minutes.

This is a simple toy example. Real fisheries models are often hundred of
times more computationally intensive as this one.

Cheers,

 Dave
~
-- 
David A. Fournier
P.O. Box 2040,
Sidney, B.C. V8l 3S3
Canada
Phone/FAX 250-655-3364
http://otter-rsch.com

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


[R] Sunflower plot error; how to deal with NA

2006-11-24 Thread Farrel Buchinsky
I suspect the problem stems from the fact that there are a couple of NA
values.

 sunflowerplot(lastoto,maxear)

Error in rep.int(i.multi, number[number  1]) :
invalid number of copies in rep.int()
So I used the subset command to get rid of the cases with NA

hell-subset(ChinOtoMayB,is.na(lastoto)==FALSE)

Then it worked perfectly

sunflowerplot(hell$lastoto,hell$maxear)

Is there a method with greater finesse to deal with NA values? In
correlations one can state use=complete.obs
But that did not work in sunflowerplot

-- 
Farrel Buchinsky

[[alternative HTML version deleted]]

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


Re: [R] Nonlinear statistical modeling -- a comparison of R and AD Model Builder

2006-11-24 Thread Mike Prager
dave fournier [EMAIL PROTECTED] wrote:

 I  think that many R users understimate the numerical challenges
 that some of the typical nonlinear statistical model used in different
 fields present. R may not be a suitable platform for development for
 such models.
 
 Around 10 years ago John Schnute, Laura Richards, and Norm Olsen
 with Canadian federal fisheries undertook an investigation
 comparing various statistical modeling packages for a simple
 age-structured statistical model of the type commonly used in
 fisheries. [...] It is possible
 to produce a working model with the present day version of R so that
 R can now be directly compared with AD Model Builder for this type of model.
 
 The results are that AD Model builder is roughly 1000 times faster than
 R for this problem. ADMB takes about 2 seconds to converge while
 R takes over 90 minutes.

Our group's experiences reflect, at least qualitatively, what
Dave says above.  We use R for analyzing results from models
written in his AD Model Builder, and a couple of years ago, we
started programming one of our models directly in R.  We quickly
abandoned that idea because of lengthy execution time under R.
That is not a judgement of either piece of software.  R and ADMB
are designed for different types of task, and it seems to me
that they complement each other well.

That experience was in part the genesis of our X2R software (now
at CRAN -- pardon the plug), which saves results from ADMB
models into a format that R can read as a list.  We feel that
now we have the best of both worlds -- fast execution with ADMB,
followed by the programming ease and excellent graphics of R for
analysis of results and projections under numerous scenarios.

-- 
Mike Prager, NOAA, Beaufort, NC
* Opinions expressed are personal and not represented otherwise.
* Any use of tradenames does not constitute a NOAA endorsement.

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


Re: [R] Nonlinear statistical modeling -- a comparison of R and AD Model Builder

2006-11-24 Thread Tony Plate
Did you try supplying gradient information to nlminb?  (I note that 
nlminb is used for the optimization, but I don't see any gradient 
information supplied to it.) I would suspect that supplying gradient 
information would greatly speed up the computation (as you note in 
comments at http://otter-rsch.ca/tresults.htm.)

I'm curious -- when you say R may not be a suitable platform for 
development for such models, what aspect of R do you feel is lacking? 
Is it the specific optimization routines available, or is it some other 
more general aspect?

Also, another optimization algorithm available in R is the L-BFGS-B 
method for optim() in the MASS package.  I've had extremely good 
experiences with using this code in S-PLUS.  It can take box 
constraints, and can use gradient information.  It is my first choice 
for most optimization problems, and I believe it is very widely used. 
Did you try using that optimization routine with this problem?

-- Tony Plate

dave fournier wrote:
 There has recently been some discussion on the list about
 AD Model builder and the suitability of R for constructing the
 types of models used in fisheries management.
 
https://stat.ethz.ch/pipermail/r-help/2006-January/086841.html
 
https://stat.ethz.ch/pipermail/r-help/2006-January/086858.html
 
 I  think that many R users understimate the numerical challenges
 that some of the typical nonlinear statistical model used in different
 fields present. R may not be a suitable platform for development for
 such models.
 
 Around 10 years ago John Schnute, Laura Richards, and Norm Olsen
 with Canadian federal fisheries undertook an investigation
 comparing various statistical modeling packages for a simple
 age-structured statistical model of the type commonly used in
 fisheries. They compared AD Mdel Builder, Gauss, Matlab, and
 Splus. Unfortunately a working model could not be produced with Splus
 so its times could not be included in the comparison. It is possible
 to produce a working model with the present day version of R so that
 R can now be directly compared with AD Model Builder for this type of model.
 
 I have put the results of the test together with the original
 Schnute and Richards paper and the working R and AD Model Builder
 codes on Otter's web site
 
  http://otter-rsch.ca/tresults.htm
 
 The results are that AD Model builder is roughly 1000 times faster than
 R for this problem. ADMB takes about 2 seconds to converge while
 R takes over 90 minutes.
 
 This is a simple toy example. Real fisheries models are often hundred of
 times more computationally intensive as this one.
 
 Cheers,
 
  Dave
 ~

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


Re: [R] axis in an image() plot

2006-11-24 Thread jim holtman
If you data is a matrix, then try:

image(1:5, 1:14, data.rect)
text(row(data.rect), col(data.rect), data.rect)


On 11/24/06, Ricardo Rodríguez - Your XEN ICT Team [EMAIL PROTECTED] wrote:
 Hi all,

 I've found quite usefull colored-grid created by image() but I'm facing a 
 doubt I am not able to solve.

 Given the following data rectangle...

 1  2  3  4  5  6  7  8  9 10 11 12 13 14
  1 12 22  0  7  2  1  0  2  0  2  6 -3  0  3
  2  0 -1  0  9  3 -4  0  0  0  0  3  0  0  0
  3 29 45  6 12 16 85 -2  0 -3 -4 89 -1 -1  1
  4  2  9  3  6 17  3 -2 -9 -2  8 -1  0  0  0
  5 44 16 -3 21 23  3  2  1  0 -2 13 18 -5  2

 I am not able to draw x and y axis labeled 1 to 14 and 1 to 5 by 1. I've 
 tried a number of options by using axis() to no avail.

 It will be also very helpfull to be able to draw the value of each x,y 
 combination within its position in the grid. Text() seems to be the answer, 
 but I am not able yet to get the correct position for each label.

 Please, could you point me in the right direction or offer some example?

 Thank you in advance,

 --
 Ricardo Rodríguez
 Your XEN ICT Team

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



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

What is the problem you are trying to solve?

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


[R] vector problem

2006-11-24 Thread bunny , lautloscrew.com
Hello out there,

i am not yet that experienced and trying to my best on a real survey.  
but i am stuck with a little matrix / vector problem.
my vector of answers could have a length of  3 or only one. i want to  
rbind all the answers into one matrix. (one vector for each participant)
answers vectors for one participant could look like:

p1:   100
p2:  20 80
p3: 40 10 50

i have the following loop which should rbind them but i get no proper  
matrix.
i´d like to have something like this

100 0 0
20 80 0
40 10 50

Here´s my loop:


for(s in 1:length(fr))
{
### ansmini is a complete survey of one participant
ansmini3=answers[relevant[s,],]

for(x in 1:length(qidsb3))
{
# outputs all answer ids out of all data which belong to the 
qids  
out of the question id vector
ansmin3=ansmini3[,1][ansmini3[,3]==qidsb3[x]]
newb3=rbind(newb3,ansmin3)
}
}

this is doing the job for question with only one answer, so i have  
only one answer id.
in this new case i have 1-3 possible answers.. that´s why i´am not  
getting a nice newb3 matrix...

i´d really be happy about some advice!
thanks in advance

matthias



[[alternative HTML version deleted]]

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


Re: [R] Random effects ANOVA

2006-11-24 Thread John Fox
Dear Luis,

I believe that lmer() doesn't print out standard errors or Wald tests for
variance-covariance components because these tests aren't reliable. You can
omit each random effect from the model in turn, and use anova() to perform a
likelihood-ratio test for the corresponding variance and covariance
components, comparing each null model with the full model.

I hope this helps,
 John


John Fox
Department of Sociology
McMaster University
Hamilton, Ontario
Canada L8S 4M4
905-525-9140x23604
http://socserv.mcmaster.ca/jfox 
 

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of Luis Salasar
 Sent: Friday, November 24, 2006 1:57 PM
 To: r-help@stat.math.ethz.ch
 Subject: [R] Random effects ANOVA
 
 Hi,
 
 I used the command lmer, in package lme4, to fit a random 
 effects ANOVA. But i didn't get the p-values of significance 
 tests of variance components. Does anyone know how  to do it? Thanks,
 
 Luis Ernesto
 
   
 -
 
   [[alternative HTML version deleted]]
 
 __
 R-help@stat.math.ethz.ch mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide 
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] vector problem

2006-11-24 Thread Gabor Grothendieck
merge.zoo can do that:

 p1 - 100
 p2 -  c(20, 80)
 p3 - c(40, 10, 50)
 library(zoo)
 t(merge(p1 = zoo(p1), p2 = zoo(p2), p3 = zoo(p3), fill = 0))
 1  2  3
p1 100  0  0
p2  20 80  0
p3  40 10 50

On 11/24/06, bunny , lautloscrew.com [EMAIL PROTECTED] wrote:
 my vector of answers could have a length of  3 or only one. i want to
 rbind all the answers into one matrix. (one vector for each participant)
 answers vectors for one participant could look like:

 p1:   100
 p2:  20 80
 p3: 40 10 50

 i have the following loop which should rbind them but i get no proper
 matrix.
 i´d like to have something like this

 100 0 0
 20 80 0
 40 10 50

 Here´s my loop:


 for(s in 1:length(fr))
{
### ansmini is a complete survey of one participant
ansmini3=answers[relevant[s,],]

for(x in 1:length(qidsb3))
{
# outputs all answer ids out of all data which belong to the 
 qids
 out of the question id vector
ansmin3=ansmini3[,1][ansmini3[,3]==qidsb3[x]]
newb3=rbind(newb3,ansmin3)
}
}

 this is doing the job for question with only one answer, so i have
 only one answer id.
 in this new case i have 1-3 possible answers.. that´s why i´am not
 getting a nice newb3 matrix...

 i´d really be happy about some advice!
 thanks in advance

 matthias

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


Re: [R] Sunflower plot error; how to deal with NA

2006-11-24 Thread Farrel Buchinsky
I only specified one of the variables and it worked fine since I was using
the subset command.

On 11/24/06, Farrel Buchinsky [EMAIL PROTECTED] wrote:

 I suspect the problem stems from the fact that there are a couple of NA
 values.

  sunflowerplot(lastoto,maxear)

 Error in rep.int(i.multi, number[number  1]) :
 invalid number of copies in rep.int()
 So I used the subset command to get rid of the cases with NA

 hell-subset(ChinOtoMayB,is.na(lastoto)==FALSE)

 Then it worked perfectly

 sunflowerplot(hell$lastoto,hell$maxear)

 Is there a method with greater finesse to deal with NA values? In
 correlations one can state use=complete.obs
 But that did not work in sunflowerplot

 --
 Farrel Buchinsky




-- 
Farrel Buchinsky
Mobile: (412) 779-1073

[[alternative HTML version deleted]]

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


Re: [R] Nonlinear statistical modeling -- a comparison of R and AD Model Builder

2006-11-24 Thread dave fournier


   Dave
  Did you try supplying gradient information to nlminb?  (I note that 
nlminb is used for the optimization, but I don't see any gradient 
information supplied to it.) I would suspect that supplying gradient 
information would greatly speed up the computation (as you note in 
comments at http://otter-rsch.ca/tresults.htm.)

Actually you should probably ask Norm Olsen these questions.
I am not proficient in R and am just using his code.

However I can say that providing derivatives for such a model is a 
highly nontrivial exercise. As I said in my posting, the  R script and 
data are available to anyone who feels that the exercise was not carried 
out properly and would like to improve on it. Also one does not need
to provide derivatives to the AD Model Builder program.

Finally suppose that you are very good at calculating derivatives and 
manage to get them right. Then someone else comes along who wants to 
modify the model. Unless they are also very good at calculating 
derivatives there will be trouble.

 
  I'm curious -- when you say R may not be a suitable platform for 
development for such models, what aspect of R do you feel is lacking? 
Is it the specific optimization routines available, or is it some other 
more general aspect?

2 seconds vs 90 minutes. For a real problem of tihs type the timings 
would probably be something like 10 minutes vs more than 2,700 minutes.

 
  Also, another optimization algorithm available in R is the L-BFGS-B 
method for optim() in the MASS package.  I've had extremely good 
experiences with using this code in S-PLUS.  It can take box 
constraints, and can use gradient information.  It is my first choice 
for most optimization problems, and I believe it is very widely used. 
Did you try using that optimization routine with this problem?
 
  -- Tony Plate
 
  dave fournier wrote:
  There has recently been some discussion on the list about
  AD Model builder and the suitability of R for constructing the
  types of models used in fisheries management.
 
 https://stat.ethz.ch/pipermail/r-help/2006-January/086841.html
 
 https://stat.ethz.ch/pipermail/r-help/2006-January/086858.html
 
  I  think that many R users understimate the numerical challenges
  that some of the typical nonlinear statistical model used in different
  fields present. R may not be a suitable platform for development for
  such models.
 
  Around 10 years ago John Schnute, Laura Richards, and Norm Olsen
  with Canadian federal fisheries undertook an investigation
  comparing various statistical modeling packages for a simple
  age-structured statistical model of the type commonly used in
  fisheries. They compared AD Mdel Builder, Gauss, Matlab, and
  Splus. Unfortunately a working model could not be produced with Splus
  so its times could not be included in the comparison. It is possible
  to produce a working model with the present day version of R so that
  R can now be directly compared with AD Model Builder for this type 
of model.
 
  I have put the results of the test together with the original
  Schnute and Richards paper and the working R and AD Model Builder
  codes on Otter's web site
 
   http://otter-rsch.ca/tresults.htm
 
  The results are that AD Model builder is roughly 1000 times faster than
  R for this problem. ADMB takes about 2 seconds to converge while
  R takes over 90 minutes.
 
  This is a simple toy example. Real fisheries models are often hundred of
  times more computationally intensive as this one.
 
  Cheers,
 
   Dave
  ~
 


-- 
David A. Fournier
P.O. Box 2040,
Sidney, B.C. V8l 3S3
Canada
Phone/FAX 250-655-3364
http://otter-rsch.com

-- 
David A. Fournier
P.O. Box 2040,
Sidney, B.C. V8l 3S3
Canada
Phone/FAX 250-655-3364
http://otter-rsch.com

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


Re: [R] Nonlinear statistical modeling -- a comparison of R and AD Model Builder

2006-11-24 Thread Spencer Graves
Hi, Mike  Dave: 

  Have you considered nonlinear mixed effects models for the types 
of problems considered in the comparison paper you cite?  Those 
benchmark trials consider T years of data ... for A age classes and 
the total number of parameters is m = T+A+5.  Without knowing more 
about the problem, I suspect that the T year parameters and the A age 
class parameters might be better modeled as random effects.  If this 
were done, the optimization problem would then involve 7 parameters, the 
5 fixed-effect parameters suggested by the computation of m plus two 
variance parameters, one for the random year effects and another for 
the random age class effect.  This would replace the problem of 
maximizing, e.g., a likelihood over T+A+5 parameters with one of 
maximizing a marginal likelihood over 2+5 parameters after integrating 
out the T and A random effects. 

  These integrations may not be easy, and I might stick with the 
fixed-effects solution if I couldn't get answers in the available time 
using a model I thought would be theoretically more appropriate.  Also, 
I might use the fixed-effects solution to get starting values for an 
attempt to maximize a more appropriate marginal likelihood.  For the 
latter, I might first try 'nlmle'.  If that failed, I might explore 
Markov Chain Monte Carlo (MCMC).  I have not done MCMC myself, but the 
MCMCpack R package looks like it might make it feasible for the types 
of problems considered in this comparison.  The CRAN summary of that 
package led me to an Adobe Acrobat version of a PPT slide presentation 
that seemed to consider just this type of problem (e.g., 
http://mcmcpack.wustl.edu/files/MartinQuinnMCMCpackslides.pdf). 

  Have you considered that? 
  Hope this helps. 
  Spencer Graves

Mike Prager wrote:
 dave fournier [EMAIL PROTECTED] wrote:

   
 I  think that many R users understimate the numerical challenges
 that some of the typical nonlinear statistical model used in different
 fields present. R may not be a suitable platform for development for
 such models.

 Around 10 years ago John Schnute, Laura Richards, and Norm Olsen
 with Canadian federal fisheries undertook an investigation
 comparing various statistical modeling packages for a simple
 age-structured statistical model of the type commonly used in
 fisheries. [...] It is possible
 to produce a working model with the present day version of R so that
 R can now be directly compared with AD Model Builder for this type of model.

 The results are that AD Model builder is roughly 1000 times faster than
 R for this problem. ADMB takes about 2 seconds to converge while
 R takes over 90 minutes.
 

 Our group's experiences reflect, at least qualitatively, what
 Dave says above.  We use R for analyzing results from models
 written in his AD Model Builder, and a couple of years ago, we
 started programming one of our models directly in R.  We quickly
 abandoned that idea because of lengthy execution time under R.
 That is not a judgement of either piece of software.  R and ADMB
 are designed for different types of task, and it seems to me
 that they complement each other well.

 That experience was in part the genesis of our X2R software (now
 at CRAN -- pardon the plug), which saves results from ADMB
 models into a format that R can read as a list.  We feel that
 now we have the best of both worlds -- fast execution with ADMB,
 followed by the programming ease and excellent graphics of R for
 analysis of results and projections under numerous scenarios.



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


Re: [R] low-variance warning in lmer

2006-11-24 Thread Gregor Gorjanc
Ben Bolker bolker at zoo.ufl.edu writes:
   For block effects with small variance, lmer will sometimes
...
   My question is: what should I suggest to students when this
 situation comes up?  Can anyone point me to appropriate references?
 (I haven't found anything relevant in Pinheiro and Bates, but
 I may not have looked in the right place ...)

Gelman has quite some similar examples. Take a look at his schools example in
the Bayesian book.

Gregor

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


[R] nlme

2006-11-24 Thread dave fournier

-- 
David A. Fournier
P.O. Box 2040,
Sidney, B.C. V8l 3S3
Canada
Phone/FAX 250-655-3364
http://otter-rsch.com

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


[R] nlme sorry about that

2006-11-24 Thread dave fournier
Sorry list I twitched and sent the wrong stuff.
Maybe I had enough fun for today.

-- 
David A. Fournier
P.O. Box 2040,
Sidney, B.C. V8l 3S3
Canada
Phone/FAX 250-655-3364
http://otter-rsch.com

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


[R] trimming plotting area in a boxplot

2006-11-24 Thread Patrick Kuss
Hello all,

I am trying to trim the plotting area in a boxplot, such that the space 
between the plot margins (left and right) are of identical size as 
between the boxes.
In the example below I want to get rid of the space outside of the 
abline().

I appreciate any suggestions.

Factor = LETTERS[1:5]
A = c(1:5); B = c(2:6); C = c(1:5); D = c(3:7); E = c(1:5); F = c(4:8)
df = data.frame(A,B,C,D,E,F)
boxplot(df)
abline(v=c(0.5,6.5))

Cheers, Patrick

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


[R] Multiple Conditional Tranformations

2006-11-24 Thread Muenchen, Robert A (Bob)
Greetings,

 

I'm learning R and I'm stuck on a basic concept: how to specify a
logical condition once and then perform multiple transformations under
that condition. The program below is simplified to demonstrate the goal.
Its results are exactly what I want, but I would like to check the
logical state of gender only once and create both (or any number of)
scores at once.

 

mystring-

(id,group,gender,q1,q2,q3,q4

01,1,f,2,2,5,4

02,2,f,2,1,4,5

03,1,f,2,2,4,4

04,2,f,1,1,5,5

05,1,m,4,5,4,

06,2,m,5,4,5,5

07,1,m,3,3,4,5

08,2,m,5,5,5,4)

 

mydata-read.table(textConnection(mystring),header=TRUE,sep=,,row.name
s=id)

mydata 

 

#Create score1 so that it differs for males and females:

mydata$score1 - ifelse( mydata$gender==f , 

   (mydata$score1 - (2*mydata$q1)+mydata$q2),

   ifelse( mydata$gender==m,

  (mydata$score1 - (3*mydata$q1)+mydata$q2), NA )

   )

mydata

 

#Create score2 so that it too differs for males and females:

mydata$score2 - ifelse( mydata$gender==f , 

   (mydata$score2 - (2.5*mydata$q1)+mydata$q2),

   ifelse( mydata$gender==m,

  (mydata$score2 - (3.5*mydata$q1)+mydata$q2), NA )

   )

mydata

 

 

Thanks!

Bob

=
Bob Muenchen (pronounced Min'-chen), Manager 
Statistical Consulting Center
U of TN Office of Information Technology
200 Stokely Management Center, Knoxville, TN 37996-0520
Voice: (865) 974-5230 
FAX: (865) 974-4810
Email: [EMAIL PROTECTED]
Web: http://oit.utk.edu/scc http://oit.utk.edu/scc , 
News: http://listserv.utk.edu/archives/statnews.html
http://listserv.utk.edu/archives/statnews.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
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Rgraphviz -404 Page not found

2006-11-24 Thread j.joshua thomas
THX Group.

On 11/25/06, Seth Falcon [EMAIL PROTECTED] wrote:

 j.joshua thomas [EMAIL PROTECTED] writes:

  Again i have problem in locating the package for clique-graphs
  I tried with BioConductor under Browse for packages, it doesn't work
 atall.
 
  Kindly guid me

 I think you already got an answer on the Bioconductor list: RBGL is
 likely the package you are looking for.




-- 
Lecturer J. Joshua Thomas
KDU College Penang Campus
Research Student,
University Sains Malaysia

[[alternative HTML version deleted]]

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


Re: [R] Multiple Conditional Tranformations

2006-11-24 Thread Muenchen, Robert A (Bob)
Mark,

Here's what I get when I try that approach.

Thanks,
Bob

 mydata$score1-numeric(mydata$q1) #just initializing.
 mydata$score2-numeric(mydata$q1)
 mydata$score1-NA
 mydata$score2-NA
 mydata
  group gender q1 q2 q3 q4 score1 score2
1 1  f  2  2  5  4 NA NA
2 2  f  2  1  4  5 NA NA
3 1  f  2  2  4  4 NA NA
4 2  f  1  1  5  5 NA NA
5 1  m  4  5  4 NA NA NA
6 2  m  5  4  5  5 NA NA
7 1  m  3  3  4  5 NA NA
8 2  m  5  5  5  4 NA NA
 mydata$score1[mydata$gender == f]-2*mydata$q1 + mydata$q2
Warning message:
number of items to replace is not a multiple of replacement length 
 mydata$score2[mydata$gender == f]-2.5*mydata$q1 + mydata$q2
Warning message:
number of items to replace is not a multiple of replacement length 
 mydata$score1[mydata$gender == m]-3*mydata$q1 + mydata$q2
Warning message:
number of items to replace is not a multiple of replacement length 
 mydata$score2[mydata$gender == m]-3.5*mydata$q1 + mydata$q2
Warning message:
number of items to replace is not a multiple of replacement length 


-Original Message-
From: Leeds, Mark (IED) [mailto:[EMAIL PROTECTED] 
Sent: Friday, November 24, 2006 8:45 PM
To: Muenchen, Robert A (Bob)
Subject: RE: [R] Multiple Conditional Tranformations

I'm not sure if I understand your question but I don't think you need
iflelse statements.

myscore-numeric(q1) ( because I'm not sure how to initialize a list so
initialize a vector with q1 elements )

myscore-NA ( I think this should set all the values in myscore to NA )
myscore[mydata$gender == f]-2*mydata$q1 + mydata$q2
myscore[mydata$gender == m]-3*mydata$q1 + mydata$q2

the above should do what you do in the first part of your code but I
don't know if that was your question ?
also, it does it making myscore a vector because I didn't know how to
initialize a list.
Someone else may goive a better solution. I'm no expert.


-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Muenchen, Robert
A (Bob)
Sent: Friday, November 24, 2006 8:27 PM
To: r-help@stat.math.ethz.ch
Subject: [R] Multiple Conditional Tranformations

Greetings,

 

I'm learning R and I'm stuck on a basic concept: how to specify a
logical condition once and then perform multiple transformations under
that condition. The program below is simplified to demonstrate the goal.
Its results are exactly what I want, but I would like to check the
logical state of gender only once and create both (or any number of)
scores at once.

 

mystring-

(id,group,gender,q1,q2,q3,q4

01,1,f,2,2,5,4

02,2,f,2,1,4,5

03,1,f,2,2,4,4

04,2,f,1,1,5,5

05,1,m,4,5,4,

06,2,m,5,4,5,5

07,1,m,3,3,4,5

08,2,m,5,5,5,4)

 

mydata-read.table(textConnection(mystring),header=TRUE,sep=,,row.name
s=id)

mydata 

 

#Create score1 so that it differs for males and females:

mydata$score1 - ifelse( mydata$gender==f , 

   (mydata$score1 - (2*mydata$q1)+mydata$q2),

   ifelse( mydata$gender==m,

  (mydata$score1 - (3*mydata$q1)+mydata$q2), NA )

   )

mydata

 

#Create score2 so that it too differs for males and females:

mydata$score2 - ifelse( mydata$gender==f , 

   (mydata$score2 - (2.5*mydata$q1)+mydata$q2),

   ifelse( mydata$gender==m,

  (mydata$score2 - (3.5*mydata$q1)+mydata$q2), NA )

   )

mydata

 

 

Thanks!

Bob

=
Bob Muenchen (pronounced Min'-chen), Manager Statistical Consulting
Center U of TN Office of Information Technology 200 Stokely Management
Center, Knoxville, TN 37996-0520
Voice: (865) 974-5230
FAX: (865) 974-4810
Email: [EMAIL PROTECTED]
Web: http://oit.utk.edu/scc http://oit.utk.edu/scc ,
News: http://listserv.utk.edu/archives/statnews.html
http://listserv.utk.edu/archives/statnews.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
and provide commented, minimal, self-contained, reproducible code.


This is not an offer (or solicitation of an offer) to buy/sell the
securities/instruments mentioned or an official confirmation.  Morgan
Stanley may deal as principal in or own or act as market maker for
securities/instruments mentioned or may advise the issuers.  This is not
research and is not from MS Research but it may refer to a research
analyst/research report.  Unless indicated, these views are the author's
and may differ from those of Morgan Stanley research or others in the
Firm.  We do not represent this is accurate or complete and we may not
update this.  Past performance is not indicative of future returns.  For
additional information, research reports and important disclosures,
contact me or see 

Re: [R] Multiple Conditional Tranformations

2006-11-24 Thread Muenchen, Robert A (Bob)
Mark,

I finally got that approach to work by spreading the logical condition
everywhere. That gets the lengths to match. Still, I can't help but
think there must be a way to specify the logic once per condition.

Thanks,
Bob

mydata$score1-numeric(mydata$q1) #just initializing.
mydata$score2-numeric(mydata$q1)
mydata$score1-NA
mydata$score2-NA
mydata

mydata$score1[mydata$gender == f]-  2*mydata$q1[mydata$gender==f] +

  mydata$q2[mydata$gender==f]
mydata$score2[mydata$gender == f]-2.5*mydata$q1[mydata$gender==f] +

  mydata$q2[mydata$gender==f]
mydata$score1[mydata$gender == m]-3*mydata$q1[mydata$gender==m] + 
  mydata$q2[mydata$gender==m]
mydata$score2[mydata$gender == m]-3.5*mydata$q1[mydata$gender==m] +

  mydata$q2[mydata$gender==m]
mydata

=
Bob Muenchen (pronounced Min'-chen), Manager 
Statistical Consulting Center
U of TN Office of Information Technology
200 Stokely Management Center, Knoxville, TN 37996-0520
Voice: (865) 974-5230 
FAX: (865) 974-4810
Email: [EMAIL PROTECTED]
Web: http://oit.utk.edu/scc, 
News: http://listserv.utk.edu/archives/statnews.html
=


-Original Message-
From: Leeds, Mark (IED) [mailto:[EMAIL PROTECTED] 
Sent: Friday, November 24, 2006 8:45 PM
To: Muenchen, Robert A (Bob)
Subject: RE: [R] Multiple Conditional Tranformations

I'm not sure if I understand your question but I don't think you need
iflelse statements.

myscore-numeric(q1) ( because I'm not sure how to initialize a list so
initialize a vector with q1 elements )

myscore-NA ( I think this should set all the values in myscore to NA )
myscore[mydata$gender == f]-2*mydata$q1 + mydata$q2
myscore[mydata$gender == m]-3*mydata$q1 + mydata$q2

the above should do what you do in the first part of your code but I
don't know if that was your question ?
also, it does it making myscore a vector because I didn't know how to
initialize a list.
Someone else may goive a better solution. I'm no expert.


-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Muenchen, Robert
A (Bob)
Sent: Friday, November 24, 2006 8:27 PM
To: r-help@stat.math.ethz.ch
Subject: [R] Multiple Conditional Tranformations

Greetings,

 

I'm learning R and I'm stuck on a basic concept: how to specify a
logical condition once and then perform multiple transformations under
that condition. The program below is simplified to demonstrate the goal.
Its results are exactly what I want, but I would like to check the
logical state of gender only once and create both (or any number of)
scores at once.

 

mystring-

(id,group,gender,q1,q2,q3,q4

01,1,f,2,2,5,4

02,2,f,2,1,4,5

03,1,f,2,2,4,4

04,2,f,1,1,5,5

05,1,m,4,5,4,

06,2,m,5,4,5,5

07,1,m,3,3,4,5

08,2,m,5,5,5,4)

 

mydata-read.table(textConnection(mystring),header=TRUE,sep=,,row.name
s=id)

mydata 

 

#Create score1 so that it differs for males and females:

mydata$score1 - ifelse( mydata$gender==f , 

   (mydata$score1 - (2*mydata$q1)+mydata$q2),

   ifelse( mydata$gender==m,

  (mydata$score1 - (3*mydata$q1)+mydata$q2), NA )

   )

mydata

 

#Create score2 so that it too differs for males and females:

mydata$score2 - ifelse( mydata$gender==f , 

   (mydata$score2 - (2.5*mydata$q1)+mydata$q2),

   ifelse( mydata$gender==m,

  (mydata$score2 - (3.5*mydata$q1)+mydata$q2), NA )

   )

mydata

 

 

Thanks!

Bob

=
Bob Muenchen (pronounced Min'-chen), Manager Statistical Consulting
Center U of TN Office of Information Technology 200 Stokely Management
Center, Knoxville, TN 37996-0520
Voice: (865) 974-5230
FAX: (865) 974-4810
Email: [EMAIL PROTECTED]
Web: http://oit.utk.edu/scc http://oit.utk.edu/scc ,
News: http://listserv.utk.edu/archives/statnews.html
http://listserv.utk.edu/archives/statnews.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
and provide commented, minimal, self-contained, reproducible code.


This is not an offer (or solicitation of an offer) to buy/sell the
securities/instruments mentioned or an official confirmation.  Morgan
Stanley may deal as principal in or own or act as market maker for
securities/instruments mentioned or may advise the issuers.  This is not
research and is not from MS Research but it may refer to a research
analyst/research report.  Unless indicated, these views are the author's
and may differ from those of Morgan Stanley research or others in the
Firm.  We do not represent this is accurate or complete and we may not
update this.  Past performance is not indicative of future returns.  For
additional information, research 

Re: [R] Multiple Conditional Tranformations

2006-11-24 Thread Muenchen, Robert A (Bob)
Good idea. I'm still getting used to how flexible R is on substitutions
like that! -Bob


-Original Message-
From: Leeds, Mark (IED) [mailto:[EMAIL PROTECTED] 
Sent: Friday, November 24, 2006 10:20 PM
To: Muenchen, Robert A (Bob)
Subject: RE: [R] Multiple Conditional Tranformations

You could set temp-which(my$gender[my$gender == f]) and then temp
will have the female indices and
Then you could just put temp everywhere instead of the statement but I
think that's the best you can do.
Definitely, someone will reply and there may be a shorter way that I am
unaware of.



-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Muenchen, Robert
A (Bob)
Sent: Friday, November 24, 2006 10:09 PM
To: r-help@stat.math.ethz.ch
Subject: Re: [R] Multiple Conditional Tranformations

Mark,

I finally got that approach to work by spreading the logical condition
everywhere. That gets the lengths to match. Still, I can't help but
think there must be a way to specify the logic once per condition.

Thanks,
Bob

mydata$score1-numeric(mydata$q1) #just initializing.
mydata$score2-numeric(mydata$q1)
mydata$score1-NA
mydata$score2-NA
mydata

mydata$score1[mydata$gender == f]-  2*mydata$q1[mydata$gender==f] +

  mydata$q2[mydata$gender==f]
mydata$score2[mydata$gender == f]-2.5*mydata$q1[mydata$gender==f] +

  mydata$q2[mydata$gender==f]
mydata$score1[mydata$gender == m]-3*mydata$q1[mydata$gender==m] +
  mydata$q2[mydata$gender==m]
mydata$score2[mydata$gender == m]-3.5*mydata$q1[mydata$gender==m] +

  mydata$q2[mydata$gender==m]
mydata

=
Bob Muenchen (pronounced Min'-chen), Manager Statistical Consulting
Center U of TN Office of Information Technology 200 Stokely Management
Center, Knoxville, TN 37996-0520
Voice: (865) 974-5230
FAX: (865) 974-4810
Email: [EMAIL PROTECTED]
Web: http://oit.utk.edu/scc,
News: http://listserv.utk.edu/archives/statnews.html
=


-Original Message-
From: Leeds, Mark (IED) [mailto:[EMAIL PROTECTED]
Sent: Friday, November 24, 2006 8:45 PM
To: Muenchen, Robert A (Bob)
Subject: RE: [R] Multiple Conditional Tranformations

I'm not sure if I understand your question but I don't think you need
iflelse statements.

myscore-numeric(q1) ( because I'm not sure how to initialize a list so
initialize a vector with q1 elements )

myscore-NA ( I think this should set all the values in myscore to NA )
myscore[mydata$gender == f]-2*mydata$q1 + mydata$q2
myscore[mydata$gender == m]-3*mydata$q1 + mydata$q2

the above should do what you do in the first part of your code but I
don't know if that was your question ?
also, it does it making myscore a vector because I didn't know how to
initialize a list.
Someone else may goive a better solution. I'm no expert.


-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of Muenchen, Robert
A (Bob)
Sent: Friday, November 24, 2006 8:27 PM
To: r-help@stat.math.ethz.ch
Subject: [R] Multiple Conditional Tranformations

Greetings,

 

I'm learning R and I'm stuck on a basic concept: how to specify a
logical condition once and then perform multiple transformations under
that condition. The program below is simplified to demonstrate the goal.
Its results are exactly what I want, but I would like to check the
logical state of gender only once and create both (or any number of)
scores at once.

 

mystring-

(id,group,gender,q1,q2,q3,q4

01,1,f,2,2,5,4

02,2,f,2,1,4,5

03,1,f,2,2,4,4

04,2,f,1,1,5,5

05,1,m,4,5,4,

06,2,m,5,4,5,5

07,1,m,3,3,4,5

08,2,m,5,5,5,4)

 

mydata-read.table(textConnection(mystring),header=TRUE,sep=,,row.name
s=id)

mydata 

 

#Create score1 so that it differs for males and females:

mydata$score1 - ifelse( mydata$gender==f , 

   (mydata$score1 - (2*mydata$q1)+mydata$q2),

   ifelse( mydata$gender==m,

  (mydata$score1 - (3*mydata$q1)+mydata$q2), NA )

   )

mydata

 

#Create score2 so that it too differs for males and females:

mydata$score2 - ifelse( mydata$gender==f , 

   (mydata$score2 - (2.5*mydata$q1)+mydata$q2),

   ifelse( mydata$gender==m,

  (mydata$score2 - (3.5*mydata$q1)+mydata$q2), NA )

   )

mydata

 

 

Thanks!

Bob

=
Bob Muenchen (pronounced Min'-chen), Manager Statistical Consulting
Center U of TN Office of Information Technology 200 Stokely Management
Center, Knoxville, TN 37996-0520
Voice: (865) 974-5230
FAX: (865) 974-4810
Email: [EMAIL PROTECTED]
Web: http://oit.utk.edu/scc http://oit.utk.edu/scc ,
News: http://listserv.utk.edu/archives/statnews.html
http://listserv.utk.edu/archives/statnews.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

[R] predict and arima

2006-11-24 Thread Franck Arnaud
Hi all,

Forecasting from an arima model is easy with predict.

But I can't manage to backcast : invent data from the model before the
begining of the sample.
The theory is easy : take your parameters, reverse your data, forecast, and
then reverse the forecast
I've tried to adapt the predict function to do that (i'm not sure that the
statistical procedure is fine (with the residuals), but that's not my point
right now) :

mav.backcast.arima-function(model,n.backcast,...)
{
if (class(model)[1]!=Arima) stop(argument model must be an object
of class 'Arima' (see ?arima))

model2-model
model$residuals-rev(model$residuals)
if (is.ts(model2$residuals))
model$residuals-ts(model$residuals,start=start(model2$residuals),
frequency=frequency(model2$residuals))
pred.before-predict(model,n.ahead=n.backcast,...)

freq-frequency(model$residuals)
startingdate-per.sub(start(model2$residuals),n.backcast,freq=freq)

pred-ts(rev(pred.before$pred),start=startingdate,freq=freq)
se-ts(rev(pred.before$se),start=startingdate,freq=freq)

return((list(pred = pred, se =se)))
}

This function does not work : it gives always the same result, it does not
depend on the residuals (i've tried to insert
a model$residuals-rep(1,100) after the definition, to check that).

Then i look at the code, with getS3method(predict,Arima). And i get even
more confused (!) :
where does data play a role in the function ? residuals are loaded into rsd,
but this variable is not used after...
I looked at KalmanForecast and at the C code of KalmanFore, but it did not
help me understand what was going on.

thanks

Franck A.

btw, it has nothing to do with it, but i've done some stuff on time series
(filtering with Hodrick prescott or Baxter King, for instance) that you can
find on http://arnaud.ensae.net

[[alternative HTML version deleted]]

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


Re: [R] Multiple Conditional Tranformations

2006-11-24 Thread Gabor Grothendieck
Try this:


transform(mydata,
   score1 = (2   + (gender == m)) * q1 + q2,
   score2 = (2.5 + (gender == m)) * q1 + q2
)


On 11/24/06, Muenchen, Robert A (Bob) [EMAIL PROTECTED] wrote:
 Mark,

 I finally got that approach to work by spreading the logical condition
 everywhere. That gets the lengths to match. Still, I can't help but
 think there must be a way to specify the logic once per condition.

 Thanks,
 Bob

 mydata$score1-numeric(mydata$q1) #just initializing.
 mydata$score2-numeric(mydata$q1)
 mydata$score1-NA
 mydata$score2-NA
 mydata

 mydata$score1[mydata$gender == f]-  2*mydata$q1[mydata$gender==f] +

  mydata$q2[mydata$gender==f]
 mydata$score2[mydata$gender == f]-2.5*mydata$q1[mydata$gender==f] +

  mydata$q2[mydata$gender==f]
 mydata$score1[mydata$gender == m]-3*mydata$q1[mydata$gender==m] +
  mydata$q2[mydata$gender==m]
 mydata$score2[mydata$gender == m]-3.5*mydata$q1[mydata$gender==m] +

  mydata$q2[mydata$gender==m]
 mydata

 =
 Bob Muenchen (pronounced Min'-chen), Manager
 Statistical Consulting Center
 U of TN Office of Information Technology
 200 Stokely Management Center, Knoxville, TN 37996-0520
 Voice: (865) 974-5230
 FAX: (865) 974-4810
 Email: [EMAIL PROTECTED]
 Web: http://oit.utk.edu/scc,
 News: http://listserv.utk.edu/archives/statnews.html
 =


 -Original Message-
 From: Leeds, Mark (IED) [mailto:[EMAIL PROTECTED]
 Sent: Friday, November 24, 2006 8:45 PM
 To: Muenchen, Robert A (Bob)
 Subject: RE: [R] Multiple Conditional Tranformations

 I'm not sure if I understand your question but I don't think you need
 iflelse statements.

 myscore-numeric(q1) ( because I'm not sure how to initialize a list so
 initialize a vector with q1 elements )

 myscore-NA ( I think this should set all the values in myscore to NA )
 myscore[mydata$gender == f]-2*mydata$q1 + mydata$q2
 myscore[mydata$gender == m]-3*mydata$q1 + mydata$q2

 the above should do what you do in the first part of your code but I
 don't know if that was your question ?
 also, it does it making myscore a vector because I didn't know how to
 initialize a list.
 Someone else may goive a better solution. I'm no expert.


 -Original Message-
 From: [EMAIL PROTECTED]
 [mailto:[EMAIL PROTECTED] On Behalf Of Muenchen, Robert
 A (Bob)
 Sent: Friday, November 24, 2006 8:27 PM
 To: r-help@stat.math.ethz.ch
 Subject: [R] Multiple Conditional Tranformations

 Greetings,



 I'm learning R and I'm stuck on a basic concept: how to specify a
 logical condition once and then perform multiple transformations under
 that condition. The program below is simplified to demonstrate the goal.
 Its results are exactly what I want, but I would like to check the
 logical state of gender only once and create both (or any number of)
 scores at once.



 mystring-

 (id,group,gender,q1,q2,q3,q4

 01,1,f,2,2,5,4

 02,2,f,2,1,4,5

 03,1,f,2,2,4,4

 04,2,f,1,1,5,5

 05,1,m,4,5,4,

 06,2,m,5,4,5,5

 07,1,m,3,3,4,5

 08,2,m,5,5,5,4)



 mydata-read.table(textConnection(mystring),header=TRUE,sep=,,row.name
 s=id)

 mydata



 #Create score1 so that it differs for males and females:

 mydata$score1 - ifelse( mydata$gender==f ,

   (mydata$score1 - (2*mydata$q1)+mydata$q2),

   ifelse( mydata$gender==m,

  (mydata$score1 - (3*mydata$q1)+mydata$q2), NA )

   )

 mydata



 #Create score2 so that it too differs for males and females:

 mydata$score2 - ifelse( mydata$gender==f ,

   (mydata$score2 - (2.5*mydata$q1)+mydata$q2),

   ifelse( mydata$gender==m,

  (mydata$score2 - (3.5*mydata$q1)+mydata$q2), NA )

   )

 mydata





 Thanks!

 Bob

 =
 Bob Muenchen (pronounced Min'-chen), Manager Statistical Consulting
 Center U of TN Office of Information Technology 200 Stokely Management
 Center, Knoxville, TN 37996-0520
 Voice: (865) 974-5230
 FAX: (865) 974-4810
 Email: [EMAIL PROTECTED]
 Web: http://oit.utk.edu/scc http://oit.utk.edu/scc ,
 News: http://listserv.utk.edu/archives/statnews.html
 http://listserv.utk.edu/archives/statnews.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
 and provide commented, minimal, self-contained, reproducible code.
 

 This is not an offer (or solicitation of an offer) to buy/sell the
 securities/instruments mentioned or an official confirmation.  Morgan
 Stanley may deal as principal in or own or act as market maker for
 securities/instruments mentioned or may advise the issuers.  This is not
 research and is not from MS Research but it may refer to a research
 analyst/research report.  Unless indicated, these 

Re: [R] Multiple Conditional Tranformations

2006-11-24 Thread Gabor Grothendieck
And here is a variation:

transform(mydata,
   score1 = (2 + (gender == m)) * q1 + q2,
   score2 = score1 + 0.5 * q1
)

or

transform(
   transform(mydata, score1 = (2 + (gender == m)) * q1 + q2),
   score2 = score1 + 0.5 * q1
)


On 11/25/06, Gabor Grothendieck [EMAIL PROTECTED] wrote:
 Try this:


 transform(mydata,
   score1 = (2   + (gender == m)) * q1 + q2,
   score2 = (2.5 + (gender == m)) * q1 + q2
 )


 On 11/24/06, Muenchen, Robert A (Bob) [EMAIL PROTECTED] wrote:
  Mark,
 
  I finally got that approach to work by spreading the logical condition
  everywhere. That gets the lengths to match. Still, I can't help but
  think there must be a way to specify the logic once per condition.
 
  Thanks,
  Bob
 
  mydata$score1-numeric(mydata$q1) #just initializing.
  mydata$score2-numeric(mydata$q1)
  mydata$score1-NA
  mydata$score2-NA
  mydata
 
  mydata$score1[mydata$gender == f]-  2*mydata$q1[mydata$gender==f] +
 
   mydata$q2[mydata$gender==f]
  mydata$score2[mydata$gender == f]-2.5*mydata$q1[mydata$gender==f] +
 
   mydata$q2[mydata$gender==f]
  mydata$score1[mydata$gender == m]-3*mydata$q1[mydata$gender==m] +
   mydata$q2[mydata$gender==m]
  mydata$score2[mydata$gender == m]-3.5*mydata$q1[mydata$gender==m] +
 
   mydata$q2[mydata$gender==m]
  mydata
 
  =
  Bob Muenchen (pronounced Min'-chen), Manager
  Statistical Consulting Center
  U of TN Office of Information Technology
  200 Stokely Management Center, Knoxville, TN 37996-0520
  Voice: (865) 974-5230
  FAX: (865) 974-4810
  Email: [EMAIL PROTECTED]
  Web: http://oit.utk.edu/scc,
  News: http://listserv.utk.edu/archives/statnews.html
  =
 
 
  -Original Message-
  From: Leeds, Mark (IED) [mailto:[EMAIL PROTECTED]
  Sent: Friday, November 24, 2006 8:45 PM
  To: Muenchen, Robert A (Bob)
  Subject: RE: [R] Multiple Conditional Tranformations
 
  I'm not sure if I understand your question but I don't think you need
  iflelse statements.
 
  myscore-numeric(q1) ( because I'm not sure how to initialize a list so
  initialize a vector with q1 elements )
 
  myscore-NA ( I think this should set all the values in myscore to NA )
  myscore[mydata$gender == f]-2*mydata$q1 + mydata$q2
  myscore[mydata$gender == m]-3*mydata$q1 + mydata$q2
 
  the above should do what you do in the first part of your code but I
  don't know if that was your question ?
  also, it does it making myscore a vector because I didn't know how to
  initialize a list.
  Someone else may goive a better solution. I'm no expert.
 
 
  -Original Message-
  From: [EMAIL PROTECTED]
  [mailto:[EMAIL PROTECTED] On Behalf Of Muenchen, Robert
  A (Bob)
  Sent: Friday, November 24, 2006 8:27 PM
  To: r-help@stat.math.ethz.ch
  Subject: [R] Multiple Conditional Tranformations
 
  Greetings,
 
 
 
  I'm learning R and I'm stuck on a basic concept: how to specify a
  logical condition once and then perform multiple transformations under
  that condition. The program below is simplified to demonstrate the goal.
  Its results are exactly what I want, but I would like to check the
  logical state of gender only once and create both (or any number of)
  scores at once.
 
 
 
  mystring-
 
  (id,group,gender,q1,q2,q3,q4
 
  01,1,f,2,2,5,4
 
  02,2,f,2,1,4,5
 
  03,1,f,2,2,4,4
 
  04,2,f,1,1,5,5
 
  05,1,m,4,5,4,
 
  06,2,m,5,4,5,5
 
  07,1,m,3,3,4,5
 
  08,2,m,5,5,5,4)
 
 
 
  mydata-read.table(textConnection(mystring),header=TRUE,sep=,,row.name
  s=id)
 
  mydata
 
 
 
  #Create score1 so that it differs for males and females:
 
  mydata$score1 - ifelse( mydata$gender==f ,
 
(mydata$score1 - (2*mydata$q1)+mydata$q2),
 
ifelse( mydata$gender==m,
 
   (mydata$score1 - (3*mydata$q1)+mydata$q2), NA )
 
)
 
  mydata
 
 
 
  #Create score2 so that it too differs for males and females:
 
  mydata$score2 - ifelse( mydata$gender==f ,
 
(mydata$score2 - (2.5*mydata$q1)+mydata$q2),
 
ifelse( mydata$gender==m,
 
   (mydata$score2 - (3.5*mydata$q1)+mydata$q2), NA )
 
)
 
  mydata
 
 
 
 
 
  Thanks!
 
  Bob
 
  =
  Bob Muenchen (pronounced Min'-chen), Manager Statistical Consulting
  Center U of TN Office of Information Technology 200 Stokely Management
  Center, Knoxville, TN 37996-0520
  Voice: (865) 974-5230
  FAX: (865) 974-4810
  Email: [EMAIL PROTECTED]
  Web: http://oit.utk.edu/scc http://oit.utk.edu/scc ,
  News: http://listserv.utk.edu/archives/statnews.html
  http://listserv.utk.edu/archives/statnews.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
  and provide commented, minimal, self-contained, reproducible