[R] Complex to polar?

2008-07-27 Thread rkevinburton
I am using fft(x) to compute the fft of a real-valued sequence. This returns 
(as expected an array of complex numbers. Is ther an R function to transform 
these complex numbers into two arrays for the magnitude and phase? Basically 
polar notation?

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.


Re: [R] Complex to polar?

2008-07-27 Thread Prof Brian Ripley

Please read ?complex.

On Sun, 27 Jul 2008, [EMAIL PROTECTED] wrote:

I am using fft(x) to compute the fft of a real-valued sequence. This 
returns (as expected an array of complex numbers. Is ther an R function 
to transform these complex numbers into two arrays for the magnitude and 
phase? Basically polar notation?


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.


PLEASE do.


--
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] S4 : setGeneric for classical methods

2008-07-27 Thread cgenolin

Martin Maechler [EMAIL PROTECTED] a écrit :


CG == Christophe Genolini [EMAIL PROTECTED]
on Sat, 26 Jul 2008 12:12:12 +0200 writes:


CG Martin Maechler [EMAIL PROTECTED] a écrit :
 CG == Christophe Genolini [EMAIL PROTECTED]
 on Tue, 22 Jul 2008 19:04:37 +0200 writes:

CG Prof Brian Ripley [EMAIL PROTECTED] a écrit :
  On Tue, 22 Jul 2008, [EMAIL PROTECTED] wrote:
 
  Hi the list (well, half of the list, only the one who
  are not on holidays...)
 
  I am trying to make an S4 package. When I run it on a
  console, everything seems ok. When I run R CMD check, I
  got an error --- 8
  --
  Error in setMethod(plot, ClassX, function(x, y) { :
  no existing definition for function plot --- 8
  --
 
  If I add
  setGeneric(plot,function(x,y,...){standardGeneric(plot)})
  in my code, then everything is OK.
 
  This is a surprise for me since I taught that we do not
  need to redefine as generic the function that are
  already generic, like plot. Am I wrong ?
 
  Yes. And do read the error message.  It says
 
  no existing definition for function plot
 
  so this is not if the function is S3 or S4 generic, but
  that no such function is visible.
 
  Looks like you forgot to declare a dependence on (or
  import) package graphics.

CG I would have forgot if I knew that I have to declare
CG such a dependence...  Do we have to declare all the
CG depence to all the package ? To base ? To stats ?

 not to base, but to all other packages (since R can be loaded
 without any packages but base).

CG I thaught that the package that are include in R when we
CG start it had not to be include.

 You thought wrongly, and all our documentation mentions that you
 need to use 'Depends:' correctly (and  'Imports(..)' in
 NAMESPACE if you make use of one).

 {and please strongly note the correct spelling of thought !}

 Martin

   CG I correct DESCRIPTION and NAMESPACE adding Depends and import.
   CG But I still get the message :

   CG checking for missing documentation entries ... WARNING
   CG Undocumented S4 methods:
   CG generic 'plot' and siglist 'ClassX'

   CG I can't find what is wrong...

well, it's  a message about *documentation*
and therefore not related to DESCRIPTION
and not much related to NAMESPACE.

It tells you that you don't have correct documentation for the
plot method for ClassX.



What if I don't want to document this method? I did not put ClassX in 
the export, I do not want the user to acces to it.

I want :

 plot,ClassX-method - private
 ClassY-class   - public
 plot,ClassY-method - public

Is it possible ?

Christophe


-- Writing R Extensions,
   section 2.1.3 'Documenting S4 classes and methods'

Regards,
Martin

   CG - DESCRIPTION - Package: packS4
   CG Type: Package
   CG Title: Toy example of S4 package
   CG Version: 0.5
   CG Date: 2008-07-22
   CG Author: Christophe Genolini / INSERM U669
   CG Maintainer: [EMAIL PROTECTED]
   CG Description: Package built to illustration package 
construction with S4

   CG License: GPL (=2)
   CG LazyLoad: yes
   CG Depends: methods, graphics
   CG Collate: global.R ClassX.R ClassY.R ClassX-ClassY.R ClassZ.R


   CG -- NAMESPACE --
   CG export(
   CG classZ,
   CG functionClassicA
   CG )
   CG exportMethods(
   CG getZ1,
   CG setZ1-,
   CG publicA,
   CG plot
   CG )
   CG exportClasses(
   CG ClassY,
   CG ClassZ
   CG )

   CG import(graphics)

   CG Christophe




   CG 
   CG Ce message a ete envoye par IMP, grace a l'Universite Paris 10 
Nanterre






__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Retain plot?

2008-07-27 Thread rkevinburton
Every time I issue a plot command it removes the current plot and replaces it 
with the plot that I just issued. I would like to keep both plots. I looked in 
the documentation and found plot.new but was unable to get it to work. I think 
I want that functionality. Any suggestions?

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] A easy way to write formula

2008-07-27 Thread Jinsong Zhao
Hi

I have a data frame, including x1, x2, x3, and y. I use lm() to fit
second-order linear model, like the following:

ft - lm(y ~ x1 + x2 + x3 + I(x1 * x1) + I(x1 * x2) + I(x1 * x3) + I(x2
* x2) + I(x2 * x3) + I(x3 * x3), mydata)

if the independent variable number is large, the formula will be very
long. Is there a easy way to write formula like the above one? I have
read the R introduction, however, I don't find a easy way.

Any hints will be appreciated. Thanks,

Jinsong

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Color of box frame in Legend (Was: Matrix barplot)

2008-07-27 Thread Andreas Tille

On Sun, 27 Jul 2008, S Ellison wrote:


Looking at the legend() source the filled box line colour is hardcoded :
   if (mfill) {
   if (plot) {
   fill - rep(fill, length.out = n.leg)
   rect2(left = xt, top = yt + ybox/2, dx = xbox, dy = ybox,
   col = fill, density = density, angle = angle,
   border = black)
   }
   xt - xt + dx.fill
   }

... so it looks like you can have any colour as long as it's black.


Ups, but this is a bug, isn't it?  Why is this black hardcoded?


However, you could copy the legend code (type legendCR to see it, then
paste it into a text editor, modify the black to, say, par(fg) or
even add your own extra parameter to the function, then just paste the
function back into your own version.


I could try this for the moment - but where can I report this as a bug?

Kind regards

 Andreas.
--
http://fam-tille.de

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Weighted variance function?

2008-07-27 Thread Arun Kumar Saha
ur prog gives following result:

weighted.var(c(1,-1), c(0.5,0.5))
[1] 2

is it ok?

On Thu, Jul 24, 2008 at 7:57 PM, Gavin Simpson [EMAIL PROTECTED]wrote:

 On Thu, 2008-07-24 at 02:25 +0530, Arun Kumar Saha wrote:
  There is a R function to calculate weighted mean : weighted.mean() under
  stats package. Is there any direct R function for calculating weighted
  variance as well?

 Here are two ways; weighted.var() is via the usual formula and
 weighted.var2() uses a running sums approach. The formulae for which are
 both on the weighted mean entry page on wikipedia for example.

 The removal of NA is as per weighted.mean, but I have not included any
 of the sanity checks that that functions contains.

 weighted.var - function(x, w, na.rm = FALSE) {
if (na.rm) {
w - w[i - !is.na(x)]
x - x[i]
}
sum.w - sum(w)
sum.w2 - sum(w^2)
mean.w - sum(x * w) / sum(w)
(sum.w / (sum.w^2 - sum.w2)) * sum(w * (x - mean.w)^2, na.rm =
 na.rm)
 }

 weighted.var2 - function(x, w, na.rm = FALSE) {
if (na.rm) {
w - w[i - !is.na(x)]
x - x[i]
}
sum.w - sum(w)
(sum(w*x^2) * sum.w - sum(w*x)^2) / (sum.w^2 - sum(w^2))
 }

 ## from example section in ?weighted.mean
 ## GPA from Siegel 1994
 wt - c(5,  5,  4,  1)/15
 x - c(3.7,3.3,3.5,2.8)
 weighted.mean(x,wt)

 weighted.var(x, wt)

 weighted.var2(x, wt)

 And some timings:

  system.time(replicate(10, weighted.var(x, wt)))
   user  system elapsed
  2.679   0.014   2.820
  system.time(replicate(10, weighted.var2(x, wt)))
   user  system elapsed
  2.224   0.010   2.315

 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/http://www.ucl.ac.uk/%7Eucfagls/
  UK. WC1E 6BT. [w] http://www.freshwaters.org.uk
 %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%



[[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] Transfer Function Modeling

2008-07-27 Thread Ashish Kumar
Hi,

 

I would like to know how to build the transfer function modeling for
bivariate time series data. I tried searching for relevant threads, but
could not find much help with it.

 

Thanks

Ashish


[[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] Lattice wireframe: How to avoid drawing lines around polygons when using shade=TRUE

2008-07-27 Thread Oliver M. Haynold
On Sun, 27 Jul 2008 05:00:48 +, Oliver M. Haynold wrote:
 I am using wireframe from the lattice package, with the shade option set
 to TRUE. When I output to PDF or Postscript, a line gets drawn around
 each polygon of my surface which causes ugly Moiré effects and doesn't
 make sense in my application (think the plot of Maunga Whau--gridlines
 don't make sense). This does not happen when I display on the screen--
 then I just get the surface as intended.

I have done some research on the problem. On inspection of the PDF file 
generated by lattice, it turns out that it does not actually contain any 
instructions to draw lines between the triangles of the surface plot. 
Yet, Evince and Acrobat Reader both show these lines, with disagreeable 
consequences for printing and screen display. Ghostscript, however, does 
not and renders the surface as a solid block, as it should be. Thus, a 
workaround for the problem is to put the file generated by lattice 
through Ghostscript and to convert it into a high-resolution image file.

The most likely reason for the problem I can think of is that the 
triangles generated by lattice don't overlap but just border each other. 
This might cause some renderers to leave a little space between them. If 
this explanation is correct, a solution would be to have lattice, or 
specifically wireframePanelCalculations, plot the triangles with a tiny 
bit of overlap.

Best,

Oliver

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


Re: [R] A easy way to write formula

2008-07-27 Thread ONKELINX, Thierry
(x1 + x2 + x3) ^2 will give you the main effects and the interactions.
So this will shorten your equation to 

ft - lm(y ~ (x1 + x2 + x3) ^2 + I(x1 * x1)+ I(x2 * x2) + I(x3 * x3),
mydata)

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 Jinsong Zhao
Verzonden: zondag 27 juli 2008 9:39
Aan: r-help@r-project.org
Onderwerp: [R] A easy way to write formula

Hi

I have a data frame, including x1, x2, x3, and y. I use lm() to fit
second-order linear model, like the following:

ft - lm(y ~ x1 + x2 + x3 + I(x1 * x1) + I(x1 * x2) + I(x1 * x3) + I(x2
* x2) + I(x2 * x3) + I(x3 * x3), mydata)

if the independent variable number is large, the formula will be very
long. Is there a easy way to write formula like the above one? I have
read the R introduction, however, I don't find a easy way.

Any hints will be appreciated. Thanks,

Jinsong

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Retain plot?

2008-07-27 Thread ONKELINX, Thierry
Kevin,

Open a new plot window before plotting the second plot.

plot(rnorm(10), runif(10))
X11()
plot(rnorm(10), runif(10))

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 [EMAIL PROTECTED]
Verzonden: zondag 27 juli 2008 9:34
Aan: r-help@r-project.org
Onderwerp: [R] Retain plot?

Every time I issue a plot command it removes the current plot and
replaces it with the plot that I just issued. I would like to keep both
plots. I looked in the documentation and found plot.new but was unable
to get it to work. I think I want that functionality. Any suggestions?

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] A easy way to write formula

2008-07-27 Thread Mark Difford

Hi Jinsong and Thierry,

 (x1 + x2 + x3) ^2 will give you the main effects and the interactions.

Although it wasn't specifically requested it is perhaps important to note
that (...)^2 doesn't expand to give _all_ interaction terms, only
interactions to the second order, so the interaction term x1:x2:x3 will be
missing, To get these, do

(x1 + x2 + x3)^3

Regards, Mark.


ONKELINX, Thierry wrote:
 
 (x1 + x2 + x3) ^2 will give you the main effects and the interactions.
 So this will shorten your equation to 
 
 ft - lm(y ~ (x1 + x2 + x3) ^2 + I(x1 * x1)+ I(x2 * x2) + I(x3 * x3),
 mydata)
 
 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 Jinsong Zhao
 Verzonden: zondag 27 juli 2008 9:39
 Aan: r-help@r-project.org
 Onderwerp: [R] A easy way to write formula
 
 Hi
 
 I have a data frame, including x1, x2, x3, and y. I use lm() to fit
 second-order linear model, like the following:
 
 ft - lm(y ~ x1 + x2 + x3 + I(x1 * x1) + I(x1 * x2) + I(x1 * x3) + I(x2
 * x2) + I(x2 * x3) + I(x3 * x3), mydata)
 
 if the independent variable number is large, the formula will be very
 long. Is there a easy way to write formula like the above one? I have
 read the R introduction, however, I don't find a easy way.
 
 Any hints will be appreciated. Thanks,
 
 Jinsong
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/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/A-easy-way-to-write-formula-tp18674075p18674632.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] Weighted variance function?

2008-07-27 Thread Arun Kumar Saha
I feel there is something wrong in definition in weighted variance.

Say for a fair dice the variance of outcome should be : 2.92 (
http://en.wikipedia.org/wiki/Variance)

However if I use cov.wt() or weighted.var() by Gavin, I get :

cov.wt(as.data.frame(1:6), rep(0.6, 6))
$cov
1:6
1:6 3.5

$center
1:6
3.5

$n.obs
[1] 6

$wt
[1] 0.167 0.167 0.167 0.167 0.167 0.167

i.e. 3.5

Therefore if I want to calculate Variance for a r.v. with different prob for
different values then should not use those formulae. Is it the case?

On Sun, Jul 27, 2008 at 2:30 PM, Patrick Burns [EMAIL PROTECTED]wrote:

 Have you seen 'cov.wt'?


 Patrick Burns
 [EMAIL PROTECTED]
 +44 (0)20 8525 0696
 http://www.burns-stat.com
 (home of S Poetry and A Guide for the Unwilling S User)

 Arun Kumar Saha wrote:

 ur prog gives following result:

 weighted.var(c(1,-1), c(0.5,0.5))
 [1] 2

 is it ok?

 On Thu, Jul 24, 2008 at 7:57 PM, Gavin Simpson [EMAIL PROTECTED]
 wrote:



 On Thu, 2008-07-24 at 02:25 +0530, Arun Kumar Saha wrote:


 There is a R function to calculate weighted mean : weighted.mean() under
 stats package. Is there any direct R function for calculating weighted
 variance as well?


 Here are two ways; weighted.var() is via the usual formula and
 weighted.var2() uses a running sums approach. The formulae for which are
 both on the weighted mean entry page on wikipedia for example.

 The removal of NA is as per weighted.mean, but I have not included any
 of the sanity checks that that functions contains.

 weighted.var - function(x, w, na.rm = FALSE) {
   if (na.rm) {
   w - w[i - !is.na(x)]
   x - x[i]
   }
   sum.w - sum(w)
   sum.w2 - sum(w^2)
   mean.w - sum(x * w) / sum(w)
   (sum.w / (sum.w^2 - sum.w2)) * sum(w * (x - mean.w)^2, na.rm =
 na.rm)
 }

 weighted.var2 - function(x, w, na.rm = FALSE) {
   if (na.rm) {
   w - w[i - !is.na(x)]
   x - x[i]
   }
   sum.w - sum(w)
   (sum(w*x^2) * sum.w - sum(w*x)^2) / (sum.w^2 - sum(w^2))
 }

 ## from example section in ?weighted.mean
 ## GPA from Siegel 1994
 wt - c(5,  5,  4,  1)/15
 x - c(3.7,3.3,3.5,2.8)
 weighted.mean(x,wt)

 weighted.var(x, wt)

 weighted.var2(x, wt)

 And some timings:



 system.time(replicate(10, weighted.var(x, wt)))


  user  system elapsed
  2.679   0.014   2.820


 system.time(replicate(10, weighted.var2(x, wt)))


  user  system elapsed
  2.224   0.010   2.315

 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/http://www.ucl.ac.uk/%7Eucfagls/
 http://www.ucl.ac.uk/%7Eucfagls/
  UK. WC1E 6BT. [w] http://www.freshwaters.org.uk
 %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%





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







--

[[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] Weighted variance function?

2008-07-27 Thread Arun Kumar Saha
oh now I got, R is taking unbiased estimator, anyways thanks

On Sun, Jul 27, 2008 at 2:54 PM, Arun Kumar Saha
[EMAIL PROTECTED]wrote:

 I feel there is something wrong in definition in weighted variance.

 Say for a fair dice the variance of outcome should be : 2.92 (
 http://en.wikipedia.org/wiki/Variance)

 However if I use cov.wt() or weighted.var() by Gavin, I get :

 cov.wt(as.data.frame(1:6), rep(0.6, 6))
 $cov
 1:6
 1:6 3.5

 $center
 1:6
 3.5

 $n.obs
 [1] 6

 $wt
 [1] 0.167 0.167 0.167 0.167 0.167 0.167

 i.e. 3.5

 Therefore if I want to calculate Variance for a r.v. with different prob
 for different values then should not use those formulae. Is it the case?


 On Sun, Jul 27, 2008 at 2:30 PM, Patrick Burns [EMAIL PROTECTED]wrote:

 Have you seen 'cov.wt'?


 Patrick Burns
 [EMAIL PROTECTED]
 +44 (0)20 8525 0696
 http://www.burns-stat.com
 (home of S Poetry and A Guide for the Unwilling S User)

 Arun Kumar Saha wrote:

 ur prog gives following result:

 weighted.var(c(1,-1), c(0.5,0.5))
 [1] 2

 is it ok?

 On Thu, Jul 24, 2008 at 7:57 PM, Gavin Simpson [EMAIL PROTECTED]
 wrote:



 On Thu, 2008-07-24 at 02:25 +0530, Arun Kumar Saha wrote:


 There is a R function to calculate weighted mean : weighted.mean()
 under
 stats package. Is there any direct R function for calculating weighted
 variance as well?


 Here are two ways; weighted.var() is via the usual formula and
 weighted.var2() uses a running sums approach. The formulae for which are
 both on the weighted mean entry page on wikipedia for example.

 The removal of NA is as per weighted.mean, but I have not included any
 of the sanity checks that that functions contains.

 weighted.var - function(x, w, na.rm = FALSE) {
   if (na.rm) {
   w - w[i - !is.na(x)]
   x - x[i]
   }
   sum.w - sum(w)
   sum.w2 - sum(w^2)
   mean.w - sum(x * w) / sum(w)
   (sum.w / (sum.w^2 - sum.w2)) * sum(w * (x - mean.w)^2, na.rm =
 na.rm)
 }

 weighted.var2 - function(x, w, na.rm = FALSE) {
   if (na.rm) {
   w - w[i - !is.na(x)]
   x - x[i]
   }
   sum.w - sum(w)
   (sum(w*x^2) * sum.w - sum(w*x)^2) / (sum.w^2 - sum(w^2))
 }

 ## from example section in ?weighted.mean
 ## GPA from Siegel 1994
 wt - c(5,  5,  4,  1)/15
 x - c(3.7,3.3,3.5,2.8)
 weighted.mean(x,wt)

 weighted.var(x, wt)

 weighted.var2(x, wt)

 And some timings:



 system.time(replicate(10, weighted.var(x, wt)))


  user  system elapsed
  2.679   0.014   2.820


 system.time(replicate(10, weighted.var2(x, wt)))


  user  system elapsed
  2.224   0.010   2.315

 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/http://www.ucl.ac.uk/%7Eucfagls/
 http://www.ucl.ac.uk/%7Eucfagls/
  UK. WC1E 6BT. [w] http://www.freshwaters.org.uk
 %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%





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







 --




--

[[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] GLS no terms component error

2008-07-27 Thread mm745
Thank you, Bill. I think I misread some text in my statistics handbook. This has
helped to clear it up. Sorry for taking up your time and thank you for replying.

M. Grace

Quoting [EMAIL PROTECTED]:

 The help file for dwtest says

 formula: a symbolic description for the model to be tested (or a fitted
 lm object).

 It is not clear from your message, but it would seem that you used fitp2
 as the argument to dwtest.  (In any case, this gives the error message
 you quote?)

 But fitp2 is neither a formula, nor is it an object ineriting from lm.

 Why did you expect dwtest to work on a gls object?


 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: Sunday, 27 July 2008 6:26 AM
 To: r-help@r-project.org
 Subject: [R] GLS no terms component error

 Hi,

 I want to test for independence in my GLS model fitp2, but when I try to
 use the
 dwtest function in the lmtest library, I get the error message Error in
 terms.default(formula) : no terms component.

 The model and data set are below. Any suggestions would be really
 helpful!
 Thanks a lot in advance,

 M. Grace


 fitp2:
 fitp2-gls(V3_total~D_total+P_total,data=pdata,weights=varPower())

 based on the data set pdata:
 pdata
   D_total V3_total   P_total
 1   75   158000   302.033
 2  221   258000   126.664
 3 1050   28   538.547
 4 8100   208000  1090.068
 522100   235000  2038.937
 6 5300   57  2626.226
 7 3250   454000  2326.473
 8 6540   169000  5613.525
 9 5248   247000  6011.532
 10   22484   542000 15140.006
 11   30841   53 13324.048
 127480   497000 10481.812
 132664   467000  5776.274
 14 432   285000  1849.763


 --
 University of St Andrews Webmail: https://webmail.st-andrews.ac.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.







--
University of St Andrews Webmail: https://webmail.st-andrews.ac.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] re sponse surface analysis

2008-07-27 Thread Tom La Bone

There appears to be a very promising response surface package being discussed
at useR-2008, but I have been unable to find the package on CRAN or contact
the authors.

www.statistik.uni-dortmund.de/useR-2008/abstracts/Sztendur+Diamond.pdf


Tom




Jinsong Zhao wrote:
 
 Hi,
 
 Is there a package that could do response surface analysis equivalent to
 SAS RSREG procedure? Thanks!
 
 Regards,
 
 Jinsong
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/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/response-surface-analysis-tp18666106p18675308.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] Object-oriented programming in R for Java programmers?

2008-07-27 Thread Werner Wernersen
Hi,

I was wondering if anybody might have a reference for
me: My R code is growing and getting more and more
confusing. Thus, I figure it's time to switch to
object-oriented again. I have done oo programming in
C++ and Java before but the first few tutorial on R oo
were a bit confusing for me. 

Is there any brief tutorial on oo programming in R
especially for people who have done oo in Java or C++
before? That would be really helpful.

Many thanks and have a great Sunday,
  Werner


  __

Dem pfiffigeren Posteingang.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] response surface analysis

2008-07-27 Thread Jinsong Zhao

Tom La Bone wrote:

There appears to be a very promising response surface package being discussed
at useR-2008, but I have been unable to find the package on CRAN or contact
the authors.

www.statistik.uni-dortmund.de/useR-2008/abstracts/Sztendur+Diamond.pdf


Tom


Thanks for pointing me to this. From the poster, it's what I try to 
find. Hope the authors could upload rsm package to CRAN as early as 
possible.


Regards,
Jinsong

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] equivalent R functions for Numerical Recipes fitxy and fitexy ?

2008-07-27 Thread Marc Fischer

Dear Folks,

We need to fit the model y~x  assuming there are random errors in both  
x and y.  Numerical Recipes (Press et al.) have two useful functions,  
fitxy and fitexy, which handle cases of unspecified and specified  
errors respectively.


Are there equivalent functions in base R or a installable package?
Alternatively, has anyone written a wrapper to provide an interface to  
a library of Numerical Recipes routines constructed for Linux?


Kind regards,

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.


[R] Link functions in SEM

2008-07-27 Thread Jeroen Ooms

Is it possible to fit a structural equation model with link functions in R? I
am trying to build a logistic-regression-like model in sem, because
incorporating the dichotomous variables linearly seems inappropriate. Mplus
can do something similar by specifying a 'link' parameter, but I would like
to be able to do it in R, ofcourse. 

I have explored the 'sem' package from John Fox, but it does not seem to be
able to fit non-linear relations. Is there some R-package or way to get this
done? I have also considered creating a seperate latent variable in the sem
model for the systematic component of the predictors, but then I still need
a way to fix a non-linear link from the systematic component to the
dichotomous Y variable. 
-- 
View this message in context: 
http://www.nabble.com/Link-functions-in-SEM-tp18679236p18679236.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] Colors in Sweave

2008-07-27 Thread Stephen Tucker
Hi list,

I was using Sweave and was wondering if anyone has had any luck changing the 
font colors of the code chunks. For instance, in my .Rnw preample I tried 
including:

===
\usepackage[usenames]{colors}
\definecolor{darkred}{rgb}{0.545,0,0}
\definecolor{midnightblue}{rgb}{0.098,0.098,0.439}
\DefineVerbatimEnvironment{Sinput}{Verbatim}{fontshape=sl,formatcom={\color{midnightblue}}}
\DefineVerbatimEnvironment{Soutput}{Verbatim}{formatcom={\color{darkred}}}
\DefineVerbatimEnvironment{Scode}{Verbatim}{fontshape=sl,formatcom={\color{blue}}}
===

which works in the sense that colors do show up in the processed pdf document, 
but extra spaces in between the input and output code chunks appear (which is 
not the case when colors are not specified), resulting in a long document with 
many blank lines.

Alternatively, in the resulting tex document I replaced all instances of 
\begin{Sinput}... \end{Sinput} with {\color{midnightblue}\begin{Sinput} ... 
\end{Sinput}} and darkred for Soutput and so on. When I do this, I get the 
following error message:

=== LaTeX excerpt ===
\begin{Schunk}
{\color{midnightblue}\begin{Sinput}
 ISOdatetime(1970, 1, 1, 0, 0, 0, ) - ISOdatetime(1970, 1, 1, 
+ 0, 0, 0, GMT)
\end{Sinput}}
{\color{darkred}\begin{Soutput}
Time difference of 8 hours
\end{Soutput}}
\end{Schunk}
===

=== error message ===
)
! FancyVerb Error:
  Extraneous input `}\end{}' between \end{Sinput} and line end
.
[EMAIL PROTECTED] ... {FancyVerb Error:
\space \space #1
}
  
l.13 \end{Sinput}}
===

I guess I don't know enough of the Schunk/Sinput/Soutput definitions to toy 
with it and was wondering if anyone had tried something similar. 

Thanks!

Stephen

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Object-oriented programming in R for Java programmers?

2008-07-27 Thread Johannes Huesing
Werner Wernersen [EMAIL PROTECTED] [Sun, Jul 27, 2008 at 01:55:26PM CEST]:
[...]
 Is there any brief tutorial on oo programming in R
 especially for people who have done oo in Java or C++
 before? That would be really helpful.

How S4 methods work highlights the differences between R's function-centric
approach to OO and the class-centric approach highlighted by languages such
as C++. This 10-page article is written by John Chambers.

-- 
Johannes Hüsing   There is something fascinating about science. 
  One gets such wholesale returns of conjecture 
mailto:[EMAIL PROTECTED]  from such a trifling investment of fact.  
  
http://derwisch.wikidot.com (Mark Twain, Life on the Mississippi)

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] melting a list: basic question

2008-07-27 Thread hadley wickham
 Is it possible to install your development version of reshape? I could not
 find it alongside of ggplot2 on github.

Not yet, but I'm working on it.

 If not, I've added ... in the
 method for the current version and it seems to work for me.

Yes, that's the exact change I made.

Thanks!

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] A easy way to write formula

2008-07-27 Thread hadley wickham
On Sun, Jul 27, 2008 at 4:19 AM, Mark Difford [EMAIL PROTECTED] wrote:

 Hi Jinsong and Thierry,

 (x1 + x2 + x3) ^2 will give you the main effects and the interactions.

 Although it wasn't specifically requested it is perhaps important to note
 that (...)^2 doesn't expand to give _all_ interaction terms, only
 interactions to the second order, so the interaction term x1:x2:x3 will be
 missing, To get these, do

 (x1 + x2 + x3)^3

Or just

x1 * x2 * x3

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] 64-bit R on Mac OS X 10.5.4

2008-07-27 Thread joseph


Hi Matt
Your method is the easiest way for me to install the 64-bit R.  I followed the 
directions on your web site and then did the following:
R --arch=x86_64
source(http://bioconductor.org/biocLite.R;) 
biocLite(type = source,lib = 
/Library/Frameworks/R.framework/Versions/2.8/Resources/RLib64)

I got many errors and warnings which I copied to the attached file.
Any suggestions on how to fix the errors will be appreciated.
Thank you for your help
Joseph


 sessionInfo()
R version 2.8.0 Under development (unstable) (2008-07-24 r46120) 
i386-apple-darwin9.4.0 
locale:
en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8
attached base packages:
[1] stats graphics  grDevices utils datasets  methods   base 
loaded via a namespace (and not attached):
[1] tools_2.8.0

This is the computer I am using:
Hardware Overview:
  Model Name:Mac Pro
  Model Identifier:MacPro1,1
  Processor Name:Dual-Core Intel Xeon
  Processor Speed:2.66 GHz
  Number Of Processors:2
  Total Number Of Cores:4
  L2 Cache (per processor):4 MB
  Memory:20 GB
  Bus Speed:1.33 GHz
  Boot ROM Version:MP11.005C.B08
  SMC Version:1.7f10
  Serial Number:G87052SGUPZ

System Software Overview:
  System Version:Mac OS X 10.5.4 (9E17)
  Kernel Version:Darwin 9.4.0
  Boot Volume:Macintosh HD
  Boot Mode:Normal


- Original Message 
From: Matthew Keller [EMAIL PROTECTED]
To: Steven McKinney [EMAIL PROTECTED]
Cc: joseph [EMAIL PROTECTED]; r-help@r-project.org
Sent: Saturday, July 26, 2008 9:15:41 AM
Subject: Re: [R] 64-bit R on Mac OS X 10.4.5

Hi Joseph,

For what it is worth (which might not be that much!), I have written
down step by step instructions on my website for getting 64 bit R
working under Leopard - it should be much different with Tiger:
http://www.matthewckeller.com/html/64_bit_r_on_mac.html. I think it'll
work but there may be some mistakes in there. I'd be happy to take
suggestions...

Matt

On Fri, Jul 25, 2008 at 7:23 PM, Steven McKinney [EMAIL PROTECTED] wrote:

 Hello

 I haven't found better instructions, it's just
 not an easy thing to do.

 You might also consider joining the r-sig-mac group
 and reviewing threads there for additional information.

 Rather than try to configure, make and install with
 one giant command, I'd suggest breaking the task down
 until you have worked out all the details.  Then you
 can build more routinely with giant commands such as
 the one you report below.

 Note that the page of the URL you posted says this:
 R on Mac OS X 10.5 (Leopard)

 so look around the r.research.att.com
 website for older 10.4 related instructions.
 OS X 10.5 is different enough from 10.4 that
 copying and pasting the 10.5 instructions
 will not work everywhere.

 You don't say what kind of computer you have
 (Intel or PowerPC), if you have a power pc
 then -arch x86_64 is not right.  So do
 report your hardware configuration as well
 as your operating system configuration.

 Try breaking up the build process into steps
 and review the messages generated at each step.

 e.g. the configure step:  Just do these bits

 cd rd64
 ../R-devel/configure SHELL='/bin/bash' \
 r_arch=x86_64 \
 CC=gcc -arch x86_64 -std=gnu99 \
 CXX=g++ -arch x86_64 \
 OBJC=gcc -arch x86_64 \
 F77=gfortran -arch x86_64 \
 FC=gfortran -arch x86_64 \
 --with-system-zlib \
 --with-blas='-framework vecLib' \
 --with-lapack 1 configure.R.txt 21

 The end bits of this version of the
 configure command redirect output and
 error messages to a file that you can then
 read, to see which bits you are missing and
 which bits cause problems.  It's tough to catch those
 as they scroll by in a terminal window.

 If configure runs without any signs of trouble,
 try the make.

 make 1 make.R.txt 21

 Then you can review all the make output for signs
 of problems and error messages.


 Here's a configure command that worked for me to
 build a 64-bit R on 10.4 a while back:


 ./configure --host=powerpc64-apple-darwin8.10.0 
 --build=powerpc64-apple-darwin8.10.0 \
 --prefix=/usr/local/lib64 'CC=gcc-4.0 -arch ppc64' 'CXX=g++ -arch ppc64' \
 'FC=gfortran-4.0 -arch ppc64' 'F77=gfortran-4.0 -arch ppc64' \
 'CFLAGS=-g -O3 -mtune=G5 -mcpu=G5' 'FFLAGS=-g -O3 -mtune=G5 -mcpu=G5' \
 'LDFLAGS=-arch ppc64 -m64 -L/usr/local/lib' 'CXXFLAGS=-g -O3 -mtune=G5 
 -mcpu=G5' \
 'FCFLAGS=-g -O3 -mtune=G5 -mcpu=G5' --disable-R-framework --enable-R-shlib \
 '--with-blas=-framework vecLib' --with-lapack --without-iconv 1 
 configure.R.txt 21


 and also note that you need to have root privileges when doing the
 make install, so you either need to run the command as root
 (not so common on Mac OS X) or run

 sudo make install

 and enter your password etc.

 HTH

 Steve McKinney




 -Original Message-
 From: [EMAIL PROTECTED] on behalf of joseph
 Sent: Fri 7/25/2008 5:07 PM
 To: r-help@r-project.org
 Cc: r-help@r-project.org
 Subject: [R] 64-bit R on Mac OS X 10.4.5

 Hello
 I have a Mac OS X 10.4.5.  I 

Re: [R] Object-oriented programming in R for Java programmers?

2008-07-27 Thread Jeroen Ooms


Werner Wernersen wrote:
 
 I have done oo programming in C++ and Java before but the first few
 tutorial on R oo were a bit confusing for me. 
 
My personal experience is that the type of OO programming that makes for
example Java code nice and easy to structure is not possible in R. The
biggest problem for me is, correct me if i'm wrong, that R does not store
memory references in its identifiers (variables), but a creates new 'object'
for every identifier. For example, if you do in R:
a - 123
b - 456
a - b
a - 789
b

If you would do a similar thing in Java, bot the variables 'a' and 'b' would
have changed to the value 789, because both the variables now refer to the
same object. So manipulating one also manipulates the other. However, in R,
variable 'b' still stores the value of 456 at the end of this example.

What this means for OO programming is that there is no easy way to reference
to objects. For example, if you have want to access a variable which is part
of a dataframe object, which is nested in a list object, you can only
reference this object every time using list$dataframe$varname, where in
Java, you could create a new variable that references to this specific
object. Once datastructures become more complex, this handicap gets bigger
and bigger.

However, please don't take this for granted since I am only a beginning R
programmer (and slightly more advanced Java programmer). I am also really
curious what others have to say because I am having the same problem of my
R-code quickly getting messy and confusing.
-- 
View this message in context: 
http://www.nabble.com/Object-oriented-programming-in-R-for-Java-programmers--tp18675688p18679872.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] Object-oriented programming in R for Java programmers?

2008-07-27 Thread Stephen Tucker
I've never used it myself but I recall reading about this R.oo package a while 
ago - might be worth considering?

http://www1.maths.lth.se/help/R/R.oo/



- Original Message 
From: Werner Wernersen [EMAIL PROTECTED]
To: [EMAIL PROTECTED]
Sent: Sunday, July 27, 2008 4:55:26 AM
Subject: [R] Object-oriented programming in R for Java programmers?

Hi,

I was wondering if anybody might have a reference for
me: My R code is growing and getting more and more
confusing. Thus, I figure it's time to switch to
object-oriented again. I have done oo programming in
C++ and Java before but the first few tutorial on R oo
were a bit confusing for me. 

Is there any brief tutorial on oo programming in R
especially for people who have done oo in Java or C++
before? That would be really helpful.

Many thanks and have a great Sunday,
  Werner


  __

Dem pfiffigeren Posteingang.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Fitting a Bivariate Poisson log-normal distribution

2008-07-27 Thread Famoye, Felix
I wrote an R program to fit a bivariate Poisson log-normal model. The
model is first proposed by Aitchison and Ho (1989) in Biometrika, Volume
76 pages 643-653. The model involves using a double integration. The
following is the program and the data that I used. My major problem was
that R was running for a long time (more than 3 hours) with no results.
I would like to place some printing at some places to see if R was doing
something. Unfortunately, there was nothing printed when I included
printing in the function f. Note that I have used the package adapt
for double integration. Can someone provide some help on this problem? I
would like to know if R is doing anything at all in fitting the model.
Thank you for your time.

#This program will be used to estimate the parameters of
Poisson-lognormal model. 
yd = read.table(file=H:\\Bivariate\\bact.txt)
 
y1 = yd[[1]]
y2 = yd[[2]]
n = length(y1)
x0= rep(1, n)
y1  = cbind(y1)
y2  = cbind(y2)
yy=cbind(y1,y2)
xx = cbind(x0,yd[[3]])
av1=mean(y1)
av2=mean(y2)
cov=var(yy)
mu=cbind(av1, av2)
rho=cov[1,2]/sqrt(cov[1,1]*cov[2,2])
sd1=sqrt(cov[1,1])
sd2=sqrt(cov[2,2]) 
oth=rbind(sd1,sd2,rho)
z = vector(length=2)  #This is a column vector with two rows

#To find the inital estimates for the regression coefficients
para1 = ginv(t(xx) %*% xx) %*% (t(xx) %*% log(y1+0.5))
para2 = ginv(t(xx) %*% xx) %*% (t(xx) %*% log(y2+0.5))
parab = rbind(para1, para2, oth)
parab

#-
# Fitting Poisson lognorma model
#f is the function that defines the (negative) log likelihood 
f - function(parab,y1,y2,xx)
{ 
vv1=(parab[5])^2
vv2=(parab[6])^2 
para.1= cbind(parab[1:2])
para.2= cbind(parab[3:4])
mu.1 = xx%*%para.1 - vv1/2
#mu.1 = exp(mu.1)
mu.2 = xx%*%para.2 - vv2/2
#mu.2 = exp(mu.2)
mu.v=cbind(mu.1,mu.2)
cov[1,1]=vv1
cov[2,2]=vv2
cov[1,2]=parab[7]*parab[5]*parab[6]
cov[2,1]=cov[1,2]
n=length(y1)
va=rep(0.0, n)
for (i in 1:n)
{
fred=function(z) 
{
tz=log(z)-mu.v[i,]
qq=(t(tz)%*%ginv(cov))%*%tz
fq=(exp(-qq/2))/(2*pi*z[1]*z[2]*sqrt(det(cov)))
p1=((z[1]^y1[i])*exp(-z[1]))/gamma(y1[i]+1)
p2=((z[2]^y2[i])*exp(-z[2]))/gamma(y2[i]+1)
int=fq*p1*p2
}
va[i]=adapt(2,lo=c(0,0),up=c(100,100),functn=fred,min=1000,eps=1.e-6,mu.
v=mu.v,cov=cov,y1=y1,y2=y2)$value
}
-sum(va)
}

I.p = nlminb(start=parab,objective=f,y1=y1,y2=y2,xx=xx)
print(Poisson Log-normal Regression Model)
I.p



   0   0   0
   0   0   0
   0   0   0
   0   0   0
   0   0   0
   0   0   0
   0   0   0
   0   0   0
   0   0   0
   0   0   0
   0   0   0
   0   0   0
   1   0   0
   1   0   0
   1   0   0
   1   0   0
   1   0   0
   1   0   0
   1   0   0
   2   0   0
   2   0   0
   2   0   0
   2   0   0
   2   0   0
   2   0   0
   3   0   0
   3   0   0
   3   0   0
   3   0   0
   3   0   0
   3   0   0
   4   0   0
   4   0   0
   4   0   0
   5   0   0
   5   0   0
   7   0   0
   8   0   0
   9   0   0
   9   0   0
  12   0   0
  13   0   0
  14   0   0
  16   0   0
  20   0   0
   0   1   0
   0   1   0
   0   1   0
   0   1   0
   0   1   0
   0   1   0
   0   1   0
   1   1   0
   1   1   0
   1   1   0
   1   1   0
   1   1   0
   1   1   0
   1   1   0
   1   1   0
   1   1   0
   1   1   0
   1   1   0
   1   1   0
   1   1   0
   2   1   0
   2   1   0
   2   1   0
   2   1   0
   2   1   0
   2   1   0
   2   1   0
   3   1   0
   3   1   0
   3   1   0
   4   1   0
   4   1   0
   5   1   0
   5   1   0
   5   1   0
   9   1   0
  10   1   0
  10   1   0
  15   1   0
   0   2   0
   0   2   0
   0   2   0
   0   2   0
   0   2   0
   0   2   0
   0   2   0
   0   2   0
   0   2   0
   1   2   0
   1   2   0
   1   2   0
   1   2   0
   1   2   0
   1   2   0
   1   2   0
   1   2   0
   1   2   0
   2   2   0
   2   2   0
   2   2   0
   3   2   0
   4   2   0
   4   2   0
   5   2   0
   6   2   0
   6   2   0
   6   2   0
   7   2   0
   7   2   0
  12   2   0
  15   2   0
  26   2   0
   0   3   0
   0   3   0
   0   3   0
   0   3   0
   0   3   0
   0   3   0
   0   3   0
   0   3   0
   0   3   0
   1   3   0
   1   3   0
   1   3   0
   1   3   0
   1   3   0
   2   3   0
   2   3   0
   4   3   0
   5   3   0
   7   3   0
  10   3   0
   0   4   0
   0   4   0
   0   4   0
   0   4   0
   1   4   0
   1   4   0
   1   4   0
   1   4   0
   3   4   0
   3   4   0
   0   5   0
   0   6   0
   0   6   0
   1   6   0
   3   6   0
   0   7   0
   0   7   0
   1   7   0
   0   8   0
   1   8   0
   0   9   0
   0   9   0
   0   9   0
   0   0   1
   0   0   1
   0   0   1
   0   0   1
   0   0   1
   0   0   1
   0   0   1
   0   0   1
   1   0   1
   1   0   1
   1   0   1
   1   0   1
   1   0   1
   1   0   1
   1   0   1
   1   0   1
   1   0   1
   1   0   1
   1   0   1
   1   0   1
   2   0   1
   2   0   1
   2   0   1
   2   0   1
   2   0   1
   2   0   1
   3   0   1
   3   0   1
   3   0   1
   3   0   1
   3   0   1
   3   0   1
   3   0   1
   3   0   1
   4   0   1
   4   0   1
   4   0   1
   4   0   1
   4   0   1
   4   0   1
   4 

Re: [R] Object-oriented programming in R for Java programmers?

2008-07-27 Thread Stephen Tucker
This page is also a brief introduction to the S3 classes:

http://www.ibm.com/developerworks/linux/library/l-r3.html

(see section on 'Object-oriented R')

and for S4, in addition to How S4 methods work:

http://developer.r-project.org/methodDefinition.html



- Original Message 
From: Johannes Huesing [EMAIL PROTECTED]
To: r-help@r-project.org
Sent: Sunday, July 27, 2008 12:21:37 PM
Subject: Re: [R] Object-oriented programming in R for Java programmers?

Werner Wernersen [EMAIL PROTECTED] [Sun, Jul 27, 2008 at 01:55:26PM CEST]:
[...]
 Is there any brief tutorial on oo programming in R
 especially for people who have done oo in Java or C++
 before? That would be really helpful.

How S4 methods work highlights the differences between R's function-centric
approach to OO and the class-centric approach highlighted by languages such
as C++. This 10-page article is written by John Chambers.

-- 
Johannes Hüsing   There is something fascinating about science. 
  One gets such wholesale returns of conjecture 
mailto:[EMAIL PROTECTED]  from such a trifling investment of fact.  
  
http://derwisch.wikidot.com (Mark Twain, Life on the Mississippi)

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





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


[R] help with durbin.watson

2008-07-27 Thread tom soyer
Hi,

I have two time series, y and x. Diff(y) and Diff(x) both show no
autocorrelation. But durbin.watson(lm(Diff(y)~lag(Diff(x),k=-4)) gives a DW
value of zero. How come the residule is autocorrelated while Diff(y) and
Diff(x) are not? Does anyone know if in my case a DW of zero indicates
serial correlation, or is it telling me that the DW statistics is not the
appropriate statistics to use here?

Thanks,

-- 
Tom

[[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] graphing regression coefficients and standard errors

2008-07-27 Thread Ritwik Sinha
Hi Murali,

The plot function works on an lm object. The example from the help page of
lm.

Ritwik

ctl - c(4.17,5.58,5.18,6.11,4.50,4.61,5.17,4.53,5.33,5.14)
trt - c(4.81,4.17,4.41,3.59,5.87,3.83,6.03,4.89,4.32,4.69)
group - gl(2,10,20, labels=c(Ctl,Trt))
weight - c(ctl, trt)
anova(lm.D9 - lm(weight ~ group))
summary(lm.D90 - lm(weight ~ group - 1))# omitting intercept
summary(resid(lm.D9) - resid(lm.D90)) #- residuals almost identical

opar - par(mfrow = c(2,2), oma = c(0, 0, 1.1, 0))
plot(lm.D9, las = 1)  # Residuals, Fitted, ...



On Fri, Jul 25, 2008 at 4:01 PM, Murali K [EMAIL PROTECTED] wrote:

 Hello,
I am interested in plotting my regression analysis
 results(regression coefficients and standard errors obtained through OLS
 and
 Tobit models)
 in the form of graphs.What is the best way to accomplish this? Thanks.


 Murali Kuchibhotla

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




-- 
Ritwik Sinha
[EMAIL PROTECTED] | +12033042111 | http://ritwik.sinha.googlepages.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.


Re: [R] Maximization under constraits

2008-07-27 Thread Ritwik Sinha
Hi Daniela,

Will the optim function with the method L-BFGS-B work for you?
Look for the lower argument in the function.

Ritwik

On Fri, Jul 25, 2008 at 9:07 AM, Daniela Garavaglia
[EMAIL PROTECTED] wrote:

 I'm looking for a R function which can maximise this logliklihood function,
 under the constraits a0 e b0



 f-function(param){

  a-param[1]

  b -param[2]

 log(prod)-(a*s2)-(b*s)-n*log(1-((0.5*b/sqrt(a))*(exp((b^2)/(4*a)))*((sqrt(pi
 ))*(1-pnorm(-b/(2*sqrt(a)), mean=0, sd=1)}



 I've tried maxlik constrOptim e donlp2 but without success.



 Thanks so much!!!






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



--
Ritwik Sinha
[EMAIL PROTECTED] | +12033042111 | http://ritwik.sinha.googlepages.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] equivalent R functions for Numerical Recipes fitxy and fitexy ?

2008-07-27 Thread Marc Fischer

Dear Folks,

We need to fit the model y~x  assuming there are random errors in both  
x and y.  Numerical Recipes (Press et al.) have two useful functions,  
fitxy and fitexy, which handle cases of unspecified and specified  
errors respectively.


Are there equivalent functions in base R or a installable package?
Alternatively, has anyone written a wrapper to provide an interface to  
a library of Numerical Recipes routines constructed for Linux?


Kind regards,

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.


[R] Interpolating a line and then summing there values for a diurnal oxygen curve (zoo object)

2008-07-27 Thread stephen sefick
#I would like to interpolate a straight line between 06/08/06  04:16:00 -
06/08/06 20:31:00 with values and then sum them.  This is an estimate of
ecosystem #respiration and I will be using this in a larger context(48 days
of these diurnal curves), but for right now I am just trying to figure out
how to do it for this one #day example.  I have some other code for
Ecosystem (stream) Metabolism that I am trying to work up into a package,
but I am still assembling the pieces.
#thanks for all of your help

library(chron)
library(zoo)
z - 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)
plot(z)

-- 
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] Interpolating a line and then summing there values for a diurnal oxygen curve (zoo object)

2008-07-27 Thread Gabor Grothendieck
Try this using your z and library(zoo):

# extract section of interest
w - window(z, 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(z, w, lin), col = 1:3, screen = 1)


On Sun, Jul 27, 2008 at 8:12 PM, stephen sefick [EMAIL PROTECTED] wrote:
 #I would like to interpolate a straight line between 06/08/06  04:16:00 -
 06/08/06 20:31:00 with values and then sum them.  This is an estimate of
 ecosystem #respiration and I will be using this in a larger context(48 days
 of these diurnal curves), but for right now I am just trying to figure out
 how to do it for this one #day example.  I have some other code for
 Ecosystem (stream) Metabolism that I am trying to work up into a package,
 but I am still assembling the pieces.
 #thanks for all of your help

 library(chron)
 library(zoo)
 z - 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)
 plot(z)

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


[R] Mixed model question.

2008-07-27 Thread Rolf Turner


I continue to struggle with mixed models.  The square zero version
of the problem that I am trying to deal with is as follows:

A number (240) of students are measured (tested; for reading  
comprehension)

on 6 separate occasions.  Initially (square zero) I want to treat the
test time as a factor (with 6 levels).  The students are of course
``random effects''.  Later I want to look at reduced models, with the
means at the 6 times following patterns:

(mu1, mu2, mu2, mu3, mu3, mu4)

or

(mu1, mu1+theta, mu1+theta, mu1+2*theta, mu1+2*theta, mu1+3*theta)

the idea begin that the tests are given ``in pairs'' at the beginning  
and

end of the school year.  Progress is expected over the school year, no
progress over the two intervening summers.  The summers fall between
times 2 and 3 and between times 4 and 5.  The second of the two models
assumes a constant amount of progress (``theta'') in each of three
school years in which the students are observed.

But the square zero model is just a trivial repeated measures model.

The estimate of the means at the 6 observation times is just
mu-hat = xbar - apply(X,2,mean) where X is the ``data matrix'',
and the estimate of the covariance structure is just given by
Sigma-hat = S - var(X).

No problem so far.  I can also fit the two reduced models (I'm pretty
sure) via maximum likelihood, using optim(), for example.  Assuming
normal data.  Rash assumption, but that's not the issue here.

But the thing is, I also want (later!) to include such things as
a school effect (there are 6 different schools), a sex effect, and
an ethnicity effect.

Things start to get complicated --- sounds like a job for lmer().

So I'd just like to get a toe in the door by fitting the trivial
(square 0) model with lmer() --- and then if I can get my head
around that, move on to the two reduced models.  I.e. I'd like to
reproduce my simple-minded computations using lmer() --- which would
give me a little bit of confidence that I'm driving lmer() correctly.

*Can* the trivial model be fitted in lmer()?  I tried using

fit - lmer(y ~ tstnum + (1|stdnt), data=schooldat)

and got estimates of the coefficients for tstnum as follows:

Fixed effects:
Estimate Std. Error t value
(Intercept)  3.229170.09743   33.14
tstnum2  0.466670.084615.52
tstnum3  0.50.084615.91
tstnum4  0.837500.084619.90
tstnum5  0.470830.084615.56
tstnum6  0.975000.08461   11.52

The mean of (the columns of) the data matrix is

3.229167 3.695833 3.729167 4.07 3.70 4.204167

which is in exact agreement with the lmer() results when converted to
the same parameterization (mu_i = mu + alpha_i, with alpha_1 = 0).

(Notice the surprizing, depressing, and so far unexplained *drop*
in the response over the second summer.)

What I *don't* understand is the correlation structure of the estimates
produced by lmer(), which is:

Correlation of Fixed Effects:
(Intr) tstnm2 tstnm3 tstnm4 tstnm5
tstnum2 -0.434
tstnum3 -0.434  0.500
tstnum4 -0.434  0.500  0.500
tstnum5 -0.434  0.500  0.500  0.500
tstnum6 -0.434  0.500  0.500  0.500  0.500

So apparently the way I called lmer() places substantial constraints
on the covariance structure.  How can I (is there any way that I can)
tell lmer() to fit the most general possible covariance structure?

As usual, advice, insight, tutelage humbly appreciated.

If anyone wishes to experiment with the real data set (it's a bit
too big to post here) I can make it available to them via email.

Thanks.

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] Coefficients of Logistic Regression from bootstrap - how to get them?

2008-07-27 Thread Tim Hesterberg
I'll address the question of whether you can use the bootstrap to
improve estimates, and whether you can use the bootstrap to virtually
increase the size of the sample.

Short answer - no, with some exceptions (bumping / Random Forests).

Longer answer:
Suppose you have data (x1, ..., xn) and a statistic ThetaHat,
that you take a number of bootstrap samples (all of size n) and
let ThetaHatBar be the average of those bootstrap statistics from
those samples.

Is ThetaHatBar better than ThetaHat?  Usually not.  Usually it
is worse.  You have not collected any new data, you are just using the
existing data in a different way, that is usually harmful:
* If the statistic is the sample mean, all this does is to add
  some noise to the estimate
* If the statistic is nonlinear, this gives an estimate that
  has roughly double the bias, without improving the variance.

What are the exceptions?  The prime example is tree models (random
forests) - taking bootstrap averages helps smooth out the
discontinuities in tree models.  For a simple example, suppose that a
simple linear regression model really holds:
y = beta x + epsilon
but that you fit a tree model; the tree model predictions are
a step function.  If you bootstrap the data, the boundaries of
the step function will differ from one sample to another, so
the average of the bootstrap samples smears out the steps, getting
closer to the smooth linear relationship.

Aside from such exceptions, the bootstrap is used for inference
(bias, standard error, confidence intervals), not improving on
ThetaHat.

Tim Hesterberg

Hi Doran,

Maybe I am wrong, but I think bootstrap is a general resampling method which
can be used for different purposes...Usually it works well when you do not
have a presentative sample set (maybe with limited number of samples).
Therefore, I am positive with Michal...

P.S., overfitting, in my opinion, is used to depict when you got a model
which is quite specific for the training dataset but cannot be generalized
with new samples..

Thanks,

--Jerry
2008/7/21 Doran, Harold [EMAIL PROTECTED]:

  I used bootstrap to virtually increase the size of my
  dataset, it should result in estimates more close to that
  from the population - isn't it the purpose of bootstrap?

 No, not really. The bootstrap is a resampling method for variance
 estimation. It is often used when there is not an easy way, or a closed
 form expression, for estimating the sampling variance of a statistic.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] product of successive rows

2008-07-27 Thread rcoder

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] Converting english words to numeric equivalents

2008-07-27 Thread oscar Linares
Hello,

I am trying to convert english words to numeric equivalents, e.g., abc to 123. 
Is there a function or library or package for doing this in R? If not, can it 
be done easily in R?

Thanks,

Oscar



  
[[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] Rolling regression - o/p selected coefficients

2008-07-27 Thread rcoder

Hi everyone,

Using the rolling regression function - rollingRegression(formula, data,
width, ...) - is there a way to output only selected coefficients to a
results matrix. For e.g. if I only wanted the slope coefficients to be
recorded, how would I go about doing this?

Thanks,

Michael
-- 
View this message in context: 
http://www.nabble.com/Rolling-regression---o-p-selected-coefficients-tp18681338p18681338.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] Converting english words to numeric equivalents

2008-07-27 Thread Rolf Turner


On 28/07/2008, at 12:32 PM, oscar Linares wrote:


Hello,

I am trying to convert english words to numeric equivalents, e.g.,  
abc to 123. Is there a function or library or package for doing  
this in R? If not, can it be done easily in R?


Your question is very ill-posed.  What ``english word'' would  
correspond to 0?

What would the ``english word'' xyz (?!?!) correspond to?

If you can phrase a coherent question, it can probably be answered.

cheers,

Rolf Turner

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

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


[R] read.table question

2008-07-27 Thread Edna Bell
Dear R list:

I have a data set that I would like to bring it.  The fourth column
shows as numeric, but I want it to be a factor.  Is there a way to do
this from the read.table statement, or should I just wait and use the
factor function 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] read.table question

2008-07-27 Thread Rolf Turner


On 28/07/2008, at 3:26 PM, Edna Bell wrote:


Dear R list:

I have a data set that I would like to bring it.  The fourth column
shows as numeric, but I want it to be a factor.  Is there a way to do
this from the read.table statement, or should I just wait and use the
factor function please?


Look at the colClasses argument to read.table().

cheers,

Rolf Turner

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

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


[R] Help with a loop

2008-07-27 Thread Sofia Martinez
HI:
I need ideas on how to make this code shorter (maybe with a second loop?).
The code as it is works, but in this case I only have 14 samples, but it
will become insane with more, so I need a way to make it more automatic. The
problem is that the output from ts1, ts2, and so on is a vector with more
than one value, so I do not know how to solve this.
   Thanks
  Prenewbie

The code is the following:

duration-c(750, 8940,17180, 8693,10100, 7990,24820, 6770,
7390,8450,18400,16252,6080,11030)

tmax-0
dt-0
for (m in 1:14)
 {
  tmax[m]-duration[m]
   dtmin-969
   dtmax-9884
   i-1
   dt[m]-round(runif(1,dtmin,(tmax[m]-1)))}
  ts-c(1)
  ts1-ts
  ts2-ts1
  ts3-ts1
  ts4-ts1
  ts5-ts1
  ts6-ts1
  ts7-ts1
  ts8-ts1
  ts9-ts1
  ts10-ts1
  ts11-ts1
  ts12-ts1
  ts13-ts1
  ts14-ts1
  for (m in 1:2){
  while(i+dt[m]tmax[m])
  {
  ts1-append(ts,i+dt[1])
  i-i+dt[1]
  ifelse(dtmin=(tmax[1]-i), dt[1]-round(runif(1,dtmin,(tmax[1]-i))),
dt[1]-(tmax[1]-i))
  ts2-append(ts,i+dt[2])
  i-i+dt[2]
  ifelse(dtmin=(tmax[2]-i), dt[2]-round(runif(1,dtmin,(tmax[2]-i))),
dt[2]-(tmax[2]-i))
   ts3-append(ts,i+dt[3])
  i-i+dt[3]
  ifelse(dtmin=(tmax[3]-i), dt[3]-round(runif(1,dtmin,(tmax[3]-i))),
dt[3]-(tmax[3]-i))
   ts4-append(ts,i+dt[4])
  i-i+dt[4]
  ifelse(dtmin=(tmax[4]-i), dt[4]-round(runif(1,dtmin,(tmax[4]-i))),
dt[4]-(tmax[4]-i))
   ts5-append(ts,i+dt[5])
  i-i+dt[5]
  ifelse(dtmin=(tmax[5]-i), dt[5]-round(runif(1,dtmin,(tmax[5]-i))),
dt[5]-(tmax[5]-i))
   ts6-append(ts,i+dt[6])
  i-i+dt[6]
  ifelse(dtmin=(tmax[6]-i), dt[6]-round(runif(1,dtmin,(tmax[6]-i))),
dt[6]-(tmax[6]-i))
   ts7-append(ts,i+dt[7])
  i-i+dt[7]
  ifelse(dtmin=(tmax[7]-i), dt[7]-round(runif(1,dtmin,(tmax[7]-i))),
dt[7]-(tmax[7]-i))
   ts8-append(ts,i+dt[8])
  i-i+dt[8]
  ifelse(dtmin=(tmax[8]-i), dt[8]-round(runif(1,dtmin,(tmax[8]-i))),
dt[8]-(tmax[8]-i))
   ts9-append(ts,i+dt[9])
  i-i+dt[9]
  ifelse(dtmin=(tmax[9]-i), dt[9]-round(runif(1,dtmin,(tmax[9]-i))),
dt[9]-(tmax[9]-i))
   ts10-append(ts,i+dt[10])
  i-i+dt[10]
  ifelse(dtmin=(tmax[10]-i),
dt[10]-round(runif(1,dtmin,(tmax[10]-i))), dt[10]-(tmax[10]-i))
   ts11-append(ts,i+dt[11])
  i-i+dt[11]
  ifelse(dtmin=(tmax[11]-i),
dt[11]-round(runif(1,dtmin,(tmax[11]-i))), dt[11]-(tmax[11]-i))
   ts12-append(ts,i+dt[12])
  i-i+dt[12]
  ifelse(dtmin=(tmax[12]-i),
dt[12]-round(runif(1,dtmin,(tmax[12]-i))), dt[12]-(tmax[12]-i))
   ts13-append(ts,i+dt[13])
  i-i+dt[13]
  ifelse(dtmin=(tmax[13]-i),
dt[13]-round(runif(1,dtmin,(tmax[13]-i))), dt[13]-(tmax[13]-i))
   ts14-append(ts,i+dt[14])
  i-i+dt[14]
  ifelse(dtmin=(tmax[14]-i),
dt[14]-round(runif(1,dtmin,(tmax[14]-i))), dt[14]-(tmax[14]-i))
  } }
 ts1-append(ts,i+dt[1])
 ts2-append(ts,i+dt[2])
 ts3-append(ts,i+dt[3])
 ts4-append(ts,i+dt[4])
 ts5-append(ts,i+dt[5])
 ts6-append(ts,i+dt[6])
 ts7-append(ts,i+dt[7])
 ts8-append(ts,i+dt[8])
 ts9-append(ts,i+dt[9])
 ts10-append(ts,i+dt[10])
 ts11-append(ts,i+dt[11])
 ts12-append(ts,i+dt[12])
 ts13-append(ts,i+dt[13])
 ts14-append(ts,i+dt[14])

[[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] Are there any packages that can process images other than pixelmap (i.e. pnm)?

2008-07-27 Thread Arthur Roberts

Hi, all,

I am having trouble getting R to take pnm images via mogrify

i.e.

mogrify -resize 320x217 -format pnm *.png

However R via pixmap says that it can't read the file.  If you have  
any ideas like a package that can read jpeg files, etc., I would  
appreciate it.


Best wishes,
Art Roberts
University of Washington

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Are there any packages that can process images other than pixelmap (i.e. pnm)?

2008-07-27 Thread Henrik Bengtsson
See EBImage on Bioconductor.org. /Henrik

On Sun, Jul 27, 2008 at 10:21 PM, Arthur Roberts
[EMAIL PROTECTED] wrote:
 Hi, all,

 I am having trouble getting R to take pnm images via mogrify

 i.e.

 mogrify -resize 320x217 -format pnm *.png

 However R via pixmap says that it can't read the file.  If you have any
 ideas like a package that can read jpeg files, etc., I would appreciate it.

 Best wishes,
 Art Roberts
 University of Washington

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/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.