Re: [R] Help interpreting density().

2008-07-29 Thread Bill.Venables
You should read the documentation more carefully.  The bw is not
essentially the sd.  To quote the documentation the bw is the
smoothing bandwidth to be used. The kernels are scaled such that this is
the standard deviation of the smoothing kernel.  That is a very
different thing.

You are confusing the standard deviation of the distribution with the
standard deviation of the gaussian smoothing kernels.  

In the second case, density(rpois(1000, 0)), you are getting the kernel
density for a sample of 1000 zeros.  So there is just one distinct
smoothing kernel and the bw is a default used for this case.  If you 

plot(density(rpois(1000, 0)))

you will see what that smoothing kernel looks like.


Bill Venables
CSIRO Laboratories
PO Box 120, Cleveland, 4163
AUSTRALIA
Office Phone (email preferred): +61 7 3826 7251
Fax (if absolutely necessary):  +61 7 3826 7304
Mobile: +61 4 8819 4402
Home Phone: +61 7 3286 7700
mailto:[EMAIL PROTECTED]
http://www.cmis.csiro.au/bill.venables/ 

-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED]
On Behalf Of [EMAIL PROTECTED]
Sent: Tuesday, 29 July 2008 2:15 PM
To: r-help@r-project.org
Subject: [R] Help interpreting density().

I issue the following:

 d - density(rnorm(1000))
 d

and get:

Call:
density.default(x = rnorm(1000))

Data: rnorm(1000) (1000 obs.);  Bandwidth 'bw' = 0.2235

   x y
 Min.   :-3.5157   Min.   :2.416e-05  
 1st Qu.:-1.6892   1st Qu.:1.129e-02  
 Median : 0.1373   Median :7.267e-02  
 Mean   : 0.1373   Mean   :1.367e-01  
 3rd Qu.: 1.9639   3rd Qu.:2.693e-01  
 Max.   : 3.7904   Max.   :4.014e-01  

The documentation indicates that the bw is essentially the sd. Yet I
have specified an sd of 1? How am I to interpret the ranges of the
values? x ranges almost from -4 to +4 and y ranges from 0 to 0.4. The
mean x is .1 which isn't too awfully close to what I would expect (0.0).
Then there is:

 d - density(rpois(1000,0))
 d

Call:
density.default(x = rpois(1000, 0))

Data: rpois(1000, 0) (1000 obs.);   Bandwidth 'bw' = 0.2261

   x y  
 Min.   :-0.6782   Min.   :0.01979  
 1st Qu.:-0.3391   1st Qu.:0.14073  
 Median : 0.   Median :0.57178  
 Mean   : 0.   Mean   :0.73454  
 3rd Qu.: 0.3391   3rd Qu.:1.32830  
 Max.   : 0.6782   Max.   :1.76436  

Here I am getting the mean that I expect from a Poisson distribuition
but y ranges from 0 to 1.75. Again I am not sure what these numbers
mean. How can I map the output to the standard distirbution description
parameters?

Thank you.

Kevin

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

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


Re: [R] Help interpreting density().

2008-07-29 Thread Mark Difford

Hi Kevin,

 The documentation indicates that the bw is essentially the sd.
  d - density(rnorm(1000))

Not so. The documentation states that the following about bw: The kernels
are scaled such that this is the standard deviation of the smoothing
kernel..., which is a very different thing.

The default bandwidth used by density is ?bw.nrd0. Read that documentation
carefully and all might be clear.

HTH, Mark.


rkevinburton wrote:
 
 I issue the following:
 
 d - density(rnorm(1000))
 d
 
 and get:
 
 Call:
 density.default(x = rnorm(1000))
 
 Data: rnorm(1000) (1000 obs.);  Bandwidth 'bw' = 0.2235
 
x y
  Min.   :-3.5157   Min.   :2.416e-05  
  1st Qu.:-1.6892   1st Qu.:1.129e-02  
  Median : 0.1373   Median :7.267e-02  
  Mean   : 0.1373   Mean   :1.367e-01  
  3rd Qu.: 1.9639   3rd Qu.:2.693e-01  
  Max.   : 3.7904   Max.   :4.014e-01  
 
 The documentation indicates that the bw is essentially the sd. Yet I have
 specified an sd of 1? How am I to interpret the ranges of the values? x
 ranges almost from -4 to +4 and y ranges from 0 to 0.4. The mean x is .1
 which isn't too awfully close to what I would expect (0.0). Then there is:
 
 d - density(rpois(1000,0))
 d
 
 Call:
 density.default(x = rpois(1000, 0))
 
 Data: rpois(1000, 0) (1000 obs.);   Bandwidth 'bw' = 0.2261
 
x y  
  Min.   :-0.6782   Min.   :0.01979  
  1st Qu.:-0.3391   1st Qu.:0.14073  
  Median : 0.   Median :0.57178  
  Mean   : 0.   Mean   :0.73454  
  3rd Qu.: 0.3391   3rd Qu.:1.32830  
  Max.   : 0.6782   Max.   :1.76436  
 
 Here I am getting the mean that I expect from a Poisson distribuition but
 y ranges from 0 to 1.75. Again I am not sure what these numbers mean. How
 can I map the output to the standard distirbution description parameters?
 
 Thank you.
 
 Kevin
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 

-- 
View this message in context: 
http://www.nabble.com/Help-interpreting-density%28%29.-tp18704955p18706154.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Rf_error crashes entire program.

2008-07-29 Thread Prof Brian Ripley
Rf_error is used so often in R itself and packages that it is almost 
certainly not the problem -- rather something else your program has done 
has damaged R's internals (e.g. overrrun an array).


'Writing R Extensions' discusses how to debug R code, including foreign 
code.  For example, I would run under valgrind to look for incorrect 
memory accesses.


This was an appropriate topic for R-devel -- please DO study the posting 
guide.


On Mon, 28 Jul 2008, Andrew Redd wrote:


I'm having a problem with the error and warning functions.  I've tried this
on multiple machine so I'm fairly sure it's not machine dependent and I've
tried it on the latest versions 2.6.0-2.7.1.  Whenever my program gets to an
error or warning it crashes the entire program rather than throwing the
error like it should.

Here are the relevant bits.

...
   char const * const ExeedsMinVarianceError = PFDA ERR: Near zero
variance encountered.  Estimation Unstable. Terminating Estimation.;

   if(debug){printf(Da:\n);printmat(DaOld,1,*ka);fflush(stdout);}
   daxpy_(ka, mOne, Da, one, DaOld, one);
   for(i=0;i*ka;i++)convergenceCriteria+=fabs(DaOld[i]);
   if(Da[*ka]  MinVariance){
   printf(PING);fflush(stdout);
   warning(ExeedsMinVarianceError);
   break;
   }

and the output from running in batch mode that I get is this:
Da:
0.18034.988e-017

PING
and here the program crashes.   I've tried this in multiple places and
sometimes the error is thrown sometimes not.  Does anyone have an idea of
what is going on.  I could not find any discussion of this error yet.
Thank you,
Andrew Redd

[[alternative HTML version deleted]]

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



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

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


[R] How to set the parameters in Trellis Graphics (by Lattice package)

2008-07-29 Thread G.H. Zuo
Dear R users

I plot the trellis graphics by using the lattice package.  Everything is
OK.  Now, I want set
some parameters of the trellis graphics.

1. The tick label site.  By default, only two tick labels had been output in
x-axis of my plot.
I want output four or five tick labels.  In the traditional graphics system,
it will be very easy.
Just not output the axis in plot by use the parameter xaxt=n, and then
used axis to add
the axis.  But it seem useless in the grid graphics system.  Then how can I
do it?

2. How to change the margin of the figure?

3. Can I change the style of the figure like par in traditional system.  For
example, I want
change all the fontface into bold.  Now, I must change the fontface for
every item by using
command trellis.par.set() like this:

trellis.par.set(list(par.xlab.text = gpar(font = 2),
 par.ylab.text = gpar(font = 2),
 axis.text = gpar(font = 2),
 add.text = gpar(font=2)
   ))

is there a method to change this parameter in one time, like par(font=2)

thanks

G.H.Zuo

[[alternative HTML version deleted]]

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


Re: [R] product of successive rows

2008-07-29 Thread Moshe Olshansky
Assuming that the number of rows is even and that your matrix is A, 
element-wise product of pairs of rows can be calculated as

A[seq(1,nrow(A),by=2),]*A{seq(2,nrow(A),by=2),]


--- On Mon, 28/7/08, rcoder [EMAIL PROTECTED] wrote:

 From: rcoder [EMAIL PROTECTED]
 Subject: [R]  product of successive rows
 To: r-help@r-project.org
 Received: Monday, 28 July, 2008, 8:20 AM
 Hi everyone,
 
 I want to perform an operation on a matrx that outputs the
 product of
 successive pairs of rows. For example: calculating the
 product between rows
 1  2; 3  4; 5  6...etc.
 
 Does anyone know of any readily available functions that
 can do this?
 
 Thanks,
 
 rcoder
 
 
 -- 
 View this message in context:
 http://www.nabble.com/product-of-successive-rows-tp18681259p18681259.html
 Sent from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained,
 reproducible code.

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


Re: [R] Help interpreting density().

2008-07-29 Thread rkevinburton
OK. Thank you for pointing out my mistake.

I still have my original question. How does the output relate to estimating the 
parameters of a given density? I read that for a gausian kernal:

bw.nrd0 implements a rule-of-thumb for choosing the bandwidth of a Gaussian 
kernel density estimator. It defaults to 0.9 times the minimum of the standard 
deviation and the interquartile range divided by 1.34 times the sample size to 
the negative one-fifth power (= Silverman's ‘rule of thumb’

But how does that relate to say a Poisson distribution or a two-parameter 
distribution like a normal, beta, or binomial distribution?

Thank you.

Kevin

 Mark Difford [EMAIL PROTECTED] wrote: 
 
 Hi Kevin,
 
  The documentation indicates that the bw is essentially the sd.
   d - density(rnorm(1000))
 
 Not so. The documentation states that the following about bw: The kernels
 are scaled such that this is the standard deviation of the smoothing
 kernel..., which is a very different thing.
 
 The default bandwidth used by density is ?bw.nrd0. Read that documentation
 carefully and all might be clear.
 
 HTH, Mark.
 
 
 rkevinburton wrote:
  
  I issue the following:
  
  d - density(rnorm(1000))
  d
  
  and get:
  
  Call:
  density.default(x = rnorm(1000))
  
  Data: rnorm(1000) (1000 obs.);  Bandwidth 'bw' = 0.2235
  
 x y
   Min.   :-3.5157   Min.   :2.416e-05  
   1st Qu.:-1.6892   1st Qu.:1.129e-02  
   Median : 0.1373   Median :7.267e-02  
   Mean   : 0.1373   Mean   :1.367e-01  
   3rd Qu.: 1.9639   3rd Qu.:2.693e-01  
   Max.   : 3.7904   Max.   :4.014e-01  
  
  The documentation indicates that the bw is essentially the sd. Yet I have
  specified an sd of 1? How am I to interpret the ranges of the values? x
  ranges almost from -4 to +4 and y ranges from 0 to 0.4. The mean x is .1
  which isn't too awfully close to what I would expect (0.0). Then there is:
  
  d - density(rpois(1000,0))
  d
  
  Call:
  density.default(x = rpois(1000, 0))
  
  Data: rpois(1000, 0) (1000 obs.);   Bandwidth 'bw' = 0.2261
  
 x y  
   Min.   :-0.6782   Min.   :0.01979  
   1st Qu.:-0.3391   1st Qu.:0.14073  
   Median : 0.   Median :0.57178  
   Mean   : 0.   Mean   :0.73454  
   3rd Qu.: 0.3391   3rd Qu.:1.32830  
   Max.   : 0.6782   Max.   :1.76436  
  
  Here I am getting the mean that I expect from a Poisson distribuition but
  y ranges from 0 to 1.75. Again I am not sure what these numbers mean. How
  can I map the output to the standard distirbution description parameters?
  
  Thank you.
  
  Kevin
  
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
  http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
  
  
 
 -- 
 View this message in context: 
 http://www.nabble.com/Help-interpreting-density%28%29.-tp18704955p18706154.html
 Sent from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] Help interpreting density().

2008-07-29 Thread rkevinburton
Sorry, poor example. I started with normal deviates and jumped without thinking 
to Poisson. The main crux of the question is how does the output of density 
relate to the parameters that describe some of the standard distributions (mean 
and std for normal, lambda for Poisson, n and p for Binomial, alpha and beta 
for Beta, etc.).

Thank you.

Kevin

 [EMAIL PROTECTED] wrote: 
 You should read the documentation more carefully.  The bw is not
 essentially the sd.  To quote the documentation the bw is the
 smoothing bandwidth to be used. The kernels are scaled such that this is
 the standard deviation of the smoothing kernel.  That is a very
 different thing.
 
 You are confusing the standard deviation of the distribution with the
 standard deviation of the gaussian smoothing kernels.  
 
 In the second case, density(rpois(1000, 0)), you are getting the kernel
 density for a sample of 1000 zeros.  So there is just one distinct
 smoothing kernel and the bw is a default used for this case.  If you 
 
 plot(density(rpois(1000, 0)))
 
 you will see what that smoothing kernel looks like.
 
 
 Bill Venables
 CSIRO Laboratories
 PO Box 120, Cleveland, 4163
 AUSTRALIA
 Office Phone (email preferred): +61 7 3826 7251
 Fax (if absolutely necessary):  +61 7 3826 7304
 Mobile: +61 4 8819 4402
 Home Phone: +61 7 3286 7700
 mailto:[EMAIL PROTECTED]
 http://www.cmis.csiro.au/bill.venables/ 
 
 -Original Message-
 From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED]
 On Behalf Of [EMAIL PROTECTED]
 Sent: Tuesday, 29 July 2008 2:15 PM
 To: r-help@r-project.org
 Subject: [R] Help interpreting density().
 
 I issue the following:
 
  d - density(rnorm(1000))
  d
 
 and get:
 
 Call:
 density.default(x = rnorm(1000))
 
 Data: rnorm(1000) (1000 obs.);  Bandwidth 'bw' = 0.2235
 
x y
  Min.   :-3.5157   Min.   :2.416e-05  
  1st Qu.:-1.6892   1st Qu.:1.129e-02  
  Median : 0.1373   Median :7.267e-02  
  Mean   : 0.1373   Mean   :1.367e-01  
  3rd Qu.: 1.9639   3rd Qu.:2.693e-01  
  Max.   : 3.7904   Max.   :4.014e-01  
 
 The documentation indicates that the bw is essentially the sd. Yet I
 have specified an sd of 1? How am I to interpret the ranges of the
 values? x ranges almost from -4 to +4 and y ranges from 0 to 0.4. The
 mean x is .1 which isn't too awfully close to what I would expect (0.0).
 Then there is:
 
  d - density(rpois(1000,0))
  d
 
 Call:
 density.default(x = rpois(1000, 0))
 
 Data: rpois(1000, 0) (1000 obs.);   Bandwidth 'bw' = 0.2261
 
x y  
  Min.   :-0.6782   Min.   :0.01979  
  1st Qu.:-0.3391   1st Qu.:0.14073  
  Median : 0.   Median :0.57178  
  Mean   : 0.   Mean   :0.73454  
  3rd Qu.: 0.3391   3rd Qu.:1.32830  
  Max.   : 0.6782   Max.   :1.76436  
 
 Here I am getting the mean that I expect from a Poisson distribuition
 but y ranges from 0 to 1.75. Again I am not sure what these numbers
 mean. How can I map the output to the standard distirbution description
 parameters?
 
 Thank you.
 
 Kevin
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 


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


Re: [R] Help interpreting density().

2008-07-29 Thread Mark Difford

Hi Kevin,

 I still have my original question. How does the output relate to
 estimating the parameters 
 of a given density? I read that for a gausian kernal:

This isn't the place for such questions: you need to do some _basic_ reading
on the subject so that you begin to understand something about the method
you are messing about with. Basically (very basically) it's a smoothed out
histogram. And you will probably (?) know that a histogram is [still used]
to show you how a set of univariate data (random variable) is distributed.

Perhaps start with http://en.wikipedia.org/wiki/Kernel_density, then go
somewhere else. But since you have access to the web you really should have
found something like this yourself.

Regards, Mark.


rkevinburton wrote:
 
 OK. Thank you for pointing out my mistake.
 
 I still have my original question. How does the output relate to
 estimating the parameters of a given density? I read that for a gausian
 kernal:
 
 bw.nrd0 implements a rule-of-thumb for choosing the bandwidth of a
 Gaussian kernel density estimator. It defaults to 0.9 times the minimum of
 the standard deviation and the interquartile range divided by 1.34 times
 the sample size to the negative one-fifth power (= Silverman's ‘rule of
 thumb’
 
 But how does that relate to say a Poisson distribution or a two-parameter
 distribution like a normal, beta, or binomial distribution?
 
 Thank you.
 
 Kevin
 
  Mark Difford [EMAIL PROTECTED] wrote: 
 
 Hi Kevin,
 
  The documentation indicates that the bw is essentially the sd.
   d - density(rnorm(1000))
 
 Not so. The documentation states that the following about bw: The
 kernels
 are scaled such that this is the standard deviation of the smoothing
 kernel..., which is a very different thing.
 
 The default bandwidth used by density is ?bw.nrd0. Read that
 documentation
 carefully and all might be clear.
 
 HTH, Mark.
 
 
 rkevinburton wrote:
  
  I issue the following:
  
  d - density(rnorm(1000))
  d
  
  and get:
  
  Call:
  density.default(x = rnorm(1000))
  
  Data: rnorm(1000) (1000 obs.);  Bandwidth 'bw' = 0.2235
  
 x y
   Min.   :-3.5157   Min.   :2.416e-05  
   1st Qu.:-1.6892   1st Qu.:1.129e-02  
   Median : 0.1373   Median :7.267e-02  
   Mean   : 0.1373   Mean   :1.367e-01  
   3rd Qu.: 1.9639   3rd Qu.:2.693e-01  
   Max.   : 3.7904   Max.   :4.014e-01  
  
  The documentation indicates that the bw is essentially the sd. Yet I
 have
  specified an sd of 1? How am I to interpret the ranges of the values? x
  ranges almost from -4 to +4 and y ranges from 0 to 0.4. The mean x is
 .1
  which isn't too awfully close to what I would expect (0.0). Then there
 is:
  
  d - density(rpois(1000,0))
  d
  
  Call:
  density.default(x = rpois(1000, 0))
  
  Data: rpois(1000, 0) (1000 obs.);   Bandwidth 'bw' = 0.2261
  
 x y  
   Min.   :-0.6782   Min.   :0.01979  
   1st Qu.:-0.3391   1st Qu.:0.14073  
   Median : 0.   Median :0.57178  
   Mean   : 0.   Mean   :0.73454  
   3rd Qu.: 0.3391   3rd Qu.:1.32830  
   Max.   : 0.6782   Max.   :1.76436  
  
  Here I am getting the mean that I expect from a Poisson distribuition
 but
  y ranges from 0 to 1.75. Again I am not sure what these numbers mean.
 How
  can I map the output to the standard distirbution description
 parameters?
  
  Thank you.
  
  Kevin
  
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
  http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
  
  
 
 -- 
 View this message in context:
 http://www.nabble.com/Help-interpreting-density%28%29.-tp18704955p18706154.html
 Sent from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 

-- 
View this message in context: 
http://www.nabble.com/Help-interpreting-density%28%29.-tp18704955p18707522.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Help interpreting density().

2008-07-29 Thread rkevinburton
Sorry I tried WikiPedia and only found:

Wikipedia does not have an article with this exact name.

I will try to find some other sources of information.

Kevin

 Mark Difford [EMAIL PROTECTED] wrote: 
 
Hi Kevin,

 I still have my original question. How does the output relate to
 estimating the parameters 
 of a given density? I read that for a gausian kernal:

This isn't the place for such questions: you need to do some _basic_ reading
on the subject so that you begin to understand something about the method
you are messing about with. Basically (very basically) it's a smoothed out
histogram. And you will probably (?) know that a histogram is [still used]
to show you how a set of univariate data (random variable) is distributed.

Perhaps start with http://en.wikipedia.org/wiki/Kernel_density, then go
somewhere else. But since you have access to the web you really should have
found something like this yourself.

Regards, Mark.


rkevinburton wrote:
 
 OK. Thank you for pointing out my mistake.
 
 I still have my original question. How does the output relate to
 estimating the parameters of a given density? I read that for a gausian
 kernal:
 
 bw.nrd0 implements a rule-of-thumb for choosing the bandwidth of a
 Gaussian kernel density estimator. It defaults to 0.9 times the minimum of
 the standard deviation and the interquartile range divided by 1.34 times
 the sample size to the negative one-fifth power (= Silverman's ‘rule of
 thumb’
 
 But how does that relate to say a Poisson distribution or a two-parameter
 distribution like a normal, beta, or binomial distribution?
 
 Thank you.
 
 Kevin
 
  Mark Difford [EMAIL PROTECTED] wrote: 
 
 Hi Kevin,
 
  The documentation indicates that the bw is essentially the sd.
   d - density(rnorm(1000))
 
 Not so. The documentation states that the following about bw: The
 kernels
 are scaled such that this is the standard deviation of the smoothing
 kernel..., which is a very different thing.
 
 The default bandwidth used by density is ?bw.nrd0. Read that
 documentation
 carefully and all might be clear.
 
 HTH, Mark.
 
 
 rkevinburton wrote:
  
  I issue the following:
  
  d - density(rnorm(1000))
  d
  
  and get:
  
  Call:
  density.default(x = rnorm(1000))
  
  Data: rnorm(1000) (1000 obs.);  Bandwidth 'bw' = 0.2235
  
 x y
   Min.   :-3.5157   Min.   :2.416e-05  
   1st Qu.:-1.6892   1st Qu.:1.129e-02  
   Median : 0.1373   Median :7.267e-02  
   Mean   : 0.1373   Mean   :1.367e-01  
   3rd Qu.: 1.9639   3rd Qu.:2.693e-01  
   Max.   : 3.7904   Max.   :4.014e-01  
  
  The documentation indicates that the bw is essentially the sd. Yet I
 have
  specified an sd of 1? How am I to interpret the ranges of the values? x
  ranges almost from -4 to +4 and y ranges from 0 to 0.4. The mean x is
 .1
  which isn't too awfully close to what I would expect (0.0). Then there
 is:
  
  d - density(rpois(1000,0))
  d
  
  Call:
  density.default(x = rpois(1000, 0))
  
  Data: rpois(1000, 0) (1000 obs.);   Bandwidth 'bw' = 0.2261
  
 x y  
   Min.   :-0.6782   Min.   :0.01979  
   1st Qu.:-0.3391   1st Qu.:0.14073  
   Median : 0.   Median :0.57178  
   Mean   : 0.   Mean   :0.73454  
   3rd Qu.: 0.3391   3rd Qu.:1.32830  
   Max.   : 0.6782   Max.   :1.76436  
  
  Here I am getting the mean that I expect from a Poisson distribuition
 but
  y ranges from 0 to 1.75. Again I am not sure what these numbers mean.
 How
  can I map the output to the standard distirbution description
 parameters?
  
  Thank you.
  
  Kevin
  
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
  http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
  
  
 
 -- 
 View this message in context:
 http://www.nabble.com/Help-interpreting-density%28%29.-tp18704955p18706154.html
 Sent from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 

-- 
View this message in context: 
http://www.nabble.com/Help-interpreting-density%28%29.-tp18704955p18707522.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] Help interpreting density().

2008-07-29 Thread Mark Difford

Hi Kevin,

Clicking on the link I sent gets me there (?), though things are pretty slow
at the moment. Perhaps try this related link, and from it get back to the
first one:

http://en.wikipedia.org/wiki/Density_estimation

You can also get to this via histogram, so search for that in Wiki, and
then...

HTH, Mark.


rkevinburton wrote:
 
 Sorry I tried WikiPedia and only found:
 
 Wikipedia does not have an article with this exact name.
 
 I will try to find some other sources of information.
 
 Kevin
 
  Mark Difford [EMAIL PROTECTED] wrote: 
 
 Hi Kevin,
 
 I still have my original question. How does the output relate to
 estimating the parameters 
 of a given density? I read that for a gausian kernal:
 
 This isn't the place for such questions: you need to do some _basic_
 reading
 on the subject so that you begin to understand something about the method
 you are messing about with. Basically (very basically) it's a smoothed out
 histogram. And you will probably (?) know that a histogram is [still used]
 to show you how a set of univariate data (random variable) is distributed.
 
 Perhaps start with http://en.wikipedia.org/wiki/Kernel_density, then go
 somewhere else. But since you have access to the web you really should
 have
 found something like this yourself.
 
 Regards, Mark.
 
 
 rkevinburton wrote:
 
 OK. Thank you for pointing out my mistake.
 
 I still have my original question. How does the output relate to
 estimating the parameters of a given density? I read that for a gausian
 kernal:
 
 bw.nrd0 implements a rule-of-thumb for choosing the bandwidth of a
 Gaussian kernel density estimator. It defaults to 0.9 times the minimum
 of
 the standard deviation and the interquartile range divided by 1.34 times
 the sample size to the negative one-fifth power (= Silverman's ‘rule of
 thumb’
 
 But how does that relate to say a Poisson distribution or a two-parameter
 distribution like a normal, beta, or binomial distribution?
 
 Thank you.
 
 Kevin
 
  Mark Difford [EMAIL PROTECTED] wrote: 
 
 Hi Kevin,
 
  The documentation indicates that the bw is essentially the sd.
   d - density(rnorm(1000))
 
 Not so. The documentation states that the following about bw: The
 kernels
 are scaled such that this is the standard deviation of the smoothing
 kernel..., which is a very different thing.
 
 The default bandwidth used by density is ?bw.nrd0. Read that
 documentation
 carefully and all might be clear.
 
 HTH, Mark.
 
 
 rkevinburton wrote:
  
  I issue the following:
  
  d - density(rnorm(1000))
  d
  
  and get:
  
  Call:
  density.default(x = rnorm(1000))
  
  Data: rnorm(1000) (1000 obs.);  Bandwidth 'bw' = 0.2235
  
 x y
   Min.   :-3.5157   Min.   :2.416e-05  
   1st Qu.:-1.6892   1st Qu.:1.129e-02  
   Median : 0.1373   Median :7.267e-02  
   Mean   : 0.1373   Mean   :1.367e-01  
   3rd Qu.: 1.9639   3rd Qu.:2.693e-01  
   Max.   : 3.7904   Max.   :4.014e-01  
  
  The documentation indicates that the bw is essentially the sd. Yet I
 have
  specified an sd of 1? How am I to interpret the ranges of the values?
 x
  ranges almost from -4 to +4 and y ranges from 0 to 0.4. The mean x is
 .1
  which isn't too awfully close to what I would expect (0.0). Then there
 is:
  
  d - density(rpois(1000,0))
  d
  
  Call:
  density.default(x = rpois(1000, 0))
  
  Data: rpois(1000, 0) (1000 obs.);   Bandwidth 'bw' = 0.2261
  
 x y  
   Min.   :-0.6782   Min.   :0.01979  
   1st Qu.:-0.3391   1st Qu.:0.14073  
   Median : 0.   Median :0.57178  
   Mean   : 0.   Mean   :0.73454  
   3rd Qu.: 0.3391   3rd Qu.:1.32830  
   Max.   : 0.6782   Max.   :1.76436  
  
  Here I am getting the mean that I expect from a Poisson distribuition
 but
  y ranges from 0 to 1.75. Again I am not sure what these numbers mean.
 How
  can I map the output to the standard distirbution description
 parameters?
  
  Thank you.
  
  Kevin
  
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
  http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
  
  
 
 -- 
 View this message in context:
 http://www.nabble.com/Help-interpreting-density%28%29.-tp18704955p18706154.html
 Sent from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, 

Re: [R] how to add notes to the graph?

2008-07-29 Thread Jim Lemon
On Mon, 2008-07-28 at 06:32 -0700, rlearner309 wrote:
 Hi,  
 I have a simple graph:
 x - c(1,2,3)
 plot(x, pch=16,type=b)
 
 I would like to add some notes just beside these 3 dots, and the notes are
 stored in a vector:
 
 a - c(12,54,84)  
 
 So the result will be: there should be a 12 below the first dot (or next
 to it, but not replacing the solid dot), a 54 next to the second dot...
 
 what if a is not numeric? a - c(a,b,c)
Hi rlearner309,
Have a look at thigmophobe.labels in the plotrix package.

Jim

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


[R] Problem reading a particular file with read.spss()

2008-07-29 Thread Chuck Cleland

Hi All:
  I have a seemingly typical SPSS data file with 219 rows and 486 
variables.  When I attempt to read this file into R with read.spss() in 
the foreign package, I consistently get a crash.  This is the sequence 
of events:


 library(foreign)
 sessionInfo()
R version 2.7.1 Patched (2008-07-24 r46120)
i386-pc-mingw32

locale:
LC_COLLATE=English_United States.1252;LC_CTYPE=English_United 
States.1252;LC_MONETARY=English_United 
States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252


attached base packages:
[1] stats graphics  grDevices utils datasets  methods   base

other attached packages:
[1] foreign_0.8-27

 read.spss(c:/chuck/mySPSSfile.sav)

R for Windows GUI front-end has encountered a problem and needs to 
close.  We are sorry for the inconvenience.


  It get a similar result trying to load the file under Rterm.exe.  The 
file opens without any trouble in SPSS version 16.0.2 (April 10, 2008). 
 For anyone willing to investigate, I placed a copy of the SPSS file here:


http://www.chuckcleland.net/mySPSSfile.sav

  I would be grateful if someone could help to figure out why this file 
causes a problem.


thanks,

Chuck

--
Chuck Cleland, Ph.D.
NDRI, Inc. (www.ndri.org)
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@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Is there way to multiple plots on gap.plot?

2008-07-29 Thread Jim Lemon
On Mon, 2008-07-28 at 12:15 -0700, Arthur Roberts wrote:
 Hi, all,
 
 Does anyone now of a way to put multiple plots on gap.plot?
 
 Much appreciated,

Hi Art,
You must have read my mind. In solving the problem you had with
gap.plot, I considered including an add argument that would allow the
user to add values to an existing plot. The bad news is that I decided
against it. As the examples show, you can plot several series of data
simultaneously with not very much trouble. However, if you can convince
me that add is a good idea, I might do it. Right now, though, I would
appreciate your testing of the new gap.plot function below.

Jim

gap.plot-function(x,y,gap,gap.axis=y,bgcol=white,breakcol=black,
 brw=0.02,xlim,ylim,xticlab,xtics=NA,yticlab,ytics=NA,lty=rep(1,length(x)),
 col=rep(par(col),length(x)),pch=rep(1,length(x)),...) {

 if(missing(y)  !missing(x)) {
  y-x
  x-1:length(y)
 }
 if(missing(gap)) stop(gap must be specified)
 gapsize-diff(gap)
 figxy - par(usr)
 xaxl-par(xlog)
 yaxl-par(ylog)
 xgw-(figxy[2]-figxy[1])*brw
 ygw-(figxy[4]-figxy[3])*brw
 if(missing(xtics)) xtics-pretty(x)
 if(missing(ytics)) ytics-pretty(y)
 if(missing(xticlab)) xticlab-xtics
 if(missing(yticlab)) yticlab-ytics
 if(length(col)  length(y)) col-rep(col,length.out=length(y))
 if(gap.axis == y) {
  littleones-which(y  gap[1])
  if(length(gapsize)  2) {
   middleones-which(y = gap[2]  y  gap[3])
   bigones-which(y = gap[4])
   lostones-sum(c(y  gap[1]  y  gap[2], y  gap[3]  y  gap[4]))
   if(missing(ylim)) ylim-c(min(y),ygw*2
+max(y)-(gapsize[1]+gapsize[3]))
   else ylim[2]-ygw*2+ylim[2]-(gapsize[1]+gapsize[3])
  }
  else {
   middleones-NA
   bigones-which(y = gap[2])
   lostones-sum(y  gap[1]  y  gap[2])
   if(missing(ylim)) ylim-c(min(y),max(y)-gapsize[1])
   else ylim[2]-ygw+ylim[2]-gapsize[1]
  }
  if(lostones) warning(some values of y will not be displayed)
  if(missing(xlim)) xlim-range(x)
 }
 else {
  littleones-which(x  gap[1])
  if(length(gapsize)  2) {
   middleones-which(x = gap[2]  x  gap[3])
   bigones-which(x = gap[4])
   lostones-sum(c(x  gap[1]  x  gap[2], x  gap[3]  x  gap[4]))
   if(missing(xlim)) xlim-c(min(x),xgw*2
+max(x)-(gapsize[1]+gapsize[3]))
   else xlim[2]-xgw*2+xlim[2]-(gapsize[1]+gapsize[3])
  }
  else {
   middleones-NA
   bigones-which(x = gap[2])
   lostones-sum(x  gap[1]  x  gap[2])
   if(missing(xlim)) xlim-c(min(x),max(x)-gapsize[1])
   else xlim[2]-xgw+xlim[2]-gapsize[1]
  }
  if(lostones) warning(some values of x will not be displayed)
  if(missing(ylim)) ylim-range(y)
 }
 if(length(lty)  length(x)) lty-rep(lty,length.out=length(x))
 if(length(col)  length(x)) col-rep(col,length.out=length(x))
 if(length(pch)  length(x)) pch-rep(pch,length.out=length(x))
 plot(x[littleones],y[littleones],xlim=xlim,ylim=ylim,axes=FALSE,
  lty=lty[littleones],col=col[littleones],pch=pch[littleones],...)
 box()
 if(gap.axis == y) {
  if(!is.na(xtics[1])) axis(1,at=xtics,labels=xticlab)
  littletics-which(ytics  gap[1])
  if(length(gapsize)  2) {
   middletics-which(ytics = gap[2]  ytics = gap[3])
   bigtics-which(ytics = gap[4])
   show.at-c(ytics[littletics],ytics[middletics]-gapsize[1],],...)
ytics[bigtics]-(gapsize[1]+gapsize[3]))
   show.labels-c(yticlab[littletics],yticlab[middletics],yticlab
   [bigtics])),
  }
  else {
   bigtics-which(ytics = gap[2])
   show.at-c(ytics[littletics],ytics[bigtics]-gapsize[1])
   show.labels-c(ytics[littletics],yticlab[bigtics])
  }
  axis(2,at=show.at,labels=show.labels)
  axis.break(2,gap
  [1]-ygw,style=gap,bgcol=bgcol,breakcol=breakcol,brw=brw)
  if(length(gapsize)  2) {
   axis.break(2,gap[3]-(gapsize[1]+ygw),style=gap,bgcol=bgcol,
breakcol=breakcol,brw=brw)
   points(x[middleones],y[middleones]-gapsize[1],
lty=lty[middleones],col=col[middleones],pch=pch[middleones],...)
   points(x[bigones],y[bigones]-(gapsize[1]+gapsize[3]),
lty=lty[bigones],col=col[bigones],pch=pch[bigones],...)
  }
  else points(x[bigones],y[bigones]-gapsize[1],
   lty=lty[bigones],col=col[bigones],pch=pch[bigones],...)
 }
 else {
  if(!is.na(ytics[1])) axis(2,at=ytics,labels=yticlab)
  littletics-which(xticsgap[1])
  if(length(gapsize)  2) {
   middletics-which(xtics = gap[2]  xtics = gap[3])
   bigtics-which(xtics = gap[4])
   show.at-c(xtics[littletics],xtics[middletics]-gapsize[1],],...)
xtics[bigtics]-(gapsize[1]+gapsize[3]))

show.labels-c(xticlab[littletics],xticlab[middletics],xticlab[bigtics])
  }
  else {
   bigtics-which(xtics = gap[2])
   show.at-c(xtics[littletics],xtics[bigtics]-gapsize[1])
   show.labels-c(xticlab[littletics],xticlab[bigtics])
  }
  axis(1,at=show.at,labels=show.labels)
  axis.break(1,gap[1]-xgw,style=gap)
  if(length(gapsize)  2) {
   axis.break(1,gap[3]-(gapsize[1]+xgw),style=gap)
   points(xgw+x[middleones]-gapsize[1],y[middleones],])
lty=lty[middleones],col=col[middleones],pch=pch[middleones],...)
   points(x[bigones]-(gapsize[1]+gapsize[3]),y[bigones],
lty=lty[bigones],col=col[bigones],pch=pch[bigones],...)
  }
  else 

[R] keywords

2008-07-29 Thread Edna Bell
Hi R Gurus!

When you build a package, you need to put in keywords in the Rd files.

Where would you find the list of keywords, please?

TIA,
Edna Bell

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


Re: [R] list

2008-07-29 Thread jim holtman
Is this what you are after: this comes from the example of wilcox.test
and shows the structure and how to reference the p-value:

 z - wilcox.test(x, y, paired = TRUE, alternative = greater)
 str(z)
List of 7
 $ statistic  : Named num 40
  ..- attr(*, names)= chr V
 $ parameter  : NULL
 $ p.value: num 0.0195
 $ null.value : Named num 0
  ..- attr(*, names)= chr location shift
 $ alternative: chr greater
 $ method : chr Wilcoxon signed rank test
 $ data.name  : chr x and y
 - attr(*, class)= chr htest
 z$p.value
[1] 0.01953125


On Tue, Jul 29, 2008 at 1:00 AM, Paul Adams [EMAIL PROTECTED] wrote:
 I have run a wilcox test on a dataframe and I have returned to me data which
 I am wanting to sort and/or pick the max.
 I have tried the following code to pick the max:
 v-wilcox.test.run
 w-max(v,na.rm=FALSE)

 I have also tried w-pmax(p-value,na.rm=FALSE) for the second line

 and this returns the following error:invalid 'type'  (list) of argument
 This error is understandable because I am trying to find the max p-value from 
 a list of
 100 with the following output:


 $g100
 Wilcoxon signed rank test
 data:newX[,i]
 V=741,  p-value=8.051e-08

 I want to be able to pull only the maximum p-value from this output and the 
 $g  that matches the p-value by using max.Any help would be appreciated.

 Thanks
 Paul





[[alternative HTML version deleted]]


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





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

What is the problem you are trying to solve?

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


[R] Bug in sd() and var() in handling vectors of NA (R version 2.7.1)?

2008-07-29 Thread Marisa Laetitia
In the previous versions of R (2.6.1), when a vector of NA was given to the 
functions 'sd' or 'var' with parameter na.rm = TRUE, it used to return NA. Now 
(2.7.1) it returns an ERROR : 

 

Example in 2.6.1: 

 sd(c(NA, NA, NA, NA), na.rm = TRUE)

[1] NA

 

Example in 2.7.1:

 sd(c(NA, NA, NA, NA), na.rm = TRUE)

Error in var(x, na.rm = na.rm) : paires d'éléments incomplètes

 

We are actually wondering if it is a bug to report, or if we have to manage it 
in our own R scripts.

 

Thanks a lot, best regards,

 

Laetitia

 








This footnote confirms that this email message has been scanned by
PineApp Mail-SeCure for the presence of malicious code, vandals  computer 
viruses.




[[alternative HTML version deleted]]

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


Re: [R] Problem reading a particular file with read.spss()

2008-07-29 Thread Prof Brian Ripley
It works for me, on Windows with that version of foreign.  I also tried 
under Linux with valgrind, and nothing untoward was reported.  The 
warnings are


Warning messages:
1: In read.spss(mySPSSfile.sav) :
  mySPSSfile.sav: File-indicated character representation code (1252) 
looks like a Windows codepage

2: In read.spss(mySPSSfile.sav) :
  mySPSSfile.sav: Unrecognized record type 7, subtype 18 encountered in 
system file


The first is innocuous for you (given your locale), and the second 
indicated that there is info that we don't know about (SPSS 16.x is not 
common).


On Tue, 29 Jul 2008, Chuck Cleland wrote:


Hi All:
 I have a seemingly typical SPSS data file with 219 rows and 486 variables. 
When I attempt to read this file into R with read.spss() in the foreign 
package, I consistently get a crash.  This is the sequence of events:



library(foreign)
sessionInfo()

R version 2.7.1 Patched (2008-07-24 r46120)
i386-pc-mingw32

locale:
LC_COLLATE=English_United States.1252;LC_CTYPE=English_United 
States.1252;LC_MONETARY=English_United 
States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252


attached base packages:
[1] stats graphics  grDevices utils datasets  methods   base

other attached packages:
[1] foreign_0.8-27


read.spss(c:/chuck/mySPSSfile.sav)


R for Windows GUI front-end has encountered a problem and needs to close. 
We are sorry for the inconvenience.


 It get a similar result trying to load the file under Rterm.exe.  The file 
opens without any trouble in SPSS version 16.0.2 (April 10, 2008).  For 
anyone willing to investigate, I placed a copy of the SPSS file here:


http://www.chuckcleland.net/mySPSSfile.sav

 I would be grateful if someone could help to figure out why this file 
causes a problem.


thanks,

Chuck

--
Chuck Cleland, Ph.D.
NDRI, Inc. (www.ndri.org)
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@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.



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

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


[R] max

2008-07-29 Thread Paul Adams
Hello everyone,
I have a vector of p-values and I am trying to find the max for the vector and 
add it to a list.I am using a loop to loop through 50 times.I have used the 
following code but it does not pick out the max and add to the list.Any help 
would be appreciated.
 
w-c(v[[1]][3],v[[2]][3],v[[3]][3],v[[4]][3],v[[5]][3],v[[6]][3],v[[7]][3],v[[8]][3],v[[9]][3],
v[[10]][3],v[[11]][3],v[[12]][3],v[[13]][3],v[[14]][3],v[[15]][3],v[[16]][3],v[[17]][3],v[[18]][3],
v[[19]][3],v[[20]][3],v[[21]][3],v[[22]][3],v[[23]][3],v[[24]][3],v[[25]][3],v[[26]][3],v[[27]][3],
v[[28]][3],v[[29]][3],v[[30]][3],v[[31]][3],v[[32]][3],v[[33]][3],v[[34]][3],v[[35]][3],v[[36]][3],
v[[37]][3],v[[38]][3],v[[39]][3],v[[40]][3],v[[41]][3],v[[42]][3],v[[43]][3],v[[44]][3],v[[45]][3],
v[[46]][3],v[[47]][3],v[[48]][3],v[[49]][3],v[[50]][3],v[[51]][3],v[[52]][3],v[[54]][3],v[[55]][3],
v[[56]][3],v[[57]][3],v[[58]][3],v[[59]][3],v[[60]][3],v[[61]][3],v[[62]][3],v[[63]][3],v[[64]][3],
v[[65]][3],v[[66]][3],v[[67]][3],v[[68]][3],v[[69]][3],v[[70]][3],v[[71]][3],v[[72]][3],v[[73]][3],
v[[74]][3],v[[75]][3],v[[76]][3],v[[77]][3],v[[78]][3],v[[79]][3],v[[80]][3],v[[81]][3],v[[82]][3],
v[[83]][3],v[[84]][3],v[[85]][3],v[[86]][3],v[[87]][3],v[[88]][3],v[[89]][3],v[[90]][3],v[[91]][3],
v[[92]][3],v[[93]][3],v[[94]][3],v[[95]][3],v[[96]][3],v[[97]][3],v[[98]][3],v[[99]][3],v[[100]][3])

w1-max(w,na.rm=FALSE)
list-c(list,w1)

This runs but it does not pick out the max of the vector nor does it add it to 
the list.The for loop shuffles the vector then iterates 50 times and there 
should be 50 values in the list.I created the list outside the for loop with 
list-NULL.When I look at the list I get the whole vector 50 times(5000 values) 
so I assume that max is not picking out the maximum value of the vector 
everytime through the loop.
Any help would be appreciated
Paul


  
[[alternative HTML version deleted]]

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


[R] FW: Installing BRugs

2008-07-29 Thread BXC (Bendix Carstensen)
A funny thing happened when I wanted a student of mine to install Brugs.
Using the InstallPackages in the windows version, firts gives an erro, but 
trying again works flawlessly.

R version is 2.7.0 on WinXP.

Any explanation?

Bendix Carstensen
__

Bendix Carstensen
Senior Statistician
Steno Diabetes Center
Niels Steensens Vej 2-4
DK-2820 Gentofte
Denmark
+45 44 43 87 38 (direct)
+45 30 75 87 38 (mobile)
[EMAIL PROTECTED]   http://www.biostat.ku.dk/~bxc

-Original Message-
From: ACJS (Anders Christian Jensen)
Sent: 29. juli 2008 11:06
To: BXC (Bendix Carstensen)
Subject:

 utils:::menuInstallPkgs()
trying URL 
'http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.7/BRugs_0.4-2.zip'
Content type 'application/zip' length 3399130 bytes (3.2 Mb) opened URL 
downloaded 3.1 Mb

Error in gzfile(file, r) : cannot open the connection In addition: Warning 
messages:
1: In download.file(url, destfile, method, mode = wb, ...) :
  downloaded length 3293996 != reported length 3399130
2: In zip.unpack(pkg, tmpDir) : error 1 in extracting from zip file
3: In gzfile(file, r) :
  cannot open compressed file 'BRugs/DESCRIPTION', probable reason 'No such 
file or directory'
 utils:::menuInstallPkgs()
trying URL 
'http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.7/BRugs_0.4-2.zip'
Content type 'application/zip' length 3399130 bytes (3.2 Mb) opened URL 
downloaded 3.2 Mb

package 'BRugs' successfully unpacked and MD5 sums checked

The downloaded packages are in
c:\temp\RtmpKoBifa\downloaded_packages
updating HTML package descriptions


___

Anders Christian Jensen

Steno Diabetes Center
Niels Steensens Vej 2-4
DK-2820 Gentofte
Denmark
+45 444 22475 (direct)
[EMAIL PROTECTED]

This e-mail (including any attachments) is intended for ...{{dropped:8}}

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


Re: [R] FW: Installing BRugs

2008-07-29 Thread Uwe Ligges



BXC (Bendix Carstensen) wrote:

A funny thing happened when I wanted a student of mine to install Brugs.
Using the InstallPackages in the windows version, firts gives an erro, but 
trying again works flawlessly.

R version is 2.7.0 on WinXP.

Any explanation?



No, particularly if it worked the second time. Perhaps a broken machine 
or a broken internte connection..

Uwe



Bendix Carstensen
__

Bendix Carstensen
Senior Statistician
Steno Diabetes Center
Niels Steensens Vej 2-4
DK-2820 Gentofte
Denmark
+45 44 43 87 38 (direct)
+45 30 75 87 38 (mobile)
[EMAIL PROTECTED]   http://www.biostat.ku.dk/~bxc

-Original Message-
From: ACJS (Anders Christian Jensen)
Sent: 29. juli 2008 11:06
To: BXC (Bendix Carstensen)
Subject:


utils:::menuInstallPkgs()

trying URL 
'http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.7/BRugs_0.4-2.zip'
Content type 'application/zip' length 3399130 bytes (3.2 Mb) opened URL 
downloaded 3.1 Mb

Error in gzfile(file, r) : cannot open the connection In addition: Warning 
messages:
1: In download.file(url, destfile, method, mode = wb, ...) :
  downloaded length 3293996 != reported length 3399130
2: In zip.unpack(pkg, tmpDir) : error 1 in extracting from zip file
3: In gzfile(file, r) :
  cannot open compressed file 'BRugs/DESCRIPTION', probable reason 'No such 
file or directory'

utils:::menuInstallPkgs()

trying URL 
'http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.7/BRugs_0.4-2.zip'
Content type 'application/zip' length 3399130 bytes (3.2 Mb) opened URL 
downloaded 3.2 Mb

package 'BRugs' successfully unpacked and MD5 sums checked

The downloaded packages are in
c:\temp\RtmpKoBifa\downloaded_packages
updating HTML package descriptions

___

Anders Christian Jensen

Steno Diabetes Center
Niels Steensens Vej 2-4
DK-2820 Gentofte
Denmark
+45 444 22475 (direct)
[EMAIL PROTECTED]

This e-mail (including any attachments) is intended for ...{{dropped:8}}

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


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


[R] mac install no font found problem in quartz display

2008-07-29 Thread Zhesi He

Dear list,

I was attempted to install 2.7.1 in my mac only recently because of  
the fantastic quartz display support.

I had no problem in building from source.
But after installation, I always got font warnings as following when  
using quartz:


Warning messages:
1: In axis(side = side, at = at, labels = labels, ...) ... :
  no font could be found for family Arial

33: In axis(side = side, at = at, labels = labels, ...) ... :
  no font could be found for family Arial
34: In title(...) ... : no font could be found for family Arial
35: In title(...) ... : no font could be found for family Arial
36: In title(...) ... : no font could be found for family Arial
37: In title(...) ... : no font could be found for family Arial


There is no such problem when using X11 as the display device.
I was guessing it might be cairo or fontconfig library? I looked into  
the configure file and found warnings:

ld: warning, duplicate dylib /opt/local/lib/libiconv.2.dylib
ld: warning, duplicate dylib /opt/local/lib/libz.1.dylib
ld: warning, duplicate dylib /opt/local/lib/libfontconfig.1.dylib
ld: warning, duplicate dylib /opt/local/lib/libXrender.1.dylib
ld: warning, duplicate dylib /usr/X11/lib/libSM.6.dylib
ld: warning, duplicate dylib /usr/X11/lib/libICE.6.dylib
ld: warning, duplicate dylib /usr/X11/lib/libX11.6.dylib
ld: warning, duplicate dylib /usr/local/lib/libfreetype.6.dylib

I used macport to install cairo and pango and gtk before, and I think  
these libraries are duplicated in /usr/lib somehow.
I don't know if R looked into the wrong directory to search for  
libraries or not.


Please advise!

Thanks a lot.


--
Zhesi He
Graham Group, CNAP
Department of Biology (Area 7)
University of York
PO Box 373, York YO10 5YW, UK
Phone:  +44-(0)1904-328774
Email:  [EMAIL PROTECTED]

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


Re: [R] product of successive rows

2008-07-29 Thread Jorge Ivan Velez
Hi rcoder,

Assuming that the number of rows of your matrix x is even, try also:

x - matrix(1:72,12)
apply(x,2, tapply, rep(1:(nrow(x)/2),each=2),prod)

# or using a function which argument x is your matrix

prod.mat=function(x) {
k=nrow(x)
g=rep(1:(k/2),each=2)
apply(x,2, tapply, g,prod)
}

prod.mat(x)


HTH,

Jorge




On Sun, Jul 27, 2008 at 6:20 PM, rcoder [EMAIL PROTECTED] wrote:


 Hi everyone,

 I want to perform an operation on a matrx that outputs the product of
 successive pairs of rows. For example: calculating the product between rows
 1  2; 3  4; 5  6...etc.

 Does anyone know of any readily available functions that can do this?

 Thanks,

 rcoder


 --
 View this message in context:
 http://www.nabble.com/product-of-successive-rows-tp18681259p18681259.html
 Sent from the R help mailing list archive at Nabble.com.

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


[[alternative HTML version deleted]]

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


[R] stringdot ?

2008-07-29 Thread renata.camargo05
Dear all,

I am using kernlab package in R, and I have amino acid sequences with different 
lenghts as input for a SVM and I need to go through this sequences using 
windows (sliding or fixed) of size X. 
Does anyone has any suggestions about which function I should use?
I thought I could use stringdot, but I am not sure whether it will do what I 
need.., I have defined my stringdot as:

mystringdot - stringdot(length = 7, lambda = 0.5, type = sequence, 
normalized = TRUE)

and my svm as:

ksvm(mydata[,1],data=mydata,type=C-svc,kernel=mystringdot,C=10)

but it doesn't work.

Thank you ever so much, Renata

-
Email sent from www.virginmedia.com/email
Virus-checked using McAfee(R) Software and scanned for spam

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


[R] optim fails when using arima

2008-07-29 Thread MAIDER MATEOS DEL PINO


Hi all,

I´m using the arima() function to study a time series but it gives me  
the following error:


Error en optim(init[mask], armafn, method = BFGS, hessian = TRUE,  
control = optim.control,  :

  non-finite finite-difference value [3]

I know that I can change the method of the arima() to CSS instead of  
ML but I'm specially interested in using maximum likelihood.


I have read that using the Nelder-Mead method in the optim() function  
could avoid this error but I think that it is not possible to change  
the method of optim from arima().


Does anyone have an idea of how could I solve this problem?

Thanks and regards,

Maider Mateos

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


[R] 'for' loop, two variables

2008-07-29 Thread Oehler, Friderike (AGPP)
Dear Rusers,
I am still an unexperienced builder of functions and loops, so my question is
very basic: Is it possible to introduce a second variable (j) into my loop.
To examplify:

# This works fine:
fn - function (x) {if (x46  x52) 1 else 0}
res -NULL 
for (i in 40:60) res -c(res,fn(i))
res

# But here, there is an error in the for expression:
fn - function (x,y) {if (x46  x52  y12) 1 else 0 }
res -NULL 
for (i in 40:60  j in 0:20) res -c(res,fn(i,j)) 
# How do I have to write the expression i in 40:60  j in 0:20? Or is there
no way to do that, i.e. I have to do the calculation in two steps?

Thanks in advance!
Friderike

[[alternative HTML version deleted]]

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


Re: [R] Is there way to multiple plots on gap.plot?

2008-07-29 Thread Jim Lemon
Hi Art,
Ignore the last email, I realized that I had already written the add
bit and it will be in plotrix 2.4-5. As far as I can see, both the
type and xlim problems are solved.

Jim

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


Re: [R] 'for' loop, two variables

2008-07-29 Thread Duncan Murdoch

On 7/29/2008 7:55 AM, Oehler, Friderike (AGPP) wrote:

Dear Rusers,
I am still an unexperienced builder of functions and loops, so my question is
very basic: Is it possible to introduce a second variable (j) into my loop.
To examplify:

# This works fine:
fn - function (x) {if (x46  x52) 1 else 0}
res -NULL 
for (i in 40:60) res -c(res,fn(i))

res

# But here, there is an error in the for expression:
fn - function (x,y) {if (x46  x52  y12) 1 else 0 }
res -NULL 
for (i in 40:60  j in 0:20) res -c(res,fn(i,j)) 
# How do I have to write the expression i in 40:60  j in 0:20? Or is there

no way to do that, i.e. I have to do the calculation in two steps?


You need two steps.  You probably want either

for (i in 40:60) for (j in 0:20) res -c(res,fn(i,j))

or

i - 40:60
j - 0:20
for (ind in seq_along(i)) res - c(res, fn(i[ind], j[ind]))

(which do quite different loops).

Duncan Murdoch

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


Re: [R] 'for' loop, two variables

2008-07-29 Thread ONKELINX, Thierry
Dear Frederike,

#Both your functions are vectorized. So you don't need loops. Working
with vectorized functions is much faster than looping.

fn - function (x) {
ifelse(x46  x52, 1, 0)
}
res - fn(40:60)

fn - function (x,y) {
ifelse(x46  x52  y12, 1, 0)
}
datagrid - expand.grid(i = 40:60, j = 0:20)
res - fn(datagrid$i, datagrid$j)

#An other option is to use the functions for the apply-family

fn - function (x) {
ifelse(x46  x52, 1, 0)
}
res - sapply(40:60, fn)

fn - function (x,y) {
ifelse(x46  x52  y12, 1, 0)
}
datagrid - expand.grid(i = 40:60, j = 0:20)
res - apply(datagrid, 1, function(z){
fn(z[i], z[j])
}) 

#or you can use a nested loop

res -NULL 
for (i in 40:60){
for(j in 0:20){
res -c(res,fn(i,j))
}
}

HTH,

Thierry



ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature
and Forest
Cel biometrie, methodologie en kwaliteitszorg / Section biometrics,
methodology and quality assurance
Gaverstraat 4
9500 Geraardsbergen
Belgium 
tel. + 32 54/436 185
[EMAIL PROTECTED] 
www.inbo.be 

To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to
say what the experiment died of.
~ Sir Ronald Aylmer Fisher

The plural of anecdote is not data.
~ Roger Brinner

The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of
data.
~ John Tukey

-Oorspronkelijk bericht-
Van: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED]
Namens Oehler, Friderike (AGPP)
Verzonden: dinsdag 29 juli 2008 13:56
Aan: Oehler, Friderike (AGPP); r-help@r-project.org
Onderwerp: [R] 'for' loop, two variables

Dear Rusers,
I am still an unexperienced builder of functions and loops, so my
question is
very basic: Is it possible to introduce a second variable (j) into my
loop.
To examplify:

# This works fine:
fn - function (x) {if (x46  x52) 1 else 0}
res -NULL 
for (i in 40:60) res -c(res,fn(i))
res

# But here, there is an error in the for expression:
fn - function (x,y) {if (x46  x52  y12) 1 else 0 }
res -NULL 
for (i in 40:60  j in 0:20) res -c(res,fn(i,j)) 
# How do I have to write the expression i in 40:60  j in 0:20? Or is
there
no way to do that, i.e. I have to do the calculation in two steps?

Thanks in advance!
Friderike

[[alternative HTML version deleted]]

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

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


Re: [R] 'for' loop, two variables

2008-07-29 Thread Benno Pütz


On 29.Jul.2008, at 14:13, ONKELINX, Thierry wrote:


Dear Frederike,

#Both your functions are vectorized. So you don't need loops. Working
with vectorized functions is much faster than looping.

fn - function (x,y) {
   ifelse(x46  x52  y12, 1, 0)
}
datagrid - expand.grid(i = 40:60, j = 0:20)
res - apply(datagrid, 1, function(z){
   fn(z[i], z[j])
})



or

outer(40:60,0:20,fn2)

which also keeps the matrix structure ...

(use as.vector(...) or as.vector(t(...)) if you need a vector)

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


Re: [R] FW: Installing BRugs

2008-07-29 Thread Prof Brian Ripley

On Tue, 29 Jul 2008, BXC (Bendix Carstensen) wrote:


A funny thing happened when I wanted a student of mine to install Brugs.
Using the InstallPackages in the windows version, firts gives an erro, but 
trying again works flawlessly.


There was a problem with the http connection.  It happens 



R version is 2.7.0 on WinXP.

Any explanation?

Bendix Carstensen
__

Bendix Carstensen
Senior Statistician
Steno Diabetes Center
Niels Steensens Vej 2-4
DK-2820 Gentofte
Denmark
+45 44 43 87 38 (direct)
+45 30 75 87 38 (mobile)
[EMAIL PROTECTED]   http://www.biostat.ku.dk/~bxc

-Original Message-
From: ACJS (Anders Christian Jensen)
Sent: 29. juli 2008 11:06
To: BXC (Bendix Carstensen)
Subject:


utils:::menuInstallPkgs()

trying URL 
'http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.7/BRugs_0.4-2.zip'
Content type 'application/zip' length 3399130 bytes (3.2 Mb) opened URL 
downloaded 3.1 Mb

Error in gzfile(file, r) : cannot open the connection In addition: Warning 
messages:
1: In download.file(url, destfile, method, mode = wb, ...) :
 downloaded length 3293996 != reported length 3399130
2: In zip.unpack(pkg, tmpDir) : error 1 in extracting from zip file
3: In gzfile(file, r) :
 cannot open compressed file 'BRugs/DESCRIPTION', probable reason 'No such file 
or directory'

utils:::menuInstallPkgs()

trying URL 
'http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.7/BRugs_0.4-2.zip'
Content type 'application/zip' length 3399130 bytes (3.2 Mb) opened URL 
downloaded 3.2 Mb

package 'BRugs' successfully unpacked and MD5 sums checked

The downloaded packages are in
   c:\temp\RtmpKoBifa\downloaded_packages
updating HTML package descriptions




___

Anders Christian Jensen

Steno Diabetes Center
Niels Steensens Vej 2-4
DK-2820 Gentofte
Denmark
+45 444 22475 (direct)
[EMAIL PROTECTED]

This e-mail (including any attachments) is intended for ...{{dropped:8}}

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



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

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


Re: [R] Bug in sd() and var() in handling vectors of NA (R version 2.7.1)?

2008-07-29 Thread Prof Brian Ripley
There was a bug in 2.6.1 which has since been corrected: there is no need 
to report corrected bugs in obsolete versions.  Two different ways to 
compute the sd of a zero-length vector gave different answers.


This is covered by the following NEWS item for 2.7.0 (version from 
R-patched)


o   co[rv](use = complete.obs) now always gives an error if there
are no complete cases: they used to give NA if
method = pearson but an error for the other two methods.
(Note that this is pretty arbitrary, but zero-length vectors
always give an error so it is at least consistent.)

Since sd(na.rm=TRUE) and var(na.rm=TRUE) both call cov(use =
complete.obs), this applies also to them.


On Tue, 29 Jul 2008, Marisa Laetitia wrote:


In the previous versions of R (2.6.1), when a vector of NA was given to the 
functions 'sd' or 'var' with parameter na.rm = TRUE, it used to return NA. Now 
(2.7.1) it returns an ERROR :



Example in 2.6.1:


sd(c(NA, NA, NA, NA), na.rm = TRUE)


[1] NA



Example in 2.7.1:


sd(c(NA, NA, NA, NA), na.rm = TRUE)


Error in var(x, na.rm = na.rm) : paires d'?l?ments incompl?tes



We are actually wondering if it is a bug to report, or if we have to manage it 
in our own R scripts.



Thanks a lot, best regards,



Laetitia










This footnote confirms that this email message has been scanned by
PineApp Mail-SeCure for the presence of malicious code, vandals  computer 
viruses.




[[alternative HTML version deleted]]




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

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


[R] Panel of pie charts

2008-07-29 Thread Amin Momin
Hi ,
  I am looking to making a panel of pie charts fo some of my dritribution data 
. I was wondering if there is a way in any R package to write a small script to 
do so.

  Thanks,
   Amin 
[[alternative HTML version deleted]]

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


Re: [R] Negative Binomial Regression

2008-07-29 Thread Ben Bolker
jcarmichael jcarmichael314 at gmail.com writes:

 
 
 Hello.
 
 I am attempting to duplicate a negative binomial regression in R.  SAS uses
 generalized estimating equations for model fitting in the GENMOD procedure.
 
 proc genmod data=mydata (where=(gender='F')); 
 by agegroup;
 class id gender type;
 model count = var1 var2 var3 /dist=NB link=log offset=lregtm;
 repeated subject=id /type=exch;  
 run;
 
 Since my dataset has several observations for each subject, I need the
 REPEATED statement in order to indicate dependence among observations with
 the same subject ID and independence amongst those with distinct subject
 IDs.  The TYPE statement goes on to specify the structure of the correlation
 matrix to be used (exchangeable in this case).

  I would try glmmPQL in the MASS package.  I don't think you
can *quite* get negative binomial regression this way, but
you can definitely get a quasipoisson model.  I think exchangeable
correlation corresponds to correlation=corCompSymm() in your
glmmPQL command.

  Ben Bolker

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


Re: [R] How to set the parameters in Trellis Graphics (by Lattice package)

2008-07-29 Thread Felix Andrews
On Tue, Jul 29, 2008 at 5:37 PM, G.H. Zuo [EMAIL PROTECTED] wrote:
 Dear R users

 I plot the trellis graphics by using the lattice package.  Everything is
 OK.  Now, I want set
 some parameters of the trellis graphics.

 1. The tick label site.  By default, only two tick labels had been output in
 x-axis of my plot.
 I want output four or five tick labels.  In the traditional graphics system,
 it will be very easy.
 Just not output the axis in plot by use the parameter xaxt=n, and then
 used axis to add
 the axis.  But it seem useless in the grid graphics system.  Then how can I
 do it?

You would generally use scales=list(ticks=5) or scales=list(x=list(ticks=5)).
See the entry for scales in ?xyplot
However, the default (suggested) number of ticks is 4 or 5, so there
may be something odd about your case. How about a reproducible
example?


 2. How to change the margin of the figure?

Why do you want to do that? Grid graphics automatically allocates
enough space for the objects to be displayed. If you really want to,
there are ways to do it, but I suggest you do not. You are more likely
looking for the aspect argument.


 3. Can I change the style of the figure like par in traditional system.  For
 example, I want
 change all the fontface into bold.  Now, I must change the fontface for
 every item by using
 command trellis.par.set() like this:

 trellis.par.set(list(par.xlab.text = gpar(font = 2),
 par.ylab.text = gpar(font = 2),
 axis.text = gpar(font = 2),
 add.text = gpar(font=2)
   ))

 is there a method to change this parameter in one time, like par(font=2)

You can try trellis.par.set(fontsize=list(text=14)) which is similar.
Otherwise you could write a function to set the font parameters;
this could be modified from latticeExtra::custom.theme



 thanks

 G.H.Zuo

[[alternative HTML version deleted]]

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




-- 
Felix Andrews / 安福立
PhD candidate
Integrated Catchment Assessment and Management Centre
The Fenner School of Environment and Society
The Australian National University (Building 48A), ACT 0200
Beijing Bag, Locked Bag 40, Kingston ACT 2604
http://www.neurofractal.org/felix/
3358 543D AAC6 22C2 D336 80D9 360B 72DD 3E4C F5D8

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


[R] combining zoo series with an overlapping index?

2008-07-29 Thread stephen sefick
day-structure(c(7.7, 7.7, 7.7, 7.7, 7.7, 7.71, 7.7, 7.71, 7.71, 7.7,
7.7, 7.7, 7.7, 7.69, 7.68, 7.68, 7.67, 7.67, 7.67, 7.66, 7.65,
7.65, 7.65, 7.64, 7.64, 7.63, 7.63, 7.63, 7.62, 7.62, 7.62, 7.62,
7.63, 7.63, 7.63, 7.63, 7.63, 7.64, 7.64, 7.65, 7.65, 7.65, 7.66,
7.66, 7.67, 7.67, 7.67, 7.68, 7.68, 7.69, 7.69, 7.69, 7.69, 7.7,
7.7, 7.7, 7.7, 7.7, 7.71, 7.7, 7.7, 7.71, 7.71, 7.7, 7.7, 7.7,
7.7, 7.7, 7.7, 7.7, 7.69, 7.68, 7.68, 7.68, 7.67, 7.67, 7.66,
7.66, 7.66, 7.66, 7.65, 7.65, 7.65, 7.65, 7.66, 7.66, 7.66, 7.67,
7.67, 7.68, 7.68, 7.68, 7.68, 7.69, 7.69, 7.69), index =
structure(c(13307.000694,
13307.01, 13307.021528, 13307.031944, 13307.042361,
13307.052778, 13307.063194, 13307.073611, 13307.084028,
13307.09, 13307.104861, 13307.115278, 13307.125694,
13307.136111, 13307.146528, 13307.156944, 13307.167361,
13307.18, 13307.188194, 13307.198611, 13307.209028,
13307.219444, 13307.229861, 13307.240278, 13307.250694,
13307.26, 13307.271528, 13307.281944, 13307.292361,
13307.302778, 13307.313194, 13307.323611, 13307.334028,
13307.34, 13307.354861, 13307.365278, 13307.375694,
13307.386111, 13307.396528, 13307.406944, 13307.417361,
13307.427778, 13307.438194, 13307.448611, 13307.459028,
13307.469444, 13307.479861, 13307.490278, 13307.500694,
13307.51, 13307.521528, 13307.531944, 13307.542361,
13307.552778, 13307.563194, 13307.573611, 13307.584028,
13307.59, 13307.604861, 13307.615278, 13307.625694,
13307.636111, 13307.646528, 13307.656944, 13307.667361,
13307.68, 13307.688194, 13307.698611, 13307.709028,
13307.719444, 13307.729861, 13307.740278, 13307.750694,
13307.76, 13307.771528, 13307.781944, 13307.792361,
13307.802778, 13307.813194, 13307.823611, 13307.834028,
13307.84, 13307.854861, 13307.865278, 13307.875694,
13307.886111, 13307.896528, 13307.906944, 13307.917361,
13307.927778, 13307.938194, 13307.948611, 13307.959028,
13307.969444, 13307.979861, 13307.990278), format =
structure(c(m/d/y,
h:m:s), .Names = c(dates, times)), origin = structure(c(1,
1, 1970), .Names = c(month, day, year)), class = c(chron,
dates, times)), class = zoo)

#Mullholland ER daytime

# extract section of interest

w - window(day, start=chron(6/8/6, 4:16:0), end=chron(6/8/6,
20:31:0))

# zap all interior points with NA's and then fill back in using linear
interpolation

lin - na.approx(replace(w, time(w)  start(w)  time(w)  end(w), NA))

# plot the 3 series on one screen

#plot(merge(day, w, lin), col = 1:3, screen = 1)

#can I do this with a command (generically?)?  I am going to use this inside
of a function, and below is the result I would like.
d- c(day[1:17,], lin, day[84:96,])
plot(d)

-- 
Let's not spend our time and resources thinking about things that are so
little or so large that all they really do for us is puff us up and make us
feel like gods. We are mammals, and have not exhausted the annoying little
problems of being mammals.

-K. Mullis

[[alternative HTML version deleted]]

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


Re: [R] keywords

2008-07-29 Thread Charilaos Skiadas

On Jul 29, 2008, at 5:24 AM, Edna Bell wrote:


Hi R Gurus!

When you build a package, you need to put in keywords in the Rd files.

Where would you find the list of keywords, please?


Simplest way is to google for r keywords. First hit is:

http://www.stat.ucl.ac.be/ISdidactique/Rhelp/doc/keywords.html

Also, if you do help.start() and follow the link Search Engine   
Keywords, you will find the list.



TIA,
Edna Bell


Haris Skiadas
Department of Mathematics and Computer Science
Hanover College

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


Re: [R] combining zoo series with an overlapping index?

2008-07-29 Thread Gabor Grothendieck
Try:

window(d, time(lin)) - coredata(lin)

On Tue, Jul 29, 2008 at 9:01 AM, stephen sefick [EMAIL PROTECTED] wrote:
 day-structure(c(7.7, 7.7, 7.7, 7.7, 7.7, 7.71, 7.7, 7.71, 7.71, 7.7,
 7.7, 7.7, 7.7, 7.69, 7.68, 7.68, 7.67, 7.67, 7.67, 7.66, 7.65,
 7.65, 7.65, 7.64, 7.64, 7.63, 7.63, 7.63, 7.62, 7.62, 7.62, 7.62,
 7.63, 7.63, 7.63, 7.63, 7.63, 7.64, 7.64, 7.65, 7.65, 7.65, 7.66,
 7.66, 7.67, 7.67, 7.67, 7.68, 7.68, 7.69, 7.69, 7.69, 7.69, 7.7,
 7.7, 7.7, 7.7, 7.7, 7.71, 7.7, 7.7, 7.71, 7.71, 7.7, 7.7, 7.7,
 7.7, 7.7, 7.7, 7.7, 7.69, 7.68, 7.68, 7.68, 7.67, 7.67, 7.66,
 7.66, 7.66, 7.66, 7.65, 7.65, 7.65, 7.65, 7.66, 7.66, 7.66, 7.67,
 7.67, 7.68, 7.68, 7.68, 7.68, 7.69, 7.69, 7.69), index =
 structure(c(13307.000694,
 13307.01, 13307.021528, 13307.031944, 13307.042361,
 13307.052778, 13307.063194, 13307.073611, 13307.084028,
 13307.09, 13307.104861, 13307.115278, 13307.125694,
 13307.136111, 13307.146528, 13307.156944, 13307.167361,
 13307.18, 13307.188194, 13307.198611, 13307.209028,
 13307.219444, 13307.229861, 13307.240278, 13307.250694,
 13307.26, 13307.271528, 13307.281944, 13307.292361,
 13307.302778, 13307.313194, 13307.323611, 13307.334028,
 13307.34, 13307.354861, 13307.365278, 13307.375694,
 13307.386111, 13307.396528, 13307.406944, 13307.417361,
 13307.427778, 13307.438194, 13307.448611, 13307.459028,
 13307.469444, 13307.479861, 13307.490278, 13307.500694,
 13307.51, 13307.521528, 13307.531944, 13307.542361,
 13307.552778, 13307.563194, 13307.573611, 13307.584028,
 13307.59, 13307.604861, 13307.615278, 13307.625694,
 13307.636111, 13307.646528, 13307.656944, 13307.667361,
 13307.68, 13307.688194, 13307.698611, 13307.709028,
 13307.719444, 13307.729861, 13307.740278, 13307.750694,
 13307.76, 13307.771528, 13307.781944, 13307.792361,
 13307.802778, 13307.813194, 13307.823611, 13307.834028,
 13307.84, 13307.854861, 13307.865278, 13307.875694,
 13307.886111, 13307.896528, 13307.906944, 13307.917361,
 13307.927778, 13307.938194, 13307.948611, 13307.959028,
 13307.969444, 13307.979861, 13307.990278), format =
 structure(c(m/d/y,
 h:m:s), .Names = c(dates, times)), origin = structure(c(1,
 1, 1970), .Names = c(month, day, year)), class = c(chron,
 dates, times)), class = zoo)

 #Mullholland ER daytime

 # extract section of interest

 w - window(day, start=chron(6/8/6, 4:16:0), end=chron(6/8/6,
 20:31:0))

 # zap all interior points with NA's and then fill back in using linear
 interpolation

 lin - na.approx(replace(w, time(w)  start(w)  time(w)  end(w), NA))

 # plot the 3 series on one screen

 #plot(merge(day, w, lin), col = 1:3, screen = 1)

 #can I do this with a command (generically?)?  I am going to use this inside
 of a function, and below is the result I would like.
 d- c(day[1:17,], lin, day[84:96,])
 plot(d)

 --
 Let's not spend our time and resources thinking about things that are so
 little or so large that all they really do for us is puff us up and make us
 feel like gods. We are mammals, and have not exhausted the annoying little
 problems of being mammals.

 -K. Mullis

[[alternative HTML version deleted]]

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


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


Re: [R] Panel of pie charts

2008-07-29 Thread John Kane
Have a look at mfcol in ?par.


--- On Tue, 7/29/08, Amin Momin [EMAIL PROTECTED] wrote:

 From: Amin Momin [EMAIL PROTECTED]
 Subject: [R] Panel of pie charts
 To: r-help@r-project.org
 Received: Tuesday, July 29, 2008, 7:34 AM
 Hi ,
   I am looking to making a panel of pie charts fo some of
 my dritribution data . I was wondering if there is a way in
 any R package to write a small script to do so.
 
   Thanks,
Amin 
   [[alternative HTML version deleted]]
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained,
 reproducible code.


  __
[[elided Yahoo spam]]

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


Re: [R] combining zoo series with an overlapping index?

2008-07-29 Thread stephen sefick
when I plot this it gives me the coredata(lin) (Ithink).  I would like
d- c(day[1:17,], lin, day[84:96,])
to be the result


On Tue, Jul 29, 2008 at 9:18 AM, Gabor Grothendieck [EMAIL PROTECTED]
 wrote:

 Try:

 window(d, time(lin)) - coredata(lin)

 On Tue, Jul 29, 2008 at 9:01 AM, stephen sefick [EMAIL PROTECTED] wrote:
  day-structure(c(7.7, 7.7, 7.7, 7.7, 7.7, 7.71, 7.7, 7.71, 7.71, 7.7,
  7.7, 7.7, 7.7, 7.69, 7.68, 7.68, 7.67, 7.67, 7.67, 7.66, 7.65,
  7.65, 7.65, 7.64, 7.64, 7.63, 7.63, 7.63, 7.62, 7.62, 7.62, 7.62,
  7.63, 7.63, 7.63, 7.63, 7.63, 7.64, 7.64, 7.65, 7.65, 7.65, 7.66,
  7.66, 7.67, 7.67, 7.67, 7.68, 7.68, 7.69, 7.69, 7.69, 7.69, 7.7,
  7.7, 7.7, 7.7, 7.7, 7.71, 7.7, 7.7, 7.71, 7.71, 7.7, 7.7, 7.7,
  7.7, 7.7, 7.7, 7.7, 7.69, 7.68, 7.68, 7.68, 7.67, 7.67, 7.66,
  7.66, 7.66, 7.66, 7.65, 7.65, 7.65, 7.65, 7.66, 7.66, 7.66, 7.67,
  7.67, 7.68, 7.68, 7.68, 7.68, 7.69, 7.69, 7.69), index =
  structure(c(13307.000694,
  13307.01, 13307.021528, 13307.031944, 13307.042361,
  13307.052778, 13307.063194, 13307.073611, 13307.084028,
  13307.09, 13307.104861, 13307.115278, 13307.125694,
  13307.136111, 13307.146528, 13307.156944, 13307.167361,
  13307.18, 13307.188194, 13307.198611, 13307.209028,
  13307.219444, 13307.229861, 13307.240278, 13307.250694,
  13307.26, 13307.271528, 13307.281944, 13307.292361,
  13307.302778, 13307.313194, 13307.323611, 13307.334028,
  13307.34, 13307.354861, 13307.365278, 13307.375694,
  13307.386111, 13307.396528, 13307.406944, 13307.417361,
  13307.427778, 13307.438194, 13307.448611, 13307.459028,
  13307.469444, 13307.479861, 13307.490278, 13307.500694,
  13307.51, 13307.521528, 13307.531944, 13307.542361,
  13307.552778, 13307.563194, 13307.573611, 13307.584028,
  13307.59, 13307.604861, 13307.615278, 13307.625694,
  13307.636111, 13307.646528, 13307.656944, 13307.667361,
  13307.68, 13307.688194, 13307.698611, 13307.709028,
  13307.719444, 13307.729861, 13307.740278, 13307.750694,
  13307.76, 13307.771528, 13307.781944, 13307.792361,
  13307.802778, 13307.813194, 13307.823611, 13307.834028,
  13307.84, 13307.854861, 13307.865278, 13307.875694,
  13307.886111, 13307.896528, 13307.906944, 13307.917361,
  13307.927778, 13307.938194, 13307.948611, 13307.959028,
  13307.969444, 13307.979861, 13307.990278), format =
  structure(c(m/d/y,
  h:m:s), .Names = c(dates, times)), origin = structure(c(1,
  1, 1970), .Names = c(month, day, year)), class = c(chron,
  dates, times)), class = zoo)
 
  #Mullholland ER daytime
 
  # extract section of interest
 
  w - window(day, start=chron(6/8/6, 4:16:0), end=chron(6/8/6,
  20:31:0))
 
  # zap all interior points with NA's and then fill back in using linear
  interpolation
 
  lin - na.approx(replace(w, time(w)  start(w)  time(w)  end(w), NA))
 
  # plot the 3 series on one screen
 
  #plot(merge(day, w, lin), col = 1:3, screen = 1)
 
  #can I do this with a command (generically?)?  I am going to use this
 inside
  of a function, and below is the result I would like.
  d- c(day[1:17,], lin, day[84:96,])
  plot(d)
 
  --
  Let's not spend our time and resources thinking about things that are so
  little or so large that all they really do for us is puff us up and make
 us
  feel like gods. We are mammals, and have not exhausted the annoying
 little
  problems of being mammals.
 
  -K. Mullis
 
 [[alternative HTML version deleted]]
 
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
 




-- 
Let's not spend our time and resources thinking about things that are so
little or so large that all they really do for us is puff us up and make us
feel like gods. We are mammals, and have not exhausted the annoying little
problems of being mammals.

-K. Mullis

[[alternative HTML version deleted]]

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


Re: [R] writing the plots

2008-07-29 Thread John Kane
?png
?write.table


--- On Mon, 7/28/08, Rajasekaramya [EMAIL PROTECTED] wrote:

 From: Rajasekaramya [EMAIL PROTECTED]
 Subject: [R]  writing the plots
 To: r-help@r-project.org
 Received: Monday, July 28, 2008, 10:54 AM
 hi there,
 
 I want to write the plots in the pdfs and the details about
 the graph  in a
 seperate notepad.
 
 plot(as.numeric(lapply(resultgenes,length)),
 main=
 Geneset.gene#.bias.test,xlab=Top.Ranked.Genesets,
 ylab=gene.number.per.geneset)
 lines(loess.smooth(c(1:1000),as.numeric(lapply(resultgenes,length)),
  span =
 2/3, degree = 1,
 family = gaussian, evaluation = 50),col=2,
 lwd=3)
 
 
 i want this graph in the pdf and some notes regarding the
 graph in a
 notepad.
 
 Kindly help me
 
 Ramya
 
 -- 
 View this message in context:
 http://www.nabble.com/writing-the-plots-tp18689539p18689539.html
 Sent from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained,
 reproducible code.


  __
[[elided Yahoo spam]]

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


[R] Howto Draw Bimodal Gamma Curve with User Supplied Parameters

2008-07-29 Thread Gundala Viswanath
Hi,

Suppose I have the following vector (data points):
 x
  [1]  36.0  57.3  73.3  92.0 300.4  80.9  19.8  31.4  85.8  44.9  24.6  48.0
 [13]  28.0  38.3  85.2 103.6 154.4 128.5  38.3  72.4 122.7 123.1  41.8  21.7
 [25] 143.6 120.2  46.6  29.2  44.8  25.0  57.3  96.4  29.4  62.9  66.4  30.0
 [37]  24.1  14.8  56.6 102.4 117.5  90.4  37.2  79.6  27.8  17.1  26.6  16.3
 [49]  41.4  48.9  24.1  23.3   9.9  11.5  15.0  23.6  29.3  27.0  19.2  18.7
 [61]   4.1  13.0   3.3   5.5  38.2   8.5  39.6  39.2  16.1  35.3  23.3  31.5
 [73]  38.8  51.5  28.4  18.8  24.1  25.4  28.8  32.8  31.0  28.8  33.3  55.5
 [85]  39.2  21.0  43.7  16.3  50.6  34.6  66.3  50.5  59.4  46.7  51.9 125.6
 [97]  69.8  43.7  86.8  50.6 132.4  56.0   6.1   4.9   7.1   7.1  12.8  12.1
[109] 164.2  69.3  15.6  11.4  34.3   9.2  17.6  21.7  19.2  30.7  61.1  35.8
[121] 185.8 118.4  13.0   9.6  19.1  45.2  94.5 248.0  56.3  24.4  13.8  12.8
[133]  35.0  31.6  22.5  50.1  18.7  22.1  28.3  39.5  48.2  33.1  43.5  35.1
[145]  37.4  30.3  15.8  13.9  15.3  16.1  12.7  11.4  13.0  13.8  31.5  25.3
[157]  65.2  39.5


And the following parameter set (2 component) for gamma function.

 comp.1comp.2
alpha (shape)  2.855444  2.152056
beta  (scale)   10.418785 39.296224

These params are predefined/precalculated by user.

My question is how can I create a bimodal gamma curve - based on the
two parameter set -
on top of the histogram of data points above?

- Gundala Viswanath
Jakarta - Indonesia

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


[R] About clustering techniques

2008-07-29 Thread pacomet
Hello R users

It's some time I am playing with a dataset to do some cluster analysis. The
data set consists of 14 columns being geographical coordinates and monthly
temperatures in annual files

latitutde - longitude - temperature 1 -. - temperature 12

I have some missing values in some cases, maybe there are 8 monthly valid
values at some points with four non valid. I don't want to supress the whole
row with 8 good/4 bad values as I wanna try annual and monthy analysis.

I first tried kmeans but found a problem with missing values. When trying
without omitting missing values kmeans gives an error and when excluding
invalid data too many values are excluded in some years of the data series.

Now I have been reading about pam, pamk and clara, I think they can handle
missing values. But can't find out the way to perform the analysis with
these functions. As I'm not an statistics nor an R expert the fpc or cluster
package documentation is not enough for me. If you know about a website or a
tutorial explaining the way to use that functions, with examples to check if
possible, please post them.

Any other help or suggestion is greatly appreciated.

Thanks in advance

Paco

-- 
_
El ponent la mou, el llevant la plou
Usuari Linux registrat: 363952
---
Fotos: http://picasaweb.google.es/pacomet

[[alternative HTML version deleted]]

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


Re: [R] Coarsening the Resolution of a Dataset

2008-07-29 Thread Steve Murray

Unfortunately, when I get to the 'myCuts' line, I receive the following error:

Error: evaluation nested too deeply: infinite recursion / options(expressions=)?

...and I also receive warnings about memory allocation being reached (even 
though I've already used memory.limit() to maximise the memory) - this is a 
fairly sizeable dataset afterall, 2160 rows by 4320 columns.

Therefore I was wondering if there are any alternative ways of coarsening a 
dataset? Or are there any packages/commands built for this sort of thing? 

Any advice would be much appreciated!

Thanks again,

Steve


_
Find the best and worst places on the planet

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


Re: [R] About clustering techniques

2008-07-29 Thread ctu

Hi Paco,
I got the same problem with you before. Thus, I just impute the missing values
For example:

newdata-as.matrix(impute(olddata, fun=random))
then I believe that you could analyze your data.

Hopefully it helps.
Chunhao


Quoting pacomet [EMAIL PROTECTED]:


Hello R users

It's some time I am playing with a dataset to do some cluster analysis. The
data set consists of 14 columns being geographical coordinates and monthly
temperatures in annual files

latitutde - longitude - temperature 1 -. - temperature 12

I have some missing values in some cases, maybe there are 8 monthly valid
values at some points with four non valid. I don't want to supress the whole
row with 8 good/4 bad values as I wanna try annual and monthy analysis.

I first tried kmeans but found a problem with missing values. When trying
without omitting missing values kmeans gives an error and when excluding
invalid data too many values are excluded in some years of the data series.

Now I have been reading about pam, pamk and clara, I think they can handle
missing values. But can't find out the way to perform the analysis with
these functions. As I'm not an statistics nor an R expert the fpc or cluster
package documentation is not enough for me. If you know about a website or a
tutorial explaining the way to use that functions, with examples to check if
possible, please post them.

Any other help or suggestion is greatly appreciated.

Thanks in advance

Paco

--
_
El ponent la mou, el llevant la plou
Usuari Linux registrat: 363952
---
Fotos: http://picasaweb.google.es/pacomet

[[alternative HTML version deleted]]

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



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


Re: [R] About clustering techniques

2008-07-29 Thread Christian Hennig

Dear Paco,

in order to use the methods in the cluster package (including pam), look up 
the help page of daisy, which is able to compute dissimilarity matrices

handling missing values appropriately (in most situations).

A good reference is the Kaufman and Rousseeuw book cited on that help page.

Christian

On Tue, 29 Jul 2008, pacomet wrote:


Hello R users

It's some time I am playing with a dataset to do some cluster analysis. The
data set consists of 14 columns being geographical coordinates and monthly
temperatures in annual files

latitutde - longitude - temperature 1 -. - temperature 12

I have some missing values in some cases, maybe there are 8 monthly valid
values at some points with four non valid. I don't want to supress the whole
row with 8 good/4 bad values as I wanna try annual and monthy analysis.

I first tried kmeans but found a problem with missing values. When trying
without omitting missing values kmeans gives an error and when excluding
invalid data too many values are excluded in some years of the data series.

Now I have been reading about pam, pamk and clara, I think they can handle
missing values. But can't find out the way to perform the analysis with
these functions. As I'm not an statistics nor an R expert the fpc or cluster
package documentation is not enough for me. If you know about a website or a
tutorial explaining the way to use that functions, with examples to check if
possible, please post them.

Any other help or suggestion is greatly appreciated.

Thanks in advance

Paco

--
_
El ponent la mou, el llevant la plou
Usuari Linux registrat: 363952
---
Fotos: http://picasaweb.google.es/pacomet

[[alternative HTML version deleted]]

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



*** --- ***
Christian Hennig
University College London, Department of Statistical Science
Gower St., London WC1E 6BT, phone +44 207 679 1698
[EMAIL PROTECTED], www.homepages.ucl.ac.uk/~ucakche

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


Re: [R] combining zoo series with an overlapping index?

2008-07-29 Thread Gabor Grothendieck
The code I gave replaces that portion of d that
overlaps with lin with lin.  Is that not what you
wanted?  (please try to  use minimal examples
as shown here)

 library(zoo)
 d - zoo(1:10) + 100
 lin - - head(d, 4)
 window(d, time(lin)) - coredata(lin)
 d
   123456789   10
-101 -102 -103 -104  105  106  107  108  109  110


On Tue, Jul 29, 2008 at 9:48 AM, stephen sefick [EMAIL PROTECTED] wrote:
 when I plot this it gives me the coredata(lin) (Ithink).  I would like
 d- c(day[1:17,], lin, day[84:96,])
 to be the result


 On Tue, Jul 29, 2008 at 9:18 AM, Gabor Grothendieck
 [EMAIL PROTECTED] wrote:

 Try:

 window(d, time(lin)) - coredata(lin)

 On Tue, Jul 29, 2008 at 9:01 AM, stephen sefick [EMAIL PROTECTED] wrote:
  day-structure(c(7.7, 7.7, 7.7, 7.7, 7.7, 7.71, 7.7, 7.71, 7.71, 7.7,
  7.7, 7.7, 7.7, 7.69, 7.68, 7.68, 7.67, 7.67, 7.67, 7.66, 7.65,
  7.65, 7.65, 7.64, 7.64, 7.63, 7.63, 7.63, 7.62, 7.62, 7.62, 7.62,
  7.63, 7.63, 7.63, 7.63, 7.63, 7.64, 7.64, 7.65, 7.65, 7.65, 7.66,
  7.66, 7.67, 7.67, 7.67, 7.68, 7.68, 7.69, 7.69, 7.69, 7.69, 7.7,
  7.7, 7.7, 7.7, 7.7, 7.71, 7.7, 7.7, 7.71, 7.71, 7.7, 7.7, 7.7,
  7.7, 7.7, 7.7, 7.7, 7.69, 7.68, 7.68, 7.68, 7.67, 7.67, 7.66,
  7.66, 7.66, 7.66, 7.65, 7.65, 7.65, 7.65, 7.66, 7.66, 7.66, 7.67,
  7.67, 7.68, 7.68, 7.68, 7.68, 7.69, 7.69, 7.69), index =
  structure(c(13307.000694,
  13307.01, 13307.021528, 13307.031944, 13307.042361,
  13307.052778, 13307.063194, 13307.073611, 13307.084028,
  13307.09, 13307.104861, 13307.115278, 13307.125694,
  13307.136111, 13307.146528, 13307.156944, 13307.167361,
  13307.18, 13307.188194, 13307.198611, 13307.209028,
  13307.219444, 13307.229861, 13307.240278, 13307.250694,
  13307.26, 13307.271528, 13307.281944, 13307.292361,
  13307.302778, 13307.313194, 13307.323611, 13307.334028,
  13307.34, 13307.354861, 13307.365278, 13307.375694,
  13307.386111, 13307.396528, 13307.406944, 13307.417361,
  13307.427778, 13307.438194, 13307.448611, 13307.459028,
  13307.469444, 13307.479861, 13307.490278, 13307.500694,
  13307.51, 13307.521528, 13307.531944, 13307.542361,
  13307.552778, 13307.563194, 13307.573611, 13307.584028,
  13307.59, 13307.604861, 13307.615278, 13307.625694,
  13307.636111, 13307.646528, 13307.656944, 13307.667361,
  13307.68, 13307.688194, 13307.698611, 13307.709028,
  13307.719444, 13307.729861, 13307.740278, 13307.750694,
  13307.76, 13307.771528, 13307.781944, 13307.792361,
  13307.802778, 13307.813194, 13307.823611, 13307.834028,
  13307.84, 13307.854861, 13307.865278, 13307.875694,
  13307.886111, 13307.896528, 13307.906944, 13307.917361,
  13307.927778, 13307.938194, 13307.948611, 13307.959028,
  13307.969444, 13307.979861, 13307.990278), format =
  structure(c(m/d/y,
  h:m:s), .Names = c(dates, times)), origin = structure(c(1,
  1, 1970), .Names = c(month, day, year)), class = c(chron,
  dates, times)), class = zoo)
 
  #Mullholland ER daytime
 
  # extract section of interest
 
  w - window(day, start=chron(6/8/6, 4:16:0), end=chron(6/8/6,
  20:31:0))
 
  # zap all interior points with NA's and then fill back in using linear
  interpolation
 
  lin - na.approx(replace(w, time(w)  start(w)  time(w)  end(w), NA))
 
  # plot the 3 series on one screen
 
  #plot(merge(day, w, lin), col = 1:3, screen = 1)
 
  #can I do this with a command (generically?)?  I am going to use this
  inside
  of a function, and below is the result I would like.
  d- c(day[1:17,], lin, day[84:96,])
  plot(d)
 
  --
  Let's not spend our time and resources thinking about things that are so
  little or so large that all they really do for us is puff us up and make
  us
  feel like gods. We are mammals, and have not exhausted the annoying
  little
  problems of being mammals.
 
  -K. Mullis
 
 [[alternative HTML version deleted]]
 
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
  http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
 



 --
 Let's not spend our time and resources thinking about things that are so
 little or so large that all they really do for us is puff us up and make us
 feel like gods. We are mammals, and have not exhausted the annoying little
 problems of being mammals.

 -K. Mullis


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

Re: [R] About clustering techniques

2008-07-29 Thread Christian Hennig
A quick comment on this: imputation is an option to make things 
technically work, but it is not 
necessarily good. Imputation always introduces some noise, ie, it fakes 
information that is not really there.


Whether it is good depends strongly on the data, the situation and the 
imputation method (random often not being a very sensible 
choice).


Christian

On Tue, 29 Jul 2008, [EMAIL PROTECTED] wrote:


Hi Paco,
I got the same problem with you before. Thus, I just impute the missing 
values

For example:

newdata-as.matrix(impute(olddata, fun=random))
then I believe that you could analyze your data.

Hopefully it helps.
Chunhao


Quoting pacomet [EMAIL PROTECTED]:


Hello R users

It's some time I am playing with a dataset to do some cluster 
analysis. The
data set consists of 14 columns being geographical coordinates and 
monthly

temperatures in annual files

latitutde - longitude - temperature 1 -. - temperature 12

I have some missing values in some cases, maybe there are 8 monthly 
valid
values at some points with four non valid. I don't want to supress the 
whole
row with 8 good/4 bad values as I wanna try annual and monthy 
analysis.


I first tried kmeans but found a problem with missing values. When 
trying
without omitting missing values kmeans gives an error and when 
excluding
invalid data too many values are excluded in some years of the data 
series.


Now I have been reading about pam, pamk and clara, I think they can 
handle
missing values. But can't find out the way to perform the analysis 
with
these functions. As I'm not an statistics nor an R expert the fpc or 
cluster
package documentation is not enough for me. If you know about a 
website or a
tutorial explaining the way to use that functions, with examples to 
check if

possible, please post them.

Any other help or suggestion is greatly appreciated.

Thanks in advance

Paco

--
_
El ponent la mou, el llevant la plou
Usuari Linux registrat: 363952
---
Fotos: http://picasaweb.google.es/pacomet

[[alternative HTML version deleted]]

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

and provide commented, minimal, self-contained, reproducible code.



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

and provide commented, minimal, self-contained, reproducible code.


*** --- ***
Christian Hennig
University College London, Department of Statistical Science
Gower St., London WC1E 6BT, phone +44 207 679 1698
[EMAIL PROTECTED], www.homepages.ucl.ac.uk/~ucakche

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


[R] Most often pairs of chars across grouping variable

2008-07-29 Thread svga
Hi list,

is there a package or function to compute the frequencies of pairs of chars in 
a variable across a grouping variable? Eg:


d - data.frame(ID=gl(2,3), F=c(A,B,C,A,C,D))
 d
  ID F
1  1 A
2  1 B
3  1 C
4  2 A
5  2 C
6  2 D


Now I want to summarize the frequencies of all pairs A-B, A-C, A-D, B-C, B-D, 
C-D across ID:

   A B C D
A  - 1 2 1
B  - - 1 0
C  - - - 1


here, the combination A-C is most frequent. The real problem behind that is 
that 'F' codes diagnoses and I search for the most often pairs of diagnoses.

Thanks, Sven

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


Re: [R] Negative Binomial Regression

2008-07-29 Thread Prof Brian Ripley

On Tue, 29 Jul 2008, Ben Bolker wrote:


jcarmichael jcarmichael314 at gmail.com writes:




Hello.

I am attempting to duplicate a negative binomial regression in R.  SAS uses
generalized estimating equations for model fitting in the GENMOD procedure.

proc genmod data=mydata (where=(gender='F'));
by agegroup;
class id gender type;
model count = var1 var2 var3 /dist=NB link=log offset=lregtm;
repeated subject=id /type=exch;
run;

Since my dataset has several observations for each subject, I need the
REPEATED statement in order to indicate dependence among observations with
the same subject ID and independence amongst those with distinct subject
IDs.  The TYPE statement goes on to specify the structure of the correlation
matrix to be used (exchangeable in this case).


 I would try glmmPQL in the MASS package.  I don't think you
can *quite* get negative binomial regression this way, but
you can definitely get a quasipoisson model.  I think exchangeable
correlation corresponds to correlation=corCompSymm() in your
glmmPQL command.


The problem here is that GLMM and GEE are not fitting the same model -- in 
one the coefficients are subject-specific and in the other 
population-average (see MASS4 or Diggle, Liang, Zeger +/- Heagarty).


There are several R packages for GEE, including gee, yags, geepack.  The 
documentation of geeglm (geepack) claims it can be used with families as 
in glm(), so you could try it with MASS's negative.binomial family.


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

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


Re: [R] Most often pairs of chars across grouping variable

2008-07-29 Thread Marc Schwartz

on 07/29/2008 09:51 AM [EMAIL PROTECTED] wrote:

Hi list,

is there a package or function to compute the frequencies of pairs of
chars in a variable across a grouping variable? Eg:


d - data.frame(ID=gl(2,3), F=c(A,B,C,A,C,D))

d

ID F 1  1 A 2  1 B 3  1 C 4  2 A 5  2 C 6  2 D


Now I want to summarize the frequencies of all pairs A-B, A-C, A-D,
B-C, B-D, C-D across ID:

A B C D A  - 1 2 1 B  - - 1 0 C  - - - 1


here, the combination A-C is most frequent. The real problem behind
that is that 'F' codes diagnoses and I search for the most often
pairs of diagnoses.

Thanks, Sven


I suspect that there might be something over in Bioconductor, but here 
is one approach:


 table(data.frame(t(do.call(cbind,
 tapply(d$F, d$ID,
function(x) combn(as.character(x), 2))
   X2
X1  B C D
  A 1 2 1
  B 0 1 0
  C 0 0 1


See ?combn to create the initial pairs from the data. This is done on a 
per ID basis using tapply. The result is transposed into a data frame 
and then table() is used to create the cross tabulation of the results.


HTH,

Marc Schwartz

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


Re: [R] Case statements in R

2008-07-29 Thread Roland Rau

Dear all,

may I suggest to include this quotation of Patrick Burns in the fortunes 
package? :-)


Best,
Roland


Patrick Burns wrote:

A good reason to use '' rather than '' is if evaluating
whatever is on the right will create an error if what is on
the left is FALSE.  '' and '||' stop if they already know
the answer from seeing the left side.  So they are likely faster
by a few nanoseconds, but it's the ability to stop before you
go off a cliff that is probably the key feature.



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


[R] subscript out of bounds error.

2008-07-29 Thread Rajasekaramya

 (mean.std.s2n.loss.gain[[1]])
  Probe.Set.ID   rho_prime  r ho_prime_sdpom   Expr1   
matchinggenes Mean std_dev 
29 SNP_A-190846347 2.47 0.75  0 PRKCZ - 
  
0.34560.  1344 30 SNP_A-213137044 2.61 
0.58 0 PRKCZ   -0.32700.1517 
31 SNP_A-220300945 2.67  0.58 0 PRKCZ   
  
-0.3242.1555
32 SNP_A-209854643 2.81  0.58 1 PRKCZ   
 
-0.3212  0.155  
33 SNP_A-226456542 2.90  0.58 1 PRKCZ   
  
0.3400.1782
34 SNP_A-178872841 2.84 0.58 1 PRKCZ
  
0.339 .1778 
1398   SNP_A-1834906   25 1.90 0.42 0 FRAP1 
   
-0.2710   0.03668
1399   SNP_A-2146504   25 1.90 0.42 0 FRAP1 
   
-0.2700.03627
1400   SNP_A-230730925 1.90 0.42 0 FRAP1
  
-0.061   0.0820
1401   SNP_A-226165827 2.07 0.42 0 FRAP1
  
-0.276 0.04101 
1402   SNP_A-427448129 2.21 0.42 0 FRAP1
  
-0.277 0.0412 
1403   SNP_A-205307029 2.21 0.42 0 FRAP1
- 0.285 0.0501 
1404   SNP_A-427638229 2.21 0.42 0 FRAP1
-0.30060.0523



mean.std.s2n.loss.gain is a l;ist object of length 1000.for sample i have
provided with 1 element in the list.

I want to calculate the mean and std dev for the unique genes.I wrote a code
for that.but it is throwing me an error subscript out of bounds for
.unique.genes[[i]]



genes.info.unique - lapply(mean.std.s2n.loss.gain,function(.unique){
.unique.genes-unique(.unique[,6])
.unique.genes-as.vector(as.matrix(.unique.genes))
snpnumber-list()
mean.genes-list()
sd.genes-list()
class(.unique.genes)
for(i in 1:length(.unique.genes))
{
snpnumber[[i]]- length(.unique[.unique[,6] %in% .unique.genes[[i]],1])
mean.genes[[i]]-mean(.unique[.unique[,6] %in%  .unique.genes[[i]],7])
sd.genes[[i]]-sd(.unique[.unique[,6] %in% .unique.genes[[i]],7])
}
cbind(.unique.genes,unlist(snpnumber),unlist(mean.genes),unlist(sd.genes))
})


kindly help me

Ramya


-- 
View this message in context: 
http://www.nabble.com/subscript-out-of-bounds-error.-tp18714563p18714563.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] lattice question

2008-07-29 Thread Troels Ring
Dear friends - I have a plot of simulated data to compare to the 
observed and want the simulated to appear accumulated so that a darker 
grey corresponds to more lines of simulated data, but on top of that I 
want the measured values (Na) to be very visible. The code below meets 
only few of the desires. Hope someone will give me help to master the 
panel formation.

Best wishes
Troels

xyplot(Na+yrep[,1:100]~time|as.factor(ID),yrepbug,type=c(g,o),
span=0.75,ylab=[Na])

--

Troels Ring - -
Department of nephrology - - 
Aalborg Hospital 9100 Aalborg, Denmark - -

+45 99326629 - -
[EMAIL PROTECTED]

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


[R] subscript out of bounds help

2008-07-29 Thread Rajasekaramya

(mean.std.s2n.loss.gain[[1]])
Probe.Set.ID   rho_prime  rho_prime_sd  pom Expr1 matchinggenes   
Meanstd_dev 
29 SNP_A-190846347 2.47 0.75 0 PRKCZ
-0.345616170.13443676  
30 SNP_A-213137044 2.61 0.58 0 PRKCZ
-0.327046920.15171799 
31 SNP_A-220300945 2.67 0.58 0 PRKCZ
-0.324239380.15554016
32 SNP_A-209854643 2.81 0.58 1 PRKCZ
-0.321298220.15939414  
33 SNP_A-226456542 2.90 0.58 1 PRKCZ 
0.340551500.17828152   
34 SNP_A-178872841 2.84 0.58 1 PRKCZ 
0.339537440.17782981
35 SNP_A-184250934 2.59 0.58 1 PRKCZ
-0.297423980.14795647  
36 SNP_A-425038933 2.60 0.50 1 PRKCZ
-0.294941910.14550642
1398   SNP_A-183490625 1.900.42 0 FRAP1
-0.271049420.03668130  
1399   SNP_A-214650425 1.900.42 0 FRAP1
-0.270215920.03623132  
1400   SNP_A-230730925 1.900.42 0 FRAP1
-0.061870950.08202400  
1401   SNP_A-226165827 2.070.42 0 FRAP1
-0.276791000.04101169  
1402   SNP_A-427448129 2.210.42 0 FRAP1
-0.277134330.0412  
1403   SNP_A-205307029 2.210.42 0 FRAP1
-0.285844830.05010136  
1404   SNP_A-427638229 2.21   0.420 FRAP1
-0.300633420.05233933  


mean.std.s2n.loss.gain is a list object of length 1000.
I want to calculate the mean and std dev for the unique genes using the
column 7 values.I wrote a code but it showing me subscript out of bounds
error for .unique.genes[[i]].Could any body help me fixing the error.

genes.info.unique - lapply(mean.std.s2n.loss.gain,function(.unique){
.unique.genes-unique(.unique[,6])
.unique.genes-as.vector(as.matrix(.unique.genes))
snpnumber-list()
mean.genes-list()
sd.genes-list()
class(.unique.genes)
for(i in 1:length(.unique.genes))
{
snpnumber[[i]]- length(.unique[.unique[,6] %in% .unique.genes[[i]],1])
mean.genes[[i]]-mean(.unique[.unique[,6] %in%  .unique.genes[[i]],7])
sd.genes[[i]]-sd(.unique[.unique[,6] %in% .unique.genes[[i]],7])
}
cbind(.unique.genes,unlist(snpnumber),unlist(mean.genes),unlist(sd.genes))
})


kindly help me

Ramya

-- 
View this message in context: 
http://www.nabble.com/subscript-out-of-bounds-help-tp18714594p18714594.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] Removing script file

2008-07-29 Thread Dennis Fisher
Colleagues,

(Running R 2.7.0)

I have a script that I want to delete as it completes execution.  The  
penultimate line of the script (before the quit command) is:
file.remove(Scriptname)
The script is executed as:
R --no-save  Scriptname

In OS X and Linux this is successful and returns:
  file.remove(x)
 [1] TRUE
and the file is deleted

In Windows XP, it returns:
 [1] FALSE

and the file is not deleted.

I realize that this is an OS issue rather than an issue with execution  
of the command.  Is there some means to delete the script from within R?

Dennis

Dennis Fisher MD
P  (The P Less Than Company)
Phone: 1-866-PLessThan (1-866-753-7784)
Fax: 1-415-564-2220
www.PLessThan.com


[[alternative HTML version deleted]]

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


[R] tensor product of equi-spaced B-splines in the unit square

2008-07-29 Thread Patrizio Frederic
Dear all,
I need to compute tensor product of B-spline defined over equi-spaced
break-points.
I wrote my own program (it works in a 2-dimensional setting)

library(splines)
# set the break-points
Knots   = seq(-1,1,length=10)
# number of splines
M   = (length(Knots)-4)^2
# short cut to splineDesign function
bspline = function(x) splineDesign(Knots,x,outer.ok = T)
# bivariate tensor product of bspline
btens   = function(x) t(bspline(x[1]))%*%bspline(x[2])
# numebr of points to plot
ng  = 51
# create vectors for plotting
xgr = seq(-1,1,length=ng)
xgr2= expand.grid(xgr,xgr)
# generate random coef. of linear combination
bet   = rnorm(M)
# create matrix for contour-type plot
Bx  = apply(xgr2,1,btens)
Bmat= matrix(t(Bx)%*%bet,ng)
# plot the result
contour(xgr,xgr,Bmat)
persp(xgr,xgr,Bmat,theta=15)

any of you have a better idea (ie more efficient)?
Thanks in advance,

Patrizio Frederic

+-
| Patrizio Frederic
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-

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


Re: [R] Negative Binomial Regression

2008-07-29 Thread Ben Bolker

-BEGIN PGP SIGNED MESSAGE-
Hash: SHA1

Prof Brian Ripley wrote:
| On Tue, 29 Jul 2008, Ben Bolker wrote:
|
| jcarmichael jcarmichael314 at gmail.com writes:
|
|
|
| Hello.
|
| I am attempting to duplicate a negative binomial regression in R.
| SAS uses
| generalized estimating equations for model fitting in the GENMOD
| procedure.
|
| proc genmod data=mydata (where=(gender='F'));
| by agegroup;
| class id gender type;
| model count = var1 var2 var3 /dist=NB link=log offset=lregtm;
| repeated subject=id /type=exch;
| run;
|
| Since my dataset has several observations for each subject, I need the
| REPEATED statement in order to indicate dependence among observations
| with
| the same subject ID and independence amongst those with distinct subject
| IDs.  The TYPE statement goes on to specify the structure of the
| correlation
| matrix to be used (exchangeable in this case).
|
|  I would try glmmPQL in the MASS package.  I don't think you
| can *quite* get negative binomial regression this way, but
| you can definitely get a quasipoisson model.  I think exchangeable
| correlation corresponds to correlation=corCompSymm() in your
| glmmPQL command.
|
| The problem here is that GLMM and GEE are not fitting the same model --
| in one the coefficients are subject-specific and in the other
| population-average (see MASS4 or Diggle, Liang, Zeger +/- Heagarty).
|
| There are several R packages for GEE, including gee, yags, geepack.  The
| documentation of geeglm (geepack) claims it can be used with families as
| in glm(), so you could try it with MASS's negative.binomial family.
|

~  Point taken (although I guess I was pointing the original poster
to a way to do a reasonable analysis, not necessarily to duplicate
the SAS analysis as requested).  Will the negative.binomial family
really work for this, since it seems to require a fixed theta
(overdispersion) parameter?

~  If I very naively do the following:

library(geepack)
data(dietox)
mf2 - formula(Weight~Cu*Time+I(Time^2)+I(Time^3))
gee2 - geeglm(mf2, data=dietox, id=Pig,
~   family=poisson(identity),corstr=ar1)
library(MASS)
gee2 -
geeglm(mf2,data=dietox,id=Pig,family=negative.binomial(theta=100),corstr=ar1)

~  gives an error variance invalid --

~  so the whole thing would seem to take a bit of troubleshooting

~  (geeglm also gives warnings about non-integer Poisson values -- I
don't know why a Poisson link is being used in this example for
a non-integer Weight value ... ?)

~  Ben Bolker

-BEGIN PGP SIGNATURE-
Version: GnuPG v1.4.6 (GNU/Linux)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iD8DBQFIj0Euc5UpGjwzenMRAhs/AJ9EOu8RpJibN9LMLbI0p1tCTls5xQCfRLnf
6hkinpQu0PMQgh4f/Sl5BaM=
=y7CX
-END PGP SIGNATURE-

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


[R] Bootstraping GAMs for confidence intervales calculation

2008-07-29 Thread Luis Reino

Dear R-Users,


I am resending this message just to reminder my question regarding the 
calculation of a bootstrap confidence intervals for a GAM plot.


I am trying to apply a bootstrap to a GAM in order to calculate the 95% 
confidence intervals for a smooth curve obtained by the “plot.gam” 
function of the mgcv package. Nonetheless, I am getting some 
difficulties in transposing the results for the graphs.


I used the following commands in R, “mgcv” and “boot” packages:

* attach(bbvc_11Jul08)*

* model.boot-function(data,indices){*

*+ sub.data-data[indices,]*

*+ 
model-gam(asin(sqrt(Cov0_30))~s(age,k=4,bs=cr),family=gaussian,data=sub.data)*


*+ coef(model)}*

* gam.boot-boot(bbvc_11Jul08,model.boot,R=1000)*

* *

* gam.boot*

ORDINARY NONPARAMETRIC BOOTSTRAP

Call:

boot(data = bbvc_11Jul08, statistic = model.boot, R = 1000)

Bootstrap Statistics :

original bias std. error

t1* 0.40370112 0.0002762674 0.02017378

t2* 0.04715501 -0.0144132280 0.07648348

t3* 0.14590936 -0.0135820501 0.05161244

t4* 0.13520098 -0.0096405733 0.04470793

I obtained 4 indexes in the “bootstrap” output. I know that for a lm 
function, with an independent variable, t1 is the index of the intercept 
and t2 is the index of the variable, but for GAMs I don't find what each 
indexes really mean!


Moreover, to calculate the confidence intervals of 95% I proceed as follows:

* *

*boot.ci(gam.boot,conf=0.95) *

BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS

Based on 1000 bootstrap replicates

CALL :

boot.ci(boot.out = gam.boot, conf = 0.95)

Intervals :

Level Normal Basic Studentized

95% ( 0.3639, 0.4430 ) ( 0.3635, 0.4434 ) ( 0.3583, 0.4507 )

Level Percentile BCa

95% ( 0.3640, 0.4439 ) ( 0.3623, 0.4434 )

Calculations and Intervals on Original Scale

Warning message:

In sqrt(tv[, 2]) : NaNs produced

My first doubt is when I make boot.ci, which of the four should I 
choose? By default, R calculates index=1:min(2,length(boot.out$t0) but 
I don't know if this is the correct form to do it.


Finally, I don’t know how to introduce the 95% confidence intervals in 
the plot.gam.


I hope you can help me. Thanks in advance.

Luís Reino


--
Luis Reino, PhD Student
Centro de Estudos Florestais
Departamento de Engenharia Florestal
Instituto Superior de Agronomia
Universidade Técnica de Lisboa
1349-017 Lisboa
Portugal

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


[R] more environment questions

2008-07-29 Thread Edna Bell
Hi R Gurus:

Here is some code that I was experimenting with, please:

 f1 - function(x) {
+ e1 - new.env(parent=.GlobalEnv)
+ environment(e1)
+ print(environment())
+ return(mean(x))
+ }
 f1(1:15)
environment: 0x02525444
[1] 8


My question:  why isn't the environment within the function set to e1, please?

Thanks,
Edna Bell

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


Re: [R] Panel of pie charts

2008-07-29 Thread Richard . Cotton
   I am looking to making a panel of pie charts fo some of my 
 dritribution data . I was wondering if there is a way in any R 
 package to write a small script to do so.
pie() will do you a one-off pie chart, but there is no equivalent using 
grid/ lattice graphics.  You could write a panel.pie function to draw 
them, but be warned: pie charts are almost never the best option for 
displaying data.  (See 
http://www.perceptualedge.com/articles/visual_business_intelligence/save_the_pies_for_dessert.pdf
 
for example.)  Have you considered using barcharts instead?

Regards,
Richie.

Mathematical Sciences Unit
HSL



ATTENTION:

This message contains privileged and confidential inform...{{dropped:20}}

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


Re: [R] more environment questions

2008-07-29 Thread Gabor Grothendieck
e1 - ...
creates a new environment e1

environment(e1)
does nothing

print(environment(e1))
print environment e1

By the way, if you are doing a lot of manipulations of environments
you might want to look at the proto package which reframes the
whole thing in terms of object oriented programming.

On Tue, Jul 29, 2008 at 12:12 PM, Edna Bell [EMAIL PROTECTED] wrote:
 Hi R Gurus:

 Here is some code that I was experimenting with, please:

 f1 - function(x) {
 + e1 - new.env(parent=.GlobalEnv)
 + environment(e1)
 + print(environment())
 + return(mean(x))
 + }
 f1(1:15)
 environment: 0x02525444
 [1] 8


 My question:  why isn't the environment within the function set to e1, please?

 Thanks,
 Edna Bell

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


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


[R] Question regarding statisticians

2008-07-29 Thread Max

Hi Everyone,

I'm apologize for asking this in the R-general list, but I'm unsure of 
where else to ask this burning question of mine:


Where do statisticians talk on the internet about professional/ career 
developement/ issues? I've found many a list (much like this one) that 
specailize in talking about statistics, but not about *being* a 
statistician. Any pointers to a message board or mailing list would be 
very helpful.


The closest I've found is a mailing list of statisticians in the UK, 
but it seems to be spammed with job offers from the UK and the odd 
question about an area of statistics. The closest thing I've found, 
which is for actuaries are actuary.com or 
http://www.actuarialoutpost.com/actuarial_discussion_forum


Thanks in Advance,

-Max

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


[R] table questions

2008-07-29 Thread Edna Bell
Hi again!


Suppose I have the following:

 xy - round(rexp(20),1)
 xy
 [1] 0.1 3.4 1.6 0.4 1.0 1.4 0.2 0.3 1.6 0.2 0.0 0.1 0.1 1.0 2.0 0.9
2.5 0.1 1.5 0.4
 table(xy)
xy
  0 0.1 0.2 0.3 0.4 0.9   1 1.4 1.5 1.6   2 2.5 3.4
  1   4   2   1   2   1   2   1   1   2   1   1   1

Is there a way to set things up to have
0 - 0.4   0.5 - 0.9  etc. please?

I know there is the cut functions, but breaks are required.  If you
don't have breaks, what should you do, please?

Would using the breaks from the hist function work appropriately, please?

thanks
Edna Bell

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


Re: [R] more environment questions

2008-07-29 Thread Edna Bell
Is there a way to set the environment within a function,,  please?


On Tue, Jul 29, 2008 at 11:25 AM, Gabor Grothendieck
[EMAIL PROTECTED] wrote:
 e1 - ...
 creates a new environment e1

 environment(e1)
 does nothing

 print(environment(e1))
 print environment e1

 By the way, if you are doing a lot of manipulations of environments
 you might want to look at the proto package which reframes the
 whole thing in terms of object oriented programming.

 On Tue, Jul 29, 2008 at 12:12 PM, Edna Bell [EMAIL PROTECTED] wrote:
 Hi R Gurus:

 Here is some code that I was experimenting with, please:

 f1 - function(x) {
 + e1 - new.env(parent=.GlobalEnv)
 + environment(e1)
 + print(environment())
 + return(mean(x))
 + }
 f1(1:15)
 environment: 0x02525444
 [1] 8


 My question:  why isn't the environment within the function set to e1, 
 please?

 Thanks,
 Edna Bell

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



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


Re: [R] more environment questions

2008-07-29 Thread Gabor Grothendieck
No but look at proto since I suspect the creation
of proto objects is basically what you are trying
to do through the back door.  Home page:

http://r-proto.googlecode.com

On Tue, Jul 29, 2008 at 12:29 PM, Edna Bell [EMAIL PROTECTED] wrote:
 Is there a way to set the environment within a function,,  please?


 On Tue, Jul 29, 2008 at 11:25 AM, Gabor Grothendieck
 [EMAIL PROTECTED] wrote:
 e1 - ...
 creates a new environment e1

 environment(e1)
 does nothing

 print(environment(e1))
 print environment e1

 By the way, if you are doing a lot of manipulations of environments
 you might want to look at the proto package which reframes the
 whole thing in terms of object oriented programming.

 On Tue, Jul 29, 2008 at 12:12 PM, Edna Bell [EMAIL PROTECTED] wrote:
 Hi R Gurus:

 Here is some code that I was experimenting with, please:

 f1 - function(x) {
 + e1 - new.env(parent=.GlobalEnv)
 + environment(e1)
 + print(environment())
 + return(mean(x))
 + }
 f1(1:15)
 environment: 0x02525444
 [1] 8


 My question:  why isn't the environment within the function set to e1, 
 please?

 Thanks,
 Edna Bell

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




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


[R] Graphics function question

2008-07-29 Thread Peter Flom
Hello

I have created a graph using the following commands:


startBReP3O1T - diffs$BReP3O1T - diffs$diff_BReP3O1T
endBReP3O1T - diffs$BReP3O1T

x - seq(47,89, length = 10)
ymin - min(min(startBReP3O1T), min(endBReP3O1T))
ymax - max(max(startBReP3O1T), max(endBReP3O1T))
y - seq(ymin, ymax, length = 10)
plot(x,y, type = 'n', xlab = 'Age', ylab = 'BReP3O1T', main = 'Age, decline and 
BReP3O1T')
segments(x0 = startage, x1 = endage, y0 = startBReP3O1T, y1 = endBReP3O1T,   
col = decline)
legend('topleft', legend = c('Stable', 'Decline'), lty = 1, col = c(1,2))


I would like to make this into a function.  The only thing that changes is 
BReP3O1T.
The problem is that sometimes this occurs as text (e.g. in ylab and main) 
sometimes after a $ (e.g. in the
first two assigment statements), and sometimes it refers to the values (e.g. in 
ymin and ymax)

Any help appreciated.

Peter




Peter L. Flom, PhD
Brainscope, Inc.
212 263 7863 (MTW)
212 845 4485 (Th)
917 488 7176 (F)



[[alternative HTML version deleted]]

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


Re: [R] combining zoo series with an overlapping index?

2008-07-29 Thread stephen sefick
plot(replace(day, is.list(index(lin)), coredata(lin)))

[[alternative HTML version deleted]]

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


[R] try question

2008-07-29 Thread Edna Bell
Hi still yet again!

I have the following code:

 try(log(rnorm(25)),silent=TRUE)
 [1] -0.26396185 NaN NaN -0.13078069 -2.44997193
-2.15603971 NaN  0.94917495  0.07244544 NaN
[11] -1.06341127 -0.42293099 -0.53769569  0.95134763  0.93403340
  NaN -0.10502078 NaN  0.30283262 NaN
[21] -0.11696872 -3.84122332 NaN NaN -0.12808690
Warning message:
In log(rnorm(25)) : NaNs produced


I thought that putting the silent = TRUE would suppress the
warnings, please.  What should I do instead, please?

Thanks,
Edna Bell

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


[R] accessing list elements

2008-07-29 Thread Paul Adams

Hello everyone,
I have a list which I am trying to calculate a max value.I have the list as 
w-c(v[[1]][1],...v[[100]][1]). The problem I am getting is that the function 
max is saying the list is an invalid
type (list) of  argument.When I show the element v[[1]][1] it shows as 
$statistic,V,736.I am only wanting to use the number 736 from v[[1]][1] but am 
not sure how to access that number only?I believe if I just use the number then 
I should be able to calculate the max.
 
Any help would be appreciated
Paul
 


  
[[alternative HTML version deleted]]

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


Re: [R] try question

2008-07-29 Thread Richard . Cotton
  try(log(rnorm(25)),silent=TRUE)
  [1] -0.26396185 NaN NaN -0.13078069 -2.44997193
 -2.15603971 NaN  0.94917495  0.07244544 NaN
 [11] -1.06341127 -0.42293099 -0.53769569  0.95134763  0.93403340
   NaN -0.10502078 NaN  0.30283262 NaN
 [21] -0.11696872 -3.84122332 NaN NaN -0.12808690
 Warning message:
 In log(rnorm(25)) : NaNs produced
 
 
 I thought that putting the silent = TRUE would suppress the
 warnings, please.  What should I do instead, please?
The 'silent' argument to try supresses errors, not warnings.  You can 
convert warnings to errors with options(warn=2), or ignore warnings with 
options(warn=-1).  See ?options for more details.

Regards,
Richie.

Mathematical Sciences Unit
HSL



ATTENTION:

This message contains privileged and confidential inform...{{dropped:20}}

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


Re: [R] Removing script file

2008-07-29 Thread Prof Brian Ripley

On Tue, 29 Jul 2008, Dennis Fisher wrote:


Colleagues,

(Running R 2.7.0)

I have a script that I want to delete as it completes execution.  The
penultimate line of the script (before the quit command) is:
file.remove(Scriptname)
The script is executed as:
R --no-save  Scriptname

In OS X and Linux this is successful and returns:

file.remove(x)

[1] TRUE

and the file is deleted

In Windows XP, it returns:

[1] FALSE


and the file is not deleted.

I realize that this is an OS issue rather than an issue with execution
of the command.  Is there some means to delete the script from within R?


Not in an OS that does not allow open files to be deleted.  (Not the whole 
truth for WIndows, but true of the mechanism used to implement redirection 
in the compiler used for R.)


You could modify Rscript to do so (it does not use re-direction).



Dennis

Dennis Fisher MD
P  (The P Less Than Company)
Phone: 1-866-PLessThan (1-866-753-7784)
Fax: 1-415-564-2220
www.PLessThan.com


[[alternative HTML version deleted]]

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



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

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


Re: [R] Panel of pie charts

2008-07-29 Thread Gavin Simpson
On Tue, 2008-07-29 at 17:23 +0100, [EMAIL PROTECTED] wrote:
I am looking to making a panel of pie charts fo some of my 
  dritribution data . I was wondering if there is a way in any R 
  package to write a small script to do so.
 pie() will do you a one-off pie chart, but there is no equivalent using 
 grid/ lattice graphics.  You could write a panel.pie function to draw 
 them, but be warned: pie charts are almost never the best option for 
 displaying data.  (See 
 http://www.perceptualedge.com/articles/visual_business_intelligence/save_the_pies_for_dessert.pdf
  
 for example.)  Have you considered using barcharts instead?

All concerns about pie charts aside, Deepayan has already done this for
you, as part of his excellent book, Lattice: Multivariate data
Visualization with R (2008) Springer, part of the UseR series.

Check out the book's website for the code to reproduce all the figures:

http://lmdvr.r-forge.r-project.org/figures/figures.html

in particular you should look at chapter 14, figure 14.5 and the
associated code, including panel functions, to reproduce this figure.

HTH

G
-- 
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
 Dr. Gavin Simpson [t] +44 (0)20 7679 0522
 ECRC, UCL Geography,  [f] +44 (0)20 7679 0565
 Pearson Building, [e] gavin.simpsonATNOSPAMucl.ac.uk
 Gower Street, London  [w] http://www.ucl.ac.uk/~ucfagls/
 UK. WC1E 6BT. [w] http://www.freshwaters.org.uk
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%

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


Re: [R] Graphics function question

2008-07-29 Thread Sarah Goslee
You can do it very easily using subsetting and a bit of paste().

I assumed you didn't need startdata and enddata  after the plot
has been made, but if you do you can change the last line of
the function to return them.

Sarah

myfunction - function(dataname)
{
# where dataname is a string, eg myfunction(BReP3O1T)

startdata - diffs[, dataname] - diffs[, paste(diff, dataname, sep=_)]
enddata - diffs[, dataname[

x - seq(47,89, length = 10)
ymin - min(min(startdata), min(enddata))
 ymax - max(max(startdata), max(enddata))
y - seq(ymin, ymax, length = 10)
plot(x,y, type = 'n', xlab = 'Age', ylab = dataname, main =
paste(Age, decline and, dataname)
segments(x0 = startage, x1 = endage, y0 = startdata, y1 = enddata,
col = decline)
legend('topleft', legend = c('Stable', 'Decline'), lty = 1, col = c(1,2))
invisible()

}

On Tue, Jul 29, 2008 at 12:38 PM, Peter Flom [EMAIL PROTECTED] wrote:
 Hello

 I have created a graph using the following commands:

 
 startBReP3O1T - diffs$BReP3O1T - diffs$diff_BReP3O1T
 endBReP3O1T - diffs$BReP3O1T

 x - seq(47,89, length = 10)
 ymin - min(min(startBReP3O1T), min(endBReP3O1T))
 ymax - max(max(startBReP3O1T), max(endBReP3O1T))
 y - seq(ymin, ymax, length = 10)
 plot(x,y, type = 'n', xlab = 'Age', ylab = 'BReP3O1T', main = 'Age, decline 
 and BReP3O1T')
 segments(x0 = startage, x1 = endage, y0 = startBReP3O1T, y1 = endBReP3O1T,   
 col = decline)
 legend('topleft', legend = c('Stable', 'Decline'), lty = 1, col = c(1,2))


 I would like to make this into a function.  The only thing that changes is 
 BReP3O1T.
 The problem is that sometimes this occurs as text (e.g. in ylab and main) 
 sometimes after a $ (e.g. in the
 first two assigment statements), and sometimes it refers to the values (e.g. 
 in ymin and ymax)

 Any help appreciated.

 Peter




-- 
Sarah Goslee
http://www.functionaldiversity.org

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


Re: [R] accessing list elements

2008-07-29 Thread Sarah Goslee
What's  v? And w? And what exactly do you want to do?

A small reproducible example would be very helpful. Without
knowing what your data look like, it's hard to make helpful
suggestions.

Sarah

On Tue, Jul 29, 2008 at 12:56 PM, Paul Adams [EMAIL PROTECTED] wrote:

 Hello everyone,
 I have a list which I am trying to calculate a max value.I have the list as 
 w-c(v[[1]][1],...v[[100]][1]). The problem I am getting is that the function 
 max is saying the list is an invalid
 type (list) of  argument.When I show the element v[[1]][1] it shows as 
 $statistic,V,736.I am only wanting to use the number 736 from v[[1]][1] but 
 am not sure how to access that number only?I believe if I just use the number 
 then I should be able to calculate the max.

 Any help would be appreciated
 Paul


-- 
Sarah Goslee
http://www.functionaldiversity.org

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


Re: [R] Graphics function question

2008-07-29 Thread Marc Schwartz

on 07/29/2008 11:38 AM Peter Flom wrote:

Hello

I have created a graph using the following commands:


startBReP3O1T - diffs$BReP3O1T - diffs$diff_BReP3O1T
endBReP3O1T - diffs$BReP3O1T

x - seq(47,89, length = 10)
ymin - min(min(startBReP3O1T), min(endBReP3O1T))
ymax - max(max(startBReP3O1T), max(endBReP3O1T))
y - seq(ymin, ymax, length = 10)
plot(x,y, type = 'n', xlab = 'Age', ylab = 'BReP3O1T', main = 'Age, decline and 
BReP3O1T')
segments(x0 = startage, x1 = endage, y0 = startBReP3O1T, y1 = endBReP3O1T,   
col = decline)
legend('topleft', legend = c('Stable', 'Decline'), lty = 1, col = c(1,2))

I would like to make this into a function.  The only thing that changes is 
BReP3O1T.
The problem is that sometimes this occurs as text (e.g. in ylab and main) 
sometimes after a $ (e.g. in the
first two assigment statements), and sometimes it refers to the values (e.g. in 
ymin and ymax)

Any help appreciated.

Peter


Hey Peter, LTNS!  Job change it looks like from the e-mail address?

Here is one approach. It is not clear from the above, where 'startage', 
'endage' and 'decline' come from, so I am passing them as arguments:


In your example above, the function call would be:

MyPlot(diffs$BReP3O1T, diffs$diff_BReP3O1T, startage, endage, decline)


MyPlot - function(x, y, startage, endage, decline)
{
  start - x - y
  end - x

  s1 - seq(47, 89, length = 10)

  # Don't need to use min/max on each vector separately
  ymin - min(start, end, na.rm = TRUE)
  ymax - max(start, end, na.rm = TRUE)

  s2 - seq(ymin, ymax, length = 10)

  # Turn 'x' into a label, stripping anything before $ if present
  label - gsub(^.+\\$, , deparse(substitute(x)))

  plot(s1, s2, type = n,  xlab = Age,  ylab = label,
   main = paste(Age, decline and, label)

  segments(startage, endage, start, end,  col = decline)

  legend(topleft, legend = c(Stable, Decline), lty = 1,
 col = c(1, 2))
}



HTH,

Marc Schwartz

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


Re: [R] Negative Binomial Regression

2008-07-29 Thread Prof Brian Ripley

On Tue, 29 Jul 2008, Ben Bolker wrote:


-BEGIN PGP SIGNED MESSAGE-
Hash: SHA1

Prof Brian Ripley wrote:
| On Tue, 29 Jul 2008, Ben Bolker wrote:
|
| jcarmichael jcarmichael314 at gmail.com writes:
|
| Hello.
|
| I am attempting to duplicate a negative binomial regression in R.
| SAS uses
| generalized estimating equations for model fitting in the GENMOD
| procedure.
|
| proc genmod data=mydata (where=(gender='F'));
| by agegroup;
| class id gender type;
| model count = var1 var2 var3 /dist=NB link=log offset=lregtm;
| repeated subject=id /type=exch;
| run;
|
| Since my dataset has several observations for each subject, I need the
| REPEATED statement in order to indicate dependence among observations
| with
| the same subject ID and independence amongst those with distinct subject
| IDs.  The TYPE statement goes on to specify the structure of the
| correlation
| matrix to be used (exchangeable in this case).
|
|  I would try glmmPQL in the MASS package.  I don't think you
| can *quite* get negative binomial regression this way, but
| you can definitely get a quasipoisson model.  I think exchangeable
| correlation corresponds to correlation=corCompSymm() in your
| glmmPQL command.
|
| The problem here is that GLMM and GEE are not fitting the same model --
| in one the coefficients are subject-specific and in the other
| population-average (see MASS4 or Diggle, Liang, Zeger +/- Heagarty).
|
| There are several R packages for GEE, including gee, yags, geepack.  The
| documentation of geeglm (geepack) claims it can be used with families as
| in glm(), so you could try it with MASS's negative.binomial family.
|

~  Point taken (although I guess I was pointing the original poster
to a way to do a reasonable analysis, not necessarily to duplicate
the SAS analysis as requested).  Will the negative.binomial family
really work for this, since it seems to require a fixed theta
(overdispersion) parameter?


I was answering 'I am trying to duplicate': I don't know how SAS estimates 
the parameter by GEE. The best guess I have is that theta is estimated 
from the initial glm fit and fixed in the GEE phase, but that is only 
interpolation from very vague descriptions.



~  If I very naively do the following:

library(geepack)
data(dietox)
mf2 - formula(Weight~Cu*Time+I(Time^2)+I(Time^3))
gee2 - geeglm(mf2, data=dietox, id=Pig,
~   family=poisson(identity),corstr=ar1)
library(MASS)
gee2 -
geeglm(mf2,data=dietox,id=Pig,family=negative.binomial(theta=100),corstr=ar1)

~  gives an error variance invalid --

~  so the whole thing would seem to take a bit of troubleshooting


I wasn't placing much faith on geeglm actually being as general as it says 
it is ('claims' ... 'could try').



~  (geeglm also gives warnings about non-integer Poisson values -- I
don't know why a Poisson link is being used in this example for
a non-integer Weight value ... ?)


'poisson' _family_, I presume?

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

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


Re: [R] try question

2008-07-29 Thread Prof Brian Ripley

On Tue, 29 Jul 2008, Edna Bell wrote:


Hi still yet again!

I have the following code:


try(log(rnorm(25)),silent=TRUE)

[1] -0.26396185 NaN NaN -0.13078069 -2.44997193
-2.15603971 NaN  0.94917495  0.07244544 NaN
[11] -1.06341127 -0.42293099 -0.53769569  0.95134763  0.93403340
 NaN -0.10502078 NaN  0.30283262 NaN
[21] -0.11696872 -3.84122332 NaN NaN -0.12808690
Warning message:
In log(rnorm(25)) : NaNs produced




I thought that putting the silent = TRUE would suppress the
warnings, please.  What should I do instead, please?


No, it supresses the error message (and you have no error here).  Use 
suppressWarnings().


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

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


Re: [R] accessing list elements

2008-07-29 Thread Erik Iverson

Hello -

Paul Adams wrote:

Hello everyone, I have a list which I am trying to calculate a max
value.I have the list as w-c(v[[1]][1],...v[[100]][1]). The problem
I am getting is that the function max is saying the list is an
invalid type (list) of  argument.When I show the element v[[1]][1]
it shows as $statistic,V,736.I am only wanting to use the number 736
from v[[1]][1] but am not sure how to access that number only?I
believe if I just use the number then I should be able to calculate
the max.

Any help would be appreciated Paul


Please see footer of this message, which states ...provide commented, 
minimal, self-contained, reproducible code.


You do not do this.

Is your problem like this?

a - list(20, 30)
max(a) ##error
max(unlist(a)) ##no error







[[alternative HTML version deleted]]






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


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


Re: [R] Negative Binomial Regression

2008-07-29 Thread Ben Bolker

Prof Brian Ripley wrote:


'poisson' _family_, I presume?



  oops, yes.

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


Re: [R] Major Cloud Computing Testbed Announced at University of Illinois

2008-07-29 Thread Ajay ohri
*R we Wishing  on a cloud-*

I would like to Test R as an cloud analytical tool.

 remote submit csv files to the web for data mining ,
get data crunching with scalable processing power and disk storage ,
display Silverlight graphs of summarization,
downloadable in pdf's.

 data mining on demand, on function.

 rather than pay costly annual analytical softwares license fees.

analytics softwares charge by cpu count on terminal servers.

r is free,open source and has the largest pool of statisticians and advanced
packages available (like Linux of data mining) www.r-project.org
Advanced stats packages limited in say weather forecasting etc because
parallel processing /super-computing PC 's are fragmented.

would do it myself, but i am the little guy.


regards,

Ajay Ohri

www.decisionstats.com
http://www.linkedin.com/in/ajohri
(Networking openly.)

(cc-ed to R help list)



On Tue, Jul 29, 2008 at 10:46 PM, Suman Chaudhuri 
[EMAIL PROTECTED] wrote:

 With our partnership with Amazon, we at Patni (www.patni.com) are thinking
 of doing a prototype test suite for building apps in the cloud that would
 address something similar to what Reuven was talking about.

 I would be interested in hearing your thoughts on what features this test
 suite in the cloud can provide to anyone who would be interested in using
 this service). If anyone here has thoughts that they would like to share,
 please do so.

 We were thinking more so from an application/website perspective
 (acceptance testing, unit testing, distributed functional testing, testing
 Ruby-on-Rails, etc), but if there are other needs, like testing out the
 actual platform, maybe we can have some discussions around that.

 Suman



  Date: Tue, 29 Jul 2008 10:06:38 -0700
  Subject: Re: Major Cloud Computing Testbed Announced at University of
 Illinois
  From: [EMAIL PROTECTED]
  To: [EMAIL PROTECTED]
 
 
   Something like this (i haven't used it, only found it via google
  search, i'm sure there are others):
 
 http://it.tmcnet.com/topics/it/articles/33908-alertsite-debuts-new-web-site-load-testing-solution.htm
 
  Seems a not very difficult task to build some EC2 service that will
  add instances to hammer your site until some specified criteria are
  met like successful requests/second drops below a threshold or 500
  errors crop up. Anyone doing that yet?
 
  - randall
 
 
 
  On Jul 29, 9:49 am, Reuven Cohen [EMAIL PROTECTED] wrote:
   We've been working with Intel for the last 18 months on several
   cloud projects and our reoccurring problem has been how to do you
   test a cloud platform designed to manage upwards of 100,000 servers
   without access to a real world test bed. This is a good first step,
   but with less then 1000 servers it still doesn't adequately address
   the problem.
  
   I wonder how others developing large scale cloud platforms are
   handling scalability testing?
  
   Reuven
  
  
  
   On Tue, Jul 29, 2008 at 12:06 PM, Ameed Taylor [EMAIL PROTECTED]
 wrote:
  
The University of Illinois announced a Cloud Computing Research
Initiative with Yahoo, HP and Intel today. View the press release
here.http://tinyurl.com/6z5cpo
  
The fact that this testbed will be Global in nature makes it notable
in my opinion. Sites will include The University of Illinois, the
National Science Foundation, the Infocomm Development Authority of
Singapore (IDA), and the Karlsruhe Institute of Technology in
Germany.
  
The press release mentions the testbed will include a 1,024-core HP
system with 200 TB of disk space and will run Apache Hadoop and the
Pig parallel programming language.
  
Hopefully someone from this group will submit a proposal to take part
in the Testbed as this highly visible project could help dispel any
lingering doubts about the viability of Cloud Computing especially as
it relates to extremely data intensive requirements.
  
Ameed Taylor
President
Applation LLC
   www.applation.com
  
Blog -www.ondemandbeat.com
  
   --
   --
  
   Reuven Cohen
   Founder  Chief Technologist, Enomaly 
   Inc.www.enomaly.comhttp://inc.www.enomaly.com/::
 416 848 6036 x 1
   skype: ruv.net // aol: ruv6
  
   blog www.elasticvapor.com
   -
   Get Linked inhttp://linkedin.com/pub/0/b72/7b4
 
  --~--~-~--~~~---~--~~
 You received this message because you are subscribed to the Google
 Groups Cloud Computing group.
 To post to this group, send email to [EMAIL PROTECTED]
 To unsubscribe from this group, send email to
 [EMAIL PROTECTED]
 To post job listing, send email to [EMAIL PROTECTED] (position title,
 employer and location in subject, description in message body) or visit
 http://www.cloudjobs.net
 For more options, visit this group at
 http://groups.google.ca/group/cloud-computing?hl=en?hl=en
 Posting guidelines:
 http://groups.google.ca/group/cloud-computing/web/frequently-asked-qu...
 To invite your colleague
 

Re: [R] Graphics function question

2008-07-29 Thread Sarah Goslee
One clarification:

I did the dirty and completely unprofessional thing, and assumed that this
function would only be run in the workspace with all your data. I use functions
like this frequently for something I need to do multiple times for a particular
project, but will never do anywhere else, but it isn't proper programming, and
will fail if it can't find diffs and the other items referenced.

If you need this to be portable, you will need to pass them all as arguments
to the function.

Sarah
-- 
Sarah Goslee
http://www.functionaldiversity.org

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


[R] List names help

2008-07-29 Thread Rajasekaramya

hi ,

Is there anyways to delet the list names once created.

i tried
rm(names(mylist))

i didnt work

kindly help me

Ramya
-- 
View this message in context: 
http://www.nabble.com/List-names-help-tp18717742p18717742.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] biplot_group_colours_and_point_symbols

2008-07-29 Thread Indermaur Lukas
hi,
i am seeking for a solution to create a biplot that shows: 
a) data points instead of labels and that
b) shows two different groups (field: site with two factor levels: 
forest/corridor) in different colours
 
i tried to include pch=19 to show points instead of labels in the syntax 
below but failed.
 
biplot(model, main=, iter=0.1,  cex=c(0.3,1.2), col=c(1,2))

i appreciate any hint.
 
best regards,
lukas
 
 
 
°°° 
Lukas Indermaur, PhD student 
eawag / Swiss Federal Institute of Aquatic Science and Technology 
ECO - Department of Aquatic Ecology
Überlandstrasse 133
CH-8600 Dübendorf
Switzerland
 
Phone: +41 (0) 71 220 38 25
Fax: +41 (0) 44 823 53 15 
Email: [EMAIL PROTECTED]
www.lukasindermaur.ch http://www.lukasindermaur.ch/ 
 
 
 

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


[R] Flexible semivariogram in R (HELP)

2008-07-29 Thread Alessandro
Hey All,

 

I am a PhD student in Forestry science and I wish use kriging to spatial
LiDAR cloud points.  The shape of my semivariogram is like hole (wave)
model, but unfortunately I don't know the code of flexible model or hole
model to fit my semivariogram in R. I Found it in google or in R help but I
did't find nothing about this model in R. 

 

Thank you for Help

 

Alessandro


[[alternative HTML version deleted]]

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


Re: [R] List names help

2008-07-29 Thread Sarah Goslee
Sure, just set them to NULL.

 mylist - c(a=1, b=2, c=3)
 mylist
a b c
1 2 3
 names(mylist) - NULL
 mylist
[1] 1 2 3


Sarah

On Tue, Jul 29, 2008 at 1:45 PM, Rajasekaramya [EMAIL PROTECTED] wrote:

 hi ,

 Is there anyways to delet the list names once created.

 i tried
 rm(names(mylist))

 i didnt work

 kindly help me

 Ramya



-- 
Sarah Goslee
http://www.functionaldiversity.org

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


Re: [R] table questions

2008-07-29 Thread Patrizio Frederic
# is that what you want?

table(cut(xy,seq(0,max(xy)+.4,by=.4)))

# or this

table(cut(xy,hist(xy)$breaks)) # not the same


regards,

PF

+-
| Patrizio Frederic
| Research associate in Statistics,
| Department of Economics,
| University of Modena and Reggio Emilia,
| Via Berengario 51,
| 41100 Modena, Italy
|
| tel:  +39 059 205 6727
| fax:  +39 059 205 6947
| mail: [EMAIL PROTECTED]
+-



2008/7/29 Edna Bell [EMAIL PROTECTED]:
 Hi again!


 Suppose I have the following:

 xy - round(rexp(20),1)
 xy
  [1] 0.1 3.4 1.6 0.4 1.0 1.4 0.2 0.3 1.6 0.2 0.0 0.1 0.1 1.0 2.0 0.9
 2.5 0.1 1.5 0.4
 table(xy)
 xy
  0 0.1 0.2 0.3 0.4 0.9   1 1.4 1.5 1.6   2 2.5 3.4
  1   4   2   1   2   1   2   1   1   2   1   1   1

 Is there a way to set things up to have
 0 - 0.4   0.5 - 0.9  etc. please?

 I know there is the cut functions, but breaks are required.  If you
 don't have breaks, what should you do, please?

 Would using the breaks from the hist function work appropriately, please?

 thanks
 Edna Bell

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


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


Re: [R] Graphics function question

2008-07-29 Thread Marc Schwartz
Just a quick follow up to my own post here. In offline exchanges, Peter 
noted a typo in the call to plot(), where I was missing the final ). 
Also, I had the order of the unnamed arguments to segments() off.


Here is the corrected code, with some additional teaks:


plotdiffs - function(x, y, startage, endage, decline)
{
  startvar - x - y
  endvar - x

  s1 - seq(47, 89, length = 10)

  ymin - min(startvar, endvar, na.rm = TRUE)
  ymax - max(startvar, endvar, na.rm = TRUE)

  s2 - seq(ymin, ymax, length = 10)

  label - gsub(^.+\\$, , deparse(substitute(x)))

  plot(s1, s2, type = n,  xlab = Age,  ylab = label,
   main = paste(Age, decline and, label))

  segments(startage, startvar, endage, endvar,  col = decline)

  legend(topleft, legend = c(Stable, Decline), lty = 1,
 col = c(1, 2))
}



Marc

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


Re: [R] combining zoo series with an overlapping index?

2008-07-29 Thread stephen sefick
My apologies it worked perfectly.

On Tue, Jul 29, 2008 at 12:40 PM, stephen sefick [EMAIL PROTECTED] wrote:

 plot(replace(day, is.list(index(lin)), coredata(lin)))




-- 
Let's not spend our time and resources thinking about things that are so
little or so large that all they really do for us is puff us up and make us
feel like gods. We are mammals, and have not exhausted the annoying little
problems of being mammals.

-K. Mullis

[[alternative HTML version deleted]]

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


Re: [R] ls() and memory question

2008-07-29 Thread jim holtman
It is memory.  You should probably not have 40% of your RAM allocated
to R so that if you do create some large objects, you will have room
for them.  How much physical memory do you have, can you put some
numbers on how many objects you need, if you have a lot of them,
have you considered using lists.

On Tue, Jul 29, 2008 at 12:36 PM, Edna Bell [EMAIL PROTECTED] wrote:
 Hi again!

 I put in ls() to check the objects in my workspace.

 Is there a limit on how many objects I can have, please?  Or does it
 depend on the memory, please?

 TIA,
 Edna Bell

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




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

What is the problem you are trying to solve?

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


[R] R command history -- can it be like Matlab's?

2008-07-29 Thread losemind

Hi all,

In R GUI window, if you use up and down key, you will be able to recall
the previous and next command that has been used and stored in the command
history cache.

But there is one inconvenience:

In Matlab, you can type the first a few characters of the command that you
previously used, 

and then press “up key, you will be able to get to that previous command
very fast,

no matter how long ago it has been used...

In R GUI, there is no such functionality.

Any thoughts on making R work my productivity?

--

Also where do I find good debugger for R?

Thanks!
-- 
View this message in context: 
http://www.nabble.com/R-command-historycan-it-be-like-Matlab%27s--tp18719530p18719530.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Coarsening the Resolution of a Dataset

2008-07-29 Thread jim holtman
I assume that you are doing this on one column of the matrix which
should only have 2160 entries in it.  can you send the actual code you
are using.  I tried it with 10,000 samples and it works fine.  So I
need to understand the data structure you are using.  Also the
infinite recursion sounds strange; do you have function like 'cut' or
'c' redefined?  So it would help if you could supply a reproducible
example.

On Tue, Jul 29, 2008 at 10:09 AM, Steve Murray [EMAIL PROTECTED] wrote:

 Unfortunately, when I get to the 'myCuts' line, I receive the following error:

 Error: evaluation nested too deeply: infinite recursion / 
 options(expressions=)?

 ...and I also receive warnings about memory allocation being reached (even 
 though I've already used memory.limit() to maximise the memory) - this is a 
 fairly sizeable dataset afterall, 2160 rows by 4320 columns.

 Therefore I was wondering if there are any alternative ways of coarsening a 
 dataset? Or are there any packages/commands built for this sort of thing?

 Any advice would be much appreciated!

 Thanks again,

 Steve


 _
 Find the best and worst places on the planet
 http://clk.atdmt.com/UKM/go/101719807/direct/01/



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

What is the problem you are trying to solve?

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


Re: [R] R command history -- can it be like Matlab's?

2008-07-29 Thread Prof Brian Ripley
A patch to do this was posted on 2007-09-29 by Glenn Davis.  Some people 
not addicted to Matlab find the behaviour very inconvenient and prefer the 
getline/readline behaviour (triggered by ^R/^S) of Rterm and R on Unixen.


On Tue, 29 Jul 2008, losemind wrote:



Hi all,

In R GUI window, if you use up and down key, you will be able to recall
the previous and next command that has been used and stored in the command
history cache.

But there is one inconvenience:


Or an inconvenience of Matlab not shared by R, depending on your 
preferences.



In Matlab, you can type the first a few characters of the command that you
previously used,

and then press “up key, you will be able to get to that previous command
very fast,

no matter how long ago it has been used...

In R GUI, there is no such functionality.

Any thoughts on making R work my productivity?

--

Also where do I find good debugger for R?

Thanks!
--
View this message in context: 
http://www.nabble.com/R-command-historycan-it-be-like-Matlab%27s--tp18719530p18719530.html
Sent from the R help mailing list archive at Nabble.com.

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



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


Re: [R] try question

2008-07-29 Thread Mike Prager
Edna Bell [EMAIL PROTECTED] wrote:

 Hi still yet again!
 
 I have the following code:
 
  try(log(rnorm(25)),silent=TRUE)
  [1] -0.26396185 NaN NaN -0.13078069 -2.44997193
 -2.15603971 NaN  0.94917495  0.07244544 NaN
 [11] -1.06341127 -0.42293099 -0.53769569  0.95134763  0.93403340
   NaN -0.10502078 NaN  0.30283262 NaN
 [21] -0.11696872 -3.84122332 NaN NaN -0.12808690
 Warning message:
 In log(rnorm(25)) : NaNs produced
 
 
 I thought that putting the silent = TRUE would suppress the
 warnings, please.  What should I do instead, please?

It's not a great idea to take logs of negative numbers. Better
than supressing the resulting messages, you might try something
like this:

 a - rnorm(25)
 b - log(ifelse(a  0, NA, a))

which gives this:

 a
 [1] -0.04816269 -0.50745059  0.15229031  0.54735811
 [5] -0.29896853  1.81854119  0.19462259 -0.48984075
 [9]  0.63489288  1.47432484  1.15295160  0.75842227
[13]  0.07918115  1.04596643  1.31722543 -0.03614219
[17]  0.44072181  0.25358843 -0.49626405 -2.10954780
[21] -0.85815654  1.38983430  0.66592947  0.64700068
[25] -1.17829527

 b
 [1]  NA  NA -1.8819 -0.60265201
 [5]  NA  0.59803463 -1.63669302  NA
 [9] -0.45429898  0.38820015  0.14232526 -0.27651496
[13] -2.53601706  0.04494127  0.27552758  NA
[17] -0.81934141 -1.37204269  NA  NA
[21]  NA  0.32918453 -0.40657152 -0.43540794
[25]  NA
 

See the help page for ifelse() for more information.

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


[R] optimize simultaneously two binomials inequalities using nlm

2008-07-29 Thread emir.toktar

Dear R users,

I´m trying to optimize simultaneously two binomials inequalities used to
acceptance sampling, which are nonlinear solution, so there is no simple
direct solution. Please, let me explain shortly the the problem and the
question as following.

The objective is to obtain the smallest value of 'n' (sample size)
satisfying both inequalities: 

(1-alpha) = pbinom(c, n, p1)  pbinom(c, n, p2) = beta

Where p1 (AQL) and p2 (LTPD) are probabilities parameters (Consumer and
Producer), the alpha and beta are conumer and producer risks, and finally,
the 'n' represents the sample size and the 'c' an acceptance number or
maximum number of defects (nonconforming) in sample size.

Considering that the 'n' and 'c' values are integer variables, it is
commonly not possible to derive an OC curve including the both (p1,1-alpha)
and (p2,beta) points. Some adjacency compromise is commonly required,
achieved by searching a more precise OC curve with respect to one of the
points. 

I´m using Mathematica 6 but it is a Trial, so I would like use R intead (or
better, I need it)!

To exemplify,  In Mathematica I call the function using NMinimize passing
the restriction and parameters:

/* function name findOpt and parameters... */

restriction = (1 - alpha) = CDF[BinomialDistribution[sample_n, p1], c] 
  betha = CDF[BinomialDistribution[sample_n, p2], c] 
   0  alpha  alphamax  0  betha  bethamax   
   1  sample_n = lot_Size  0 = c  lot_size
  p1  p2  p2max ;

fcost = sample_n/lot_Size; 

result = NMinimize[{fcost, restriction}, {sample_n, c, alpha, betha, p2max},
Method - NelderMead, AccuracyGoal - 10];

/* Calling the function findOpt */
findOpt[p1=0.005, lot_size=1000, alphamax=0.05, bethamax =0.05, p2max =
0.04]

/* and I got the return of values of; minimal n, c, alpha, betha and
the p2 or (LTPD) were computed */ {0.514573, {alpha$74 - 0.0218683,
sample_n$74 - 155.231, betha$74 - 0.05,
c$74 - 2, p2$74 - 0.04}}

Now, using R, I would define the pbinom(c, n, prob) given only the minimum
and maximum values to n and c and limits to p1 and p2 probabilities
(Consumer and Producer), similar on the example above.  

I found some examples to optimize equations in R and some tips, but I not be
able to define the sintaxe to use with that functions. Among the functions
that could be used to resolve the problem presented, I found the function
optim() that it is used for unconstrained optimization and the nlm() which
is used for solving nonlinear unconstrained minimization problems. May I
wrong, but the nlm() function would be appropriate to solve this problem, is
it right?

Can I get a pointer to solve this problem using the nlm() function or where
could I get some tips/example to help me, please?


Thank you very much.

EToktar

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


[R] finding a faster way to do an iterative computation

2008-07-29 Thread dxc13

useR's,

I am trying trying to find out if there is a faster way to do a certain
computation.  I have successfully used FOR loops and the apply function to
do this, but it can take some time to fully compute, but I was wondering if
anyone may know of a different function or way to do this:
 x
[1] 1 2 3 4 5
 xk
 [1] 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 6.0

I want to do: abs(x-xk[i]) where i = 1 to length(xk)=13.  It should result
in a 13 by 5 matrix or data frame.  Does anyone have a quicker solution
than FOR loops or apply()?

Much appreciation and thanks,
dxc13
-- 
View this message in context: 
http://www.nabble.com/finding-a-faster-way-to-do-an-iterative-computation-tp18718233p18718233.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] finding a faster way to do an iterative computation

2008-07-29 Thread Christos Hatzis
matrix(rep(x, each=13) - xk, nrow=13) 

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of dxc13
 Sent: Tuesday, July 29, 2008 2:13 PM
 To: r-help@r-project.org
 Subject: [R] finding a faster way to do an iterative computation
 
 
 useR's,
 
 I am trying trying to find out if there is a faster way to do 
 a certain computation.  I have successfully used FOR loops 
 and the apply function to do this, but it can take some time 
 to fully compute, but I was wondering if anyone may know of a 
 different function or way to do this:
  x
 [1] 1 2 3 4 5
  xk
  [1] 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 6.0
 
 I want to do: abs(x-xk[i]) where i = 1 to length(xk)=13.  It 
 should result in a 13 by 5 matrix or data frame.  Does anyone 
 have a quicker solution than FOR loops or apply()?
 
 Much appreciation and thanks,
 dxc13
 --
 View this message in context: 
 http://www.nabble.com/finding-a-faster-way-to-do-an-iterative-
 computation-tp18718233p18718233.html
 Sent from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide 
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 


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


Re: [R] finding a faster way to do an iterative computation

2008-07-29 Thread Rolf Turner


On 30/07/2008, at 6:12 AM, dxc13 wrote:



useR's,

I am trying trying to find out if there is a faster way to do a  
certain
computation.  I have successfully used FOR loops and the apply  
function to
do this, but it can take some time to fully compute, but I was  
wondering if

anyone may know of a different function or way to do this:

x

[1] 1 2 3 4 5

xk

 [1] 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 6.0

I want to do: abs(x-xk[i]) where i = 1 to length(xk)=13.  It should  
result
in a 13 by 5 matrix or data frame.  Does anyone have a quicker  
solution

than FOR loops or apply()?


outer(xk,x,function(a,b){abs(a-b)})

##
Attention:\ This e-mail message is privileged and confid...{{dropped:9}}

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


Re: [R] R command history -- can it be like Matlab's?

2008-07-29 Thread hadley wickham
On Tue, Jul 29, 2008 at 2:45 PM, Prof Brian Ripley
[EMAIL PROTECTED] wrote:
 A patch to do this was posted on 2007-09-29 by Glenn Davis.  Some people not
 addicted to Matlab find the behaviour very inconvenient and prefer the
 getline/readline behaviour (triggered by ^R/^S) of Rterm and R on Unixen.

On unixen you can redefine your up/down arrows in your .inputrc:

\e[A: history-search-backward
\e[B: history-search-forward

which I find really useful, but it definitely takes a few weeks of
getting used to.  I suspect there maybe an equivalent for Rterm.

Hadley



-- 
http://had.co.nz/

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


Re: [R] R command history -- can it be like Matlab's?

2008-07-29 Thread Rolf Turner


On 30/07/2008, at 9:16 AM, hadley wickham wrote:


On Tue, Jul 29, 2008 at 2:45 PM, Prof Brian Ripley
[EMAIL PROTECTED] wrote:
A patch to do this was posted on 2007-09-29 by Glenn Davis.  Some  
people not
addicted to Matlab find the behaviour very inconvenient and prefer  
the
getline/readline behaviour (triggered by ^R/^S) of Rterm and R on  
Unixen.


On unixen you can redefine your up/down arrows in your .inputrc:

\e[A: history-search-backward
\e[B: history-search-forward

which I find really useful, but it definitely takes a few weeks of
getting used to.  I suspect there maybe an equivalent for Rterm.


Those of us who (sensibly! :-) ) use vi have:

set editing-mode vi

in our .inputrc files.

Then esc puts you into vi editing mode whence ``/'' searches
``forwards'' and ``?'' searches ``backwards'' --- with the convention
that the most recent command is the ``top'' of the file.

And then when you've found the command that you want, you can edit it
(with vi syntax) before pressing return to re-issue the command.

No home should be without one.

cheers,

Rolf Turner

##
Attention:\ This e-mail message is privileged and confid...{{dropped:9}}

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


Re: [R] R command history -- can it be like Matlab's?

2008-07-29 Thread hadley wickham
On Tue, Jul 29, 2008 at 4:39 PM, Rolf Turner [EMAIL PROTECTED] wrote:

 On 30/07/2008, at 9:16 AM, hadley wickham wrote:

 On Tue, Jul 29, 2008 at 2:45 PM, Prof Brian Ripley
 [EMAIL PROTECTED] wrote:

 A patch to do this was posted on 2007-09-29 by Glenn Davis.  Some people
 not
 addicted to Matlab find the behaviour very inconvenient and prefer the
 getline/readline behaviour (triggered by ^R/^S) of Rterm and R on Unixen.

 On unixen you can redefine your up/down arrows in your .inputrc:

 \e[A: history-search-backward
 \e[B: history-search-forward

 which I find really useful, but it definitely takes a few weeks of
 getting used to.  I suspect there maybe an equivalent for Rterm.

Those of us who (sensibly! :-) ) use vi have:

set editing-mode vi

in our .inputrc files.

Then esc puts you into vi editing mode whence ``/'' searches
``forwards'' and ``?'' searches ``backwards'' --- with the convention
that the most recent command is the ``top'' of the file.

And then when you've found the command that you want, you can edit it
(with vi syntax) before pressing return to re-issue the command.

And you might also want to beef up the size of your history file (and
automatically remove duplicates). I have the following in my
bash_profile (and I assume similar functionality exists for other
shells):

export HISTSIZE=1000
export HISTFILESIZE=5000
export HISTCONTROL=erasedups

Hadley



-- 
http://had.co.nz/

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


  1   2   >