Re: [R] getting started on Bayesian analysis

2004-09-15 Thread Jari Oksanen
On Wed, 2004-09-15 at 03:27, HALL, MARK E wrote:
   I've found 
   
 Bayesian Methods: A Social and Behavioral Sciences Approach
 by Jeff Gill 
 
 useful as an introduction.  The examples are written in R and S with generalized 
 scripts for doing 
 a variety of problems.  (Though I never got change-point analysis to successfully in 
 R.)
 
Change point analysis? I haven't seen the book, but I read lecture
handouts of one Bayesian course over here in Finland (Antti Penttinen,
Jyväskylä), and translated his example to R during one (rare) warm
summer day in a garden. So do you mean this (binary case):

 source(/mnt/flash/cb.update.R)
 cb.update
function (y, A=1, B=1, C=1, D=1, N=1200, burnin=200)
{
n - length(y)
lambda - numeric(N)
mu - numeric(N)
k - numeric(N)
lambda[1] - A/(A+B)
mu[1] - C/(C+D)
k[1] - n/2
sn - sum(y)
 
for (i in 2:N) {
kold - k[i-1]
sk - sum(y[1:kold])
lambda[i] - rbeta(1, A+sk, B + kold - sk)
mu[i] - rbeta(1, C + sn - sk, D + n - sn + sk - kold  )
knew - sample(n-1, 1)
sknew - sum(y[1:knew])
r - (sknew - sk) *
(log(lambda[i])-log(mu[i]))-(knew-kold)*(lambda[i]-mu[i])
if(min(0,r)  log(runif(1))) k[i] - knew
else k[i] - k[i-1]
}
out - cbind(lambda, mu, k)
out[(burnin+1):N, ]
}
 y - c(rbinom(60, 1, 0.8), rbinom(40, 1, 0.3))
 uh - cb.update(y, N=5200)
 colMeans(uh)
lambda mu  k
 0.8189303  0.4169367 59.077
 mean(y[1:60])
[1] 0.783
 mean(y[41:100])
[1] 0.45
 plot(density(uh[,1]))
 plot(density(uh[,2]))
 plot(table(uh[,3]), type=h)

This was off-topic. So something about business: isn't the (Win)BUGS
author working with a R port?

cheers, jari oksanen
-- 
Jari Oksanen -- Dept Biology, Univ Oulu, 90014 Oulu, Finland
email [EMAIL PROTECTED], homepage http://cc.oulu.fi/~jarioksa/

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


Re: [R] getting started on Bayesian analysis

2004-09-15 Thread Anon.
Jari Oksanen wrote:
On Wed, 2004-09-15 at 03:27, HALL, MARK E wrote:
 I've found 
 
Bayesian Methods: A Social and Behavioral Sciences Approach
by Jeff Gill 

useful as an introduction.  The examples are written in R and S with generalized scripts for doing 
a variety of problems.  (Though I never got change-point analysis to successfully in R.)

Change point analysis? I haven't seen the book, but I read lecture
handouts of one Bayesian course over here in Finland (Antti Penttinen,
Jyväskylä), and translated his example to R during one (rare) warm
summer day in a garden. So do you mean this (binary case):

snip

This was off-topic. So something about business: isn't the (Win)BUGS
author working with a R port?
Yes, he is, and he's got it working in Windows.  But if anyone wants to 
discuss how R does memory management with Andrew, he'll be all ears.

In the mean time, there is the R2WinBUGS package on CRAN for you to enjoy.
Bob
--
Bob O'Hara
Dept. of Mathematics and Statistics
P.O. Box 68 (Gustaf Hällströmin katu 2b)
FIN-00014 University of Helsinki
Finland
Telephone: +358-9-191 51479
Mobile: +358 50 599 0540
Fax:  +358-9-191 51400
WWW:  http://www.RNI.Helsinki.FI/~boh/
Journal of Negative Results - EEB: http://www.jnr-eeb.org
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Cancor

2004-09-15 Thread Irena Komprej
Dear Gabor,
thank you for your answer related to the normalization of cancor
coefficients. If you want to interpret the coefficients in terms of
variables' contributions to canonical variables, loadings, redundancy
measure, etc. ,  you have to normalize the results so that the canonical
variables have identity variance matrix. Multiplying cancor coefficients
xcoef and ycoef by as.numeric(sqrt(nrow(x)-1)) does the job.
(I sometimes miss such information in the R help.)
Best regards,
Irena Komprej

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


[R] Bessel function

2004-09-15 Thread beat.huggler


 Dear all

Currently, I'm implementing the generalized hyperbolic distribution into
Splus. Unfortunately the Bessel function is not implemented in Splus. In
R the Bessel function does exist but it is an internal function and I'm
not able to look at the code.

Is there any possibility to see the code of the Bessel function in R or
does anybody has an implementation of the Bessel function in Splus?

Thanks a lot for your help.

Beat

[[alternative HTML version deleted]]

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


Re: [R] Bessel function

2004-09-15 Thread Dimitris Rizopoulos
Hi Beat,

take a look at this link, http://www.statsci.org/s/besseli0.html for
the modified Bessel function.

I hope it helps.

Best,
Dimitris


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

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


- Original Message - 
From: [EMAIL PROTECTED]
To: [EMAIL PROTECTED]
Sent: Wednesday, September 15, 2004 10:24 AM
Subject: [R] Bessel function




  Dear all

 Currently, I'm implementing the generalized hyperbolic distribution
into
 Splus. Unfortunately the Bessel function is not implemented in
Splus. In
 R the Bessel function does exist but it is an internal function and
I'm
 not able to look at the code.

 Is there any possibility to see the code of the Bessel function in R
or
 does anybody has an implementation of the Bessel function in Splus?

 Thanks a lot for your help.

 Beat

 [[alternative HTML version deleted]]

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


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


Re: [R] Cancor

2004-09-15 Thread Prof Brian Ripley
On (she claimed) Thu, 1 Jan 1998, Irena Komprej wrote:

 Dear Gabor,

This is R-help, not `Gabor', although you keep sending mail addressed to
`Gabor' to this list.

 thank you for your answer related to the normalization of cancor
 coefficients. If you want to interpret the coefficients in terms of
 variables' contributions to canonical variables, loadings, redundancy
 measure, etc. ,  you have to normalize the results so that the canonical
 variables have identity variance matrix. Multiplying cancor coefficients
 xcoef and ycoef by as.numeric(sqrt(nrow(x)-1)) does the job.
 (I sometimes miss such information in the R help.)

I think you miss it in the references given on the R help pages.  Please
do read them.  (Any good book will tell you that the scaling of canonical 
variates is arbitrary.)

Please also ask your local IT advisors to set your computer to a 
sensible time: the world is 6.7 years ahead of you.

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

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


Re: [R] Bessel function

2004-09-15 Thread Prof Brian Ripley
On Wed, 15 Sep 2004 [EMAIL PROTECTED] wrote:

 Currently, I'm implementing the generalized hyperbolic distribution into
 Splus. Unfortunately the Bessel function is not implemented in Splus. In
 R the Bessel function does exist but it is an internal function and I'm
 not able to look at the code.
 
 Is there any possibility to see the code of the Bessel function in R or

Of course.  It *is* in the sources.  What stops you looking at the 
sources?  (Where they are is in the FAQ, for example.)

BTW, there are four Bessel functions in R, and others do exist.

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

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


[R] testing goodness of fit of linear model

2004-09-15 Thread Dewez Thomas
Dear R-users,

I've been reading a bunch of things on linear models but cannot quite find a
clear answer. How can one determine whether a linear model is significant or
not?

For background info, I am modelling the response of topographic slope to the
distance of a catchment's outlet. Some guys have shown that if there is a
significant fit to a linear model, one can deduce the dynamic state of the
basin, that is, whether erosion is as strong as rock uplift, erosion is
smaller than rock uplift, or erosion is greater than rock uplift. I am thus
to test 4 situations:

Situation 1: a linear model is inappropriate for describing the data, the
scatter is too large, and thus a linear model is unfit to explain the data.

Situation 2: the linear model of the kind y = b0 + b1 * x is fit to
describe the data, ie data points lie close to a straight line.

Situation 2a: the relationship between slope and distance is significantly
positive
Situation 2b: the relationship between slope and distance is significantly
null (ie data is clustered around a line with b1 non-significantly different
from 0)
Situation 2c: the relationship between slope and distance is significantly
negative

I am confused as to what test I should use for discriminating these
situations.

The glm offers an indication about the significance of regression
parameters. So in the case where b1 is significantly different from 0 (p
value =0.05 for a test where H0: b1=0; H1: b1 != 0), it is straightforward.
But I don't know how to discriminate between situation 1 and situation 2 (ie
whether a linear model is significant).

Any suggestion are welcome

Cheers,

Thomas
***
Le contenu de cet e-mail et de ses pièces jointes est destin...{{dropped}}

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


Re: [R] Bessel function

2004-09-15 Thread Robin Hankin
Hi
you might find package(hyperbolicDist) interesting.
best wishes
rksh

 Dear all
Currently, I'm implementing the generalized hyperbolic distribution into
Splus. Unfortunately the Bessel function is not implemented in Splus. In
R the Bessel function does exist but it is an internal function and I'm
not able to look at the code.
Is there any possibility to see the code of the Bessel function in R or
does anybody has an implementation of the Bessel function in Splus?
Thanks a lot for your help.
Beat
[[alternative HTML version deleted]]
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

--
Robin Hankin
Uncertainty Analyst
Southampton Oceanography Centre
SO14 3ZH
tel +44(0)23-8059-7743
[EMAIL PROTECTED] (edit in obvious way; spam precaution)
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Bessel function

2004-09-15 Thread Robin Hankin
Hi
you might find package(hyperbolicDist) interesting.
best wishes
robin

 Dear all
Currently, I'm implementing the generalized hyperbolic distribution into
Splus. Unfortunately the Bessel function is not implemented in Splus. In
R the Bessel function does exist but it is an internal function and I'm
not able to look at the code.
Is there any possibility to see the code of the Bessel function in R or
does anybody has an implementation of the Bessel function in Splus?
Thanks a lot for your help.
Beat
[[alternative HTML version deleted]]
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

--
Robin Hankin
Uncertainty Analyst
Southampton Oceanography Centre
SO14 3ZH
tel +44(0)23-8059-7743
[EMAIL PROTECTED] (edit in obvious way; spam precaution)
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] loading error of the Rcmdr library on Debian Sid

2004-09-15 Thread Christian Schulz
Hi,

it seems you have to install the rgl package and this demand 
the gl4 java lib what is necessary, too. 

http://www.jausoft.com/products/gl4java/gl4java_install.html

christian


Am Mittwoch, 15. September 2004 12:01 schrieb Thomas Schönhoff:
 Hello,

 I just tried to get Rcmdr package working, resulting in:


 --

   library(Rcmdr)

 Loading required package: tcltk
 Loading required package: lattice
 Loading required package: foreign
 Loading required package: abind
 Loading required package: lmtest
 Loading required package: multcomp
 Loading required package: relimp
 Loading required package: effects
 Loading required package: rgl
 Error in dyn.load(x, as.logical(local), as.logical(now)) :
  unable to load shared library
 /usr/lib/R/site-library/rgl/libs/rgl.so:
libnvidia-tls.so.1: cannot handle TLS data
 Error in .C(symbol.C(rgl_quit), success = FALSE, PACKAGE = rgl) :
  C function name not in DLL for package rgl
 Loading required package: mgcv
 This is mgcv 1.1-1
 Loading required package: car
 Error: Missing packages: rgl
 Error: .onLoad failed in loadNamespace
 Error in library(Rcmdr) : package/namespace load failed

 Missing rgl-package ?


 ---


 thomas dpkg -l |grep r-

 r-cran-abind   1.1.0-1
 r-cran-car 1.0.13-1
 r-cran-foreign 0.7-1
 r-cran-lattice 0.9.16-1
 r-cran-mgcv1.1.1.1-1
 r-cran-rgl 0.64.13-1
 r-cran-relimp  0.8.4-1
 r-cran-rcmdr   0.9.11-1
 r-cran-lmtest  0.9.6-2
 r-cran-effects 1.0.5-1
 r-cran-multcom 0.4.7-1
 r-cran-mvtnorm 0.6.8-1


 amongst other R related (basic and specific) packages.


 Am I still missing some required packages to run RCommander smoothly?


 My system:

 Debian Sid (unstable)
 GNU R 1.9.1
 Xfree 4.3

 Anyone else noticed this on Debian Sid whilst trying to run the R
 Commander !? So far I didn't find a related bug report for Linux. Maybe
 this is a Debian related problem, I really have no clue at the moment.

 Regards
 Thomas

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

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


Re: [R] loading error of the Rcmdr library on Debian Sid

2004-09-15 Thread Thomas Schönhoff
Hello Christian,
Christian Schulz schrieb:
Hi,
it seems you have to install the rgl package and this demand 
the gl4 java lib what is necessary, too. 

http://www.jausoft.com/products/gl4java/gl4java_install.html
 

Seems like that this dependency is missing in the r-cran-rgl package!? 
So I am going to write a bug report to the debian maintainer!

Thanks
Thomas
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Density Estimation

2004-09-15 Thread Brian Mac Namee
Hi there,

Sorry if this is a rather loing post. I have a simple list of single
feature data points from which I would like to generate a probability
that an unseen point comes from the same distribution. To do this I am
trying to estimate the probability density of the list of points and
use this to generate a probability for the new unseen points. I have
managed to use the R density function to generate the density estimate
but have not been able to do anything with this - i.e. generate a
rpobability that a new point comes from the same distribution. Is
there a function to do this, or am I way off the mark using the
density function at all?

Thanks in advance,

Brian.

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


[R] Read.fwf

2004-09-15 Thread Doran, Harold
Dear List

I have a fixed width file with variables of varying width. The help is
pretty transparent for this feature, but I can't seem to figure out how
I can make effective use of the package with my data.

In my dataset, the first 80 columns are of width 1 followed by other
variables with width larger than 1.

I think the correct way to do this, by brute force, would be

 read.fwf(dataset, 1,1,1,...,1,5,...).

I'm sure I do not need to actually enter 80 1s, correct? How might I
do this more effectively?

Harold

[[alternative HTML version deleted]]

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


[R] Density Estimation

2004-09-15 Thread Vito Ricci
Dear Brian,

I can suggest you to use density() function to get an
estimate of the pdf you're finding (I believe it's
unknown). Then you can plot the point you got by
density() using plot(). In this way you have a graphic
representation of you unknown pdf. According its shape
and helping by the graphic you could try to understand
what kind of pdf it would be (normal, gamma, weibul,
etc.)
After you can estimate parameters of pdf using your
data with LS or ML methods.
Then you can calculate the goodness of fit for each
model of pdf and use the best one.

I hope I get you a little help.

Cordially
Vito Ricci

[EMAIL PROTECTED]  wrote:

Hi there,

Sorry if this is a rather loing post. I have a simple
list of single
feature data points from which I would like to
generate a probability
that an unseen point comes from the same distribution.
To do this I am
trying to estimate the probability density of the list
of points and
use this to generate a probability for the new unseen
points. I have
managed to use the R density function to generate the
density estimate
but have not been able to do anything with this - i.e.
generate a
rpobability that a new point comes from the same
distribution. Is
there a function to do this, or am I way off the mark
using the
density function at all?

Thanks in advance,

Brian.

=
Diventare costruttori di soluzioni

Visitate il portale http://www.modugno.it/
e in particolare la sezione su Palese http://www.modugno.it/archivio/cat_palese.shtml



___

http://it.seriea.fantasysports.yahoo.com/

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


Re: [R] Read.fwf

2004-09-15 Thread Chuck Cleland
  You can use rep() in specifying the widths argument to avoid typing 
all those ones.  For example:

read.fwf(file=c:\\mydata.dat, widths = c(rep(1,80), 5, 4, 6))
hope this helps,
Chuck
Doran, Harold wrote:
I have a fixed width file with variables of varying width. The help is
pretty transparent for this feature, but I can't seem to figure out how
I can make effective use of the package with my data.
In my dataset, the first 80 columns are of width 1 followed by other
variables with width larger than 1.
I think the correct way to do this, by brute force, would be
read.fwf(dataset, 1,1,1,...,1,5,...).
I'm sure I do not need to actually enter 80 1s, correct? How might I
do this more effectively?

--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 452-1424 (M, W, F)
fax: (917) 438-0894
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Read.fwf

2004-09-15 Thread Uwe Ligges
Doran, Harold wrote:
Dear List
I have a fixed width file with variables of varying width. The help is
pretty transparent for this feature, but I can't seem to figure out how
I can make effective use of the package with my data.
In my dataset, the first 80 columns are of width 1 followed by other
variables with width larger than 1.
I think the correct way to do this, by brute force, would be

read.fwf(dataset, 1,1,1,...,1,5,...).

I'm sure I do not need to actually enter 80 1s, correct? How might I
do this more effectively?
?read.fwf:
widths  integer vector, giving the widths of the fixed-width fields (of 
one line).

Hence:
read.fwf(dataset, c(rep(1, 80), 5, ...))
Uwe Ligges


Harold
[[alternative HTML version deleted]]
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


RE: [R] Bessel function

2004-09-15 Thread Ravi Varadhan
Dear Beat,
 
You might also want to look at the book by Zang and Jin -  Computation of Special 
Functions, John Wiley.  They have Fortran sources for all the special functions 
covered in there.
 
Ravi.
 
Ravi Varadhan, Ph.D.
Assistant Professor,  The Center on Aging and Health
Division of Geriatric Medicine and Gerontology,
Johns Hopkins University,
2024 E. Monument Street, Suite 2-700
Baltimore, MD 21205.
(410) 502 - 7806.



From: [EMAIL PROTECTED] on behalf of [EMAIL PROTECTED]
Sent: Wed 9/15/2004 4:24 AM
To: [EMAIL PROTECTED]
Subject: [R] Bessel function





 Dear all

Currently, I'm implementing the generalized hyperbolic distribution into
Splus. Unfortunately the Bessel function is not implemented in Splus. In
R the Bessel function does exist but it is an internal function and I'm
not able to look at the code.

Is there any possibility to see the code of the Bessel function in R or
does anybody has an implementation of the Bessel function in Splus?

Thanks a lot for your help.

Beat

[[alternative HTML version deleted]]

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



[[alternative HTML version deleted]]

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


[R] control of font size colour for title, subtitles, axis, and tick marks in LATTICE graph

2004-09-15 Thread Jens Hainmueller
Hi,

I very much appreciate any help on this fine tuning problem in a lattice
graph (I am new to LATTICE and could not find an example in the help files
that worked for me. My apologies if I missed it there).

I am running the following box plots to compare conditional distributions of
x at different levels of y under two treatment conditions ID=1 (upper
panel )  ID=0 (lower panel of the plot).

bwplot(HF.ELECYEAR ~ stparvotech | ID ,
data=data, aspect=1,
layout=c(1,2),
xlab=Changes in Party Vote Shares,
xlim=(-20:20),
ylab=Periods Following Last Federal Election,
main=Divided Government,
panel = function(x,y)
{
panel.bwplot(x,y)
panel.abline(v=0, col=red)
}
)

How can I:
1. Control the font size of the main title, the panel titles, the axis, and
the tick marks? The usual cex.main=1/3, cex.xlab etc. do not work.
2. Change the color of the boxes where the panel titles (ID=1  ID=0) are
located?

Thank you very much.

Best,
Jens

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


Re: [R] pairs correlations colors

2004-09-15 Thread Uwe Ligges
Valeria Edefonti wrote:
I have the following problem.
I want to use pairs function and get a matrix of scatterplots with the  
correlations in the upper panel and the ordinary scatterplots in the  
lower panel.
Moreover, I want to have points colored in five differet ways in the  
lower panel, because I have five subgroups.
In order to do that I tried to combine examples on pairs function help.
I got a colored matrix using hints on iris dataset.
I got a black and white matrix with correlations using function  
panel.cor, exactly as it is in the example.
Unfortunately, the line:

jpeg(filename=/home/valeria/Thesis/lung/fig/scatterplotcolnames.jpg)
pairs(aggiunta[,1: 
6],labels=c(ALCAM,ITGB5,MSN,CSTB,DHCR24,TRIM29), main =  
Scatterplots selected genes,pch=21,
bg = c(red, green3, blue,  
brown,orange)[aggiunta[,7]],upper.panel=panel.cor)
dev.off()

doesn't allow me to get the desidered matrix with colors and  correlations.
I also tried to create a function panel.col for the lower.panel:
## put colors on the lower panels
 panel.col - function(datiepheno)
 {
usr - par(usr); on.exit(par(usr))
par(bg = c(red, green3, blue,  
brown,orange)[datiepheno[,7]],pch=21, usr = c(0, 1, 0, 1))
}

but it doesn't work as well.
Any idea?
I hope I'll be precise but not too much precise!
Thank you very much
Valeria

Looks like there was no answer yet:
Please specify am example with data available in R, so that we can 
reproduce your example. This will save much time for those people who 
are going to help.
Rearranging your examples with data we do not have available requires 
quite a lot of effort.

What I guess is that you are going to specify
lower.panel = function(x, y)
  points(x, y,
col = sapply(aggiunta[,7], switch, 1=red, 2=green3, .))
Uwe Ligges


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


Re: [R] loading error of the Rcmdr library on Debian Sid

2004-09-15 Thread A.J. Rossini

You can apt-get RGL in sid.  

Thomas Schönhoff [EMAIL PROTECTED] writes:

 Hello,

 I just tried to get Rcmdr package working, resulting in:


 --
   library(Rcmdr)
 Loading required package: tcltk
 Loading required package: lattice
 Loading required package: foreign
 Loading required package: abind
 Loading required package: lmtest
 Loading required package: multcomp
 Loading required package: relimp
 Loading required package: effects
 Loading required package: rgl
 Error in dyn.load(x, as.logical(local), as.logical(now)) :
  unable to load shared library
  /usr/lib/R/site-library/rgl/libs/rgl.so:
libnvidia-tls.so.1: cannot handle TLS data
 Error in .C(symbol.C(rgl_quit), success = FALSE, PACKAGE = rgl) :
  C function name not in DLL for package rgl
 Loading required package: mgcv
 This is mgcv 1.1-1
 Loading required package: car
 Error: Missing packages: rgl
 Error: .onLoad failed in loadNamespace
 Error in library(Rcmdr) : package/namespace load failed

 Missing rgl-package ?


 ---


 thomas dpkg -l |grep r-

 r-cran-abind   1.1.0-1
 r-cran-car 1.0.13-1
 r-cran-foreign 0.7-1
 r-cran-lattice 0.9.16-1
 r-cran-mgcv1.1.1.1-1
 r-cran-rgl 0.64.13-1
 r-cran-relimp  0.8.4-1
 r-cran-rcmdr   0.9.11-1
 r-cran-lmtest  0.9.6-2
 r-cran-effects 1.0.5-1
 r-cran-multcom 0.4.7-1
 r-cran-mvtnorm 0.6.8-1


 amongst other R related (basic and specific) packages.


 Am I still missing some required packages to run RCommander smoothly?


 My system:

 Debian Sid (unstable)
 GNU R 1.9.1
 Xfree 4.3

 Anyone else noticed this on Debian Sid whilst trying to run the R
 Commander !? So far I didn't find a related bug report for
 Linux. Maybe this is a Debian related problem, I really have no clue
 at the moment.

 Regards
 Thomas

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


-- 
Anthony Rossini Research Associate Professor
[EMAIL PROTECTED]http://www.analytics.washington.edu/ 
Biomedical and Health Informatics   University of Washington
Biostatistics, SCHARP/HVTN  Fred Hutchinson Cancer Research Center
UW (Tu/Th/F): 206-616-7630 FAX=206-543-3461 | Voicemail is unreliable
FHCRC  (M/W): 206-667-7025 FAX=206-667-4812 | use Email

CONFIDENTIALITY NOTICE: This e-mail message and any attachme...{{dropped}}

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


RE: [R] Bessel function

2004-09-15 Thread Kahra Hannu
Currently, I'm implementing the generalized hyperbolic distribution into
Splus. Unfortunately the Bessel function is not implemented in Splus. In
R the Bessel function does exist but it is an internal function and I'm
not able to look at the code.

Is there any possibility to see the code of the Bessel function in R or
does anybody has an implementation of the Bessel function in Splus?

Have a look at Press-Teukolsky-Vetterling-Flannery: Numerical Recipes in C, Cambridge 
University Press, Chapter 6 (2nd edition). I guess you are looking for modified Bessel 
functions (see, 6.6).


Hannu Kahra 
Progetti Speciali 
Monte Paschi Asset Management SGR S.p.A. 
Via San Vittore, 37
IT-20123 Milano, Italia 

Tel.: +39 02 43828 754 
Mobile: +39 333 876 1558 
Fax: +39 02 43828 247 
E-mail: [EMAIL PROTECTED] 
Web: www.mpsam.it

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


Re: [R] Density Estimation

2004-09-15 Thread Bob Wheeler
Try fitting it with a Johnson function -- see SuppDists. If you can fit 
it you will then be able to use the functions in SuppDists just as you 
can for any other distribution supported by R.

Brian Mac Namee wrote:
Hi there,
Sorry if this is a rather loing post. I have a simple list of single
feature data points from which I would like to generate a probability
that an unseen point comes from the same distribution. To do this I am
trying to estimate the probability density of the list of points and
use this to generate a probability for the new unseen points. I have
managed to use the R density function to generate the density estimate
but have not been able to do anything with this - i.e. generate a
rpobability that a new point comes from the same distribution. Is
there a function to do this, or am I way off the mark using the
density function at all?
Thanks in advance,
Brian.
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

--
Bob Wheeler --- http://www.bobwheeler.com/
ECHIP, Inc. ---
Randomness comes in bunches.
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


RE: [R] Read.fwf

2004-09-15 Thread Liaw, Andy
It might help to read that help file again:

 d - paste(c(rep(0, 80), 1), collapse=)
 d
[1]
000
01
 f - file(try.dat, w)
 writeLines(d, f)
 writeLines(d, f)
 writeLines(d, f)
 close(f)
 fw - c(rep(1, 80), 5)
 x - read.fwf(try.dat, fw)
 x
  V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 V11 V12 V13 V14 V15 V16 V17 V18 V19 V20 V21
V22
1  0  0  0  0  0  0  0  0  0   0   0   0   0   0   0   0   0   0   0   0   0
0
2  0  0  0  0  0  0  0  0  0   0   0   0   0   0   0   0   0   0   0   0   0
0
3  0  0  0  0  0  0  0  0  0   0   0   0   0   0   0   0   0   0   0   0   0
0
  V23 V24 V25 V26 V27 V28 V29 V30 V31 V32 V33 V34 V35 V36 V37 V38 V39 V40
V41 V42
1   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
0   0
2   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
0   0
3   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
0   0
  V43 V44 V45 V46 V47 V48 V49 V50 V51 V52 V53 V54 V55 V56 V57 V58 V59 V60
V61 V62
1   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
0   0
2   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
0   0
3   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
0   0
  V63 V64 V65 V66 V67 V68 V69 V70 V71 V72 V73 V74 V75 V76 V77 V78 V79 V80
V81
1   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
1
2   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
1
3   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0
1

Andy
 

 From: Doran, Harold
 
 Dear List
 
 I have a fixed width file with variables of varying width. The help is
 pretty transparent for this feature, but I can't seem to 
 figure out how
 I can make effective use of the package with my data.
 
 In my dataset, the first 80 columns are of width 1 followed by other
 variables with width larger than 1.
 
 I think the correct way to do this, by brute force, would be
 
  read.fwf(dataset, 1,1,1,...,1,5,...).
 
 I'm sure I do not need to actually enter 80 1s, correct? How might I
 do this more effectively?
 
 Harold
 
   [[alternative HTML version deleted]]
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 


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


Re: [R] Density Estimation

2004-09-15 Thread Wolski
Hi!

The function density returns you a object of class density.
This object has an x and an y attribute which you can access by x y,
Hi!

Use approx and runif.

eg.:

dd-density(rnorm(100,3,5))
plot(dd)

Using the function ?approx you can compute the density value for any x.
#the x is a dummy here.
mydist-function(x,dd)
{

while(1)
{
tmp - runif(1,min=min(dd$x),max=max(dd$x))
lev - approx(dd$x,dd$y,tmp)$y
if(runif(1,c(0,1)) = lev)
{
return(tmp)
}
}
}

x - 0
mydist(x,dd)

res-rep(0,500)
res-sapply(res,mydist,dd)
lines(density(res),col=2)


/E.



*** REPLY SEPARATOR  ***

On 9/15/2004 at 12:36 PM Brian Mac Namee wrote:

Hi there,

Sorry if this is a rather loing post. I have a simple list of single
feature data points from which I would like to generate a probability
that an unseen point comes from the same distribution. To do this I am
trying to estimate the probability density of the list of points and
use this to generate a probability for the new unseen points. I have
managed to use the R density function to generate the density estimate
but have not been able to do anything with this - i.e. generate a
rpobability that a new point comes from the same distribution. Is
there a function to do this, or am I way off the mark using the
density function at all?

Thanks in advance,

Brian.

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



Dipl. bio-chem. Witold Eryk Wolski @ MPI-Moleculare Genetic   
Ihnestrasse 63-73 14195 Berlin'v'
tel: 0049-30-83875219/   \   
mail: [EMAIL PROTECTED]---W-Whttp://www.molgen.mpg.de/~wolski 
  [EMAIL PROTECTED]

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


Re: [R] loading error of the Rcmdr library on Debian Sid

2004-09-15 Thread Thomas Schnhoff
Hello,
A.J. Rossini schrieb:
You can apt-get RGL in sid.  
An apt-cache search or (Synaptic search) RGL gives me r-cran-rgl 
0.64.13-1. This package is already installed!

regards
Thomas
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] fractal dimension of an image

2004-09-15 Thread Rajarshi Guha
Hello, I have a problem that I think can be solved in R but I'm not sure
how to tie things together.

I have a digital image of a crystal growth in 2 dimensions. And my aim
is to calculate the fractal dimension of the crystal. I was planning to
use the box counting method.

So I need to read in the image in R (for which I can use the pixmap or
rimage packages) and then draw a grid over at a series of resolutions
and for each grid, count how many cells of the grid are occupied by the
crystal. Since the image can be converted to grayscale, I figured that
specifying an intensity threshold would allow me to differentiate
between crystal and surrounding solution.

I know this can be done in Matlab, but since I'm more comfortable in R
can it be done in R? I found the fdim package but that calculates the
dimension for a data.frame - but I'm not sure whether this would be
feasible for my image.

Any pointers would be appreciated
---
Rajarshi Guha [EMAIL PROTECTED] http://jijo.cjb.net
GPG Fingerprint: 0CCA 8EE2 2EEB 25E2 AB04 06F7 1BB9 E634 9B87 56EE
---
So the Zen master asked the hot-dog vendor, 
Can you make me one with everything?
- TauZero on Slashdot

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


RE: [R] Density Estimation

2004-09-15 Thread Ted Harding
On 15-Sep-04 Brian Mac Namee wrote:
 Sorry if this is a rather loing post. I have a simple list of single
 feature data points from which I would like to generate a probability
 that an unseen point comes from the same distribution. To do this I am
 trying to estimate the probability density of the list of points and
 use this to generate a probability for the new unseen points. I have
 managed to use the R density function to generate the density estimate
 but have not been able to do anything with this - i.e. generate a
 rpobability that a new point comes from the same distribution. Is
 there a function to do this, or am I way off the mark using the
 density function at all?

It's not clear what you're really after, but it looks as though you
may be wanting to sample from the distribution estimated by 'density'.

A possible approach, which you could refine, is exemplified by

  x-rnorm(1000)
  d-density(x,n=4096)
  y-sample(d$x,size=1000,prob=d$y)

Check performance with

  hist(y)

Looks OK to me! See ?density and ?sample.

On an alternative interpretation, perhaps you want to first estimate
the density based on data you already have, and then when you have
got further data (but these would then be seen and not unseen)
come to a judgement about whether these new points are compatible
with coming from the distributikon you have estimated.

A possible approach to this question (again susceptible to refinement)
would be as follows.

1. Use a fine-grained grid for 'density', i.e. a large value for n.

2. Replace each of the points in the new data by the nearest point
   in this grid. Call these values z1, z2, ... , zk corresponding
   to index values i1, i2, ... , ik in d$x.

3. Evaluate the probability P(z1,...,zk) from the density as the
   product of d$y[i] where i-c(i1,...,ik).
   Better still, evaluated the logarithm of this. Call the result L.

4. Now simulate a large number of draws of k values from d on the
   lines of sample(d$x,size=k,prob=d$y) as above, and evaluate L
   for each  of these. Where is the value of L from (3) situated in
   the distribution of these values of L from (4)? If (say) only
   1 per cent of the simulated values of L from d are less than
   the value of L from (3), then you have a basis for a test that
   your new data did not come from the distribution you have estimated
   from your old data, in that the new data are from the low-density
   part of the estimated distribution.

There are of course alternative ways to view this question. The
value of k is relevant. In particular, if k is small (say 3
or 4) then the suggestion in (4) is probably the best way to
approach it. However, if k is large then you can use a test on
the lines of Kolmogorov-Smirnov with the reference distribution
estimated as the cumulative distribution of d$y and the distribution
being tested as the empirical cumulative distribution of your new
data.

Even sharper focus is available if you are in a position to make
a paramatric model for your data, but your description does not
suggest that this is the case.

Best wishes,
Ted.



E-Mail: (Ted Harding) [EMAIL PROTECTED]
Fax-to-email: +44 (0)870 167 1972
Date: 15-Sep-04   Time: 15:07:33
-- XFMail --

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


Re: [R] control of font size colour for title, subtitles, axis, and tick marks in LATTICE graph

2004-09-15 Thread Deepayan Sarkar
On Wednesday 15 September 2004 08:24, Jens Hainmueller wrote:
 Hi,

 I very much appreciate any help on this fine tuning problem in a
 lattice graph (I am new to LATTICE and could not find an example in
 the help files that worked for me. My apologies if I missed it
 there).

Well, if the examples were to cover all possible uses, there would be so 
many that they would become practically useless as reference. You are 
also expected to read the documentation. At least for your first 
question, the solution is described pretty much where you would expect 
it to be (see below).

 I am running the following box plots to compare conditional
 distributions of x at different levels of y under two treatment
 conditions ID=1 (upper panel )  ID=0 (lower panel of the plot).

 bwplot(HF.ELECYEAR ~ stparvotech | ID ,
  data=data, aspect=1,
  layout=c(1,2),
 xlab=Changes in Party Vote Shares,
  xlim=(-20:20),
 ylab=Periods Following Last Federal Election,
 main=Divided Government,
  panel = function(x,y)
   {
   panel.bwplot(x,y)
   panel.abline(v=0, col=red)
   }
  )

 How can I:
 1. Control the font size of the main title, the panel titles, the
 axis, and the tick marks? The usual cex.main=1/3, cex.xlab etc. do
 not work. 

Why do you expect them to? In help(bwplot), the entry for 'main' says:

main: character string or expression or list describing main title
  to be placed on top of each page. Defaults to 'NULL'. Can be
  a character string or expression, or a list with components
  'label, cex, col, font'. The 'label' tag can be omitted if it
  is the first element of the list. Expressions are treated as
  specification of LaTeX-like markup as in 'plotmath'

which suggests that 

main = list(label = Divided Government, cex = 1/3)

is what you want. 


 2. Change the color of the boxes where the panel titles 
 (ID=1  ID=0) are located?

See help(strip.default). You need to use something like

strip = function(..., bg) strip.default(..., bg = white)

(The next version of lattice will have a slightly easier way to do this, 
along with an example.)

All this can also be done by changing the settings globally (see 
help(trellis.par.get)), but that's probably not what you want.

Hope that helps,

Deepayan

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


[R] adding observations to lm for fast recursive residuals?

2004-09-15 Thread ivo_welch-rstat8783
dear R community:  i have been looking but failed to find the 
following:  is there a function in R that updates a plain OLS lm() 
model with one additional observation, so that I can write a function 
that computes recursive residuals *quickly*?

PS: (I looked at package strucchange, but if I am not mistaken, the 
recresid function there takes longer than iterating over the models 
fresh from start to end.)  I know the two functions do not do the same 
thing, but the main part (OLS) is the same:
   handrecurse.test - function( y, x ) { z- rep(NA, T); for (i in 
2:T)  { z[i] - coef(lm(y[1:i] ~ x[1:i]))[2]; }; return(z); }
  system.time(handrecurse.test(y,x))
   [1] 0.69 0.00 0.70 0.00 0.00
  system.time(length(recresid( y~x )))
[1] 1.44 0.07 1.59 0.00 0.00

pointers appreciated.  regards, /iaw
---
ivo welch
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] adding observations to lm for fast recursive residuals?

2004-09-15 Thread roger koenker
In my quantreg package there is a function called lm.fit.recursive() 
that, as the .Rd file
says:

Description:
 This function fits a linear model by recursive least squares.  It
 is a utility routine for the 'khmaladzize' function of the
 quantile regression package.
Usage:
 lm.fit.recursive(X, y, int=TRUE)
Arguments:
   X: Design Matrix
   y: Response Variable
 int: if TRUE then append intercept to X
Value:
 return p by n matrix of fitted parameters, where p. The ith column
 gives the solution up to time i.
It is written in fortran so it should be reasonably quick.
HTH
url:www.econ.uiuc.edu/~rogerRoger Koenker
email   [EMAIL PROTECTED]   Department of Economics
vox:217-333-4558University of Illinois
fax:217-244-6678Champaign, IL 61820
On Sep 15, 2004, at 9:53 AM, [EMAIL PROTECTED] wrote:
dear R community:  i have been looking but failed to find the 
following:  is there a function in R that updates a plain OLS lm() 
model with one additional observation, so that I can write a function 
that computes recursive residuals *quickly*?

PS: (I looked at package strucchange, but if I am not mistaken, the 
recresid function there takes longer than iterating over the models 
fresh from start to end.)  I know the two functions do not do the same 
thing, but the main part (OLS) is the same:
   handrecurse.test - function( y, x ) { z- rep(NA, T); for (i in 
2:T)  { z[i] - coef(lm(y[1:i] ~ x[1:i]))[2]; }; return(z); }
  system.time(handrecurse.test(y,x))
   [1] 0.69 0.00 0.70 0.00 0.00
  system.time(length(recresid( y~x )))
[1] 1.44 0.07 1.59 0.00 0.00

pointers appreciated.  regards, /iaw
---
ivo welch
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! 
http://www.R-project.org/posting-guide.html
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] rolling n-step predictions for ARIMA models

2004-09-15 Thread Paul Gilbert
There are functions in the dse bundle that do this. (See 
featherForecasts and horizonForecasts.) You might look through the 
users' guide to get an idea if they are exactly what you want.

Paul Gilbert
Michael Roberts wrote:
Hello:
I would like to generate rolling, multiperiod forecasts from an 
estimated ARIMA model, but the function predict.Arima seems 
only to generate forecasts from the last observation in the data 
set.  To implement this, I was looking for an argument like 
'newdata=' in predict.lm.  

I can write some code that does this for my particular problem,
but might there exist a package/function that does this that I 
cannot find?

Thanks,
-Michael

Michael J. Roberts
Resource Economics Division
Production, Management, and Technology
USDA-ERS
(202) 694-5557 (phone)
(202) 694-5775 (fax)
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
 

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


[R] [R-pkgs] RODBC 1.1-1

2004-09-15 Thread Prof Brian Ripley
The first non-maintenance update of RODBC since January 2003 is now on 
CRAN and will soon propagate to mirrors.  From the ChangeLog:

* Select the decimal point from Sys.localeconv.

* Add an external reference and finalizer so open channels get
  closed at the end of the session or when there is no R object
  referring to them.

* There is no longer a restriction to 16 channels.

* Add NAMESPACE.

* odbcConnect{Access,Dbase,Excel} allow a missing file name
  (and will bring up a dialog box to search for it).

* odbcGetInfo returns more information in a 8-element character
  vector (based on an idea of Matthew Dowle).

* The C code calls SQLExecuteDirect rather than SQLExecute and
  does not call SQLCloseCursor, based on a problem report from
  Matthew Dowle using MS SQLServer.  Repeated calls to
  sqlGetResults now work.

* New function sqlFetchMore.

* Table names in Access with embedded spaces are mapped to the
  [name space] form which Access requires.

* Table creation no longer removes _ from column names.

* New functions get/setSqlTypeInfo and the typeInfo argument to
  sqlSave allow users to specify the mapping from R types to DBMS
  datatypes.  sqlSave also allows the specification of DBMS
  datatypes by column.

* It is now possible to write more than 255 chars to a field with
  sqlSave and sqlUpdate.

* Dates and timestamps are now read as 'Date' and 'POSIXct'
  columns by sqlGetResults (unless as.is = TRUE for the column).

I have been able to test SQL Server reasonably extensively this time 
around, as well as MySQL, PostgreSQL, Access, Excel and SQLite.

-- 
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-packages mailing list
[EMAIL PROTECTED]
https://stat.ethz.ch/mailman/listinfo/r-packages

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


Re: [R] cluster analysis and null hypothesis testing

2004-09-15 Thread Christian Hennig
Hi,

testing the randomness of a cluster analysis is not a well defined
problem, because it depends crucially on your null model. In fpc, there is
nothing like this. Function prabtest in package prabclus performs such a
test, but this is for a particular data structure, namely presence-absence
data in biogeography. 

In principle, a Monte Carlo test can be constructed (and thus implemented in
R) as follows:

1) You need a null model H_0, from which you generate data.
2) You need a test statistic T.
3) Compute T on your data (call it T_0).
4) Repeat k times:
 a) Generate data from H_0
 b) Compute T on the generated data.
5) The p-value is (K+1)/(k+1), where K is the number of generated datasets
   for which T=T_0 (given that T small indicates the tendency of
   clustering). 

Standard choices for H_0 will be a normal or uniform distribution. (In
prabtest, it is a complicated distribution on presence-absence data.)
There are lots of possible choices of T. prabtest uses the ratio between  
the 25% smallest distances in the dataset and the 25% largest distances.
This should be reasonable in fairly general settings. For a discussion of
this and alternative choices (and references on them), you may take a look
into 

C. Hennig and B. Hausdorf:  Distance-based parametric bootstrap tests for
clustering of species ranges,  Computational
Statistics and Data Analysis 45 (2004), 875-896.

A preprint of this can be obtained from my web page.

If you want to test the significance of a solution from a particular cluster
analysis method, you should think about choosing T so that it is somehow
connected to the method. (In the Hennig and Hausdorf paper, there are for
example two alternatives discussed that are connected to Single Linkage.)

Best,
Christian 

On Wed, 15 Sep 2004, Patrick Giraudoux wrote:

 Hi,
 
 I am wondering if a Monte Carlo method (or equivalent) exist permitting to test the 
 randomness of a cluster analysis (eg got by
 hclust(). I went through the package fpc (maybe too superficially) but dit not 
 find such method.
 
 Thanks for any hint,
 
 Patrick Giraudoux
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
 

***
Christian Hennig
Fachbereich Mathematik-SPST/ZMS, Universitaet Hamburg
[EMAIL PROTECTED], http://www.math.uni-hamburg.de/home/hennig/
###
ich empfehle www.boag-online.de

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


Re: [R] fractal dimension of an image

2004-09-15 Thread Spencer Graves
 From www.r-project.org - search - R site search - fractal 
dimension, I got 17 hits, the third of which discussed package 
RandomFields, which includes a function fractal.  Are you familiar 
with this? 

 hope this helps. 
 spencer graves
p.s.  Have you read the posting guide! 
http://www.R-project.org/posting-guide.html;? 

Rajarshi Guha wrote:
Hello, I have a problem that I think can be solved in R but I'm not sure
how to tie things together.
I have a digital image of a crystal growth in 2 dimensions. And my aim
is to calculate the fractal dimension of the crystal. I was planning to
use the box counting method.
So I need to read in the image in R (for which I can use the pixmap or
rimage packages) and then draw a grid over at a series of resolutions
and for each grid, count how many cells of the grid are occupied by the
crystal. Since the image can be converted to grayscale, I figured that
specifying an intensity threshold would allow me to differentiate
between crystal and surrounding solution.
I know this can be done in Matlab, but since I'm more comfortable in R
can it be done in R? I found the fdim package but that calculates the
dimension for a data.frame - but I'm not sure whether this would be
feasible for my image.
Any pointers would be appreciated
---
Rajarshi Guha [EMAIL PROTECTED] http://jijo.cjb.net
GPG Fingerprint: 0CCA 8EE2 2EEB 25E2 AB04 06F7 1BB9 E634 9B87 56EE
---
So the Zen master asked the hot-dog vendor, 
Can you make me one with everything?
- TauZero on Slashdot

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

--
Spencer Graves, PhD, Senior Development Engineer
O:  (408)938-4420;  mobile:  (408)655-4567
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] heatmap2

2004-09-15 Thread Paul Lepp
Can anybody tell me where to find a copy of heatmap2?  I've seen it in my
travels across the web but didn't bookmark it and can't find it again.
Thanks in advance.

Paul

 `-:-.   ,-;`-:-.   ,-;`-:-.   ,-;`-:-.   ,-;`-:-.   ,-;`-:-.
   `=`,'=/ `=`,'=/ `=`,'=/ `=`,'=/ `=`,'=/ `=`
 ==/==/==/==/==/
   ,=,-=`.,=,-=`.,=,-=`.,=,-=`.,=,-=`.,=,
,-'-'   `-=_,-'-'   `-=_,-'-'   `-=_,-'-'   `-=_,-'-'   `-=_,-'-'
Paul Lepp, Ph.D.   Stanford School of Medicine

VAPAHCS, 154T   Dept. of Microbiology  Immunology
3801 Miranda Ave   Stanford University
Palo Alto, CA 94304   Stanford, CA
(650) 493-5000 x66762  fax: (650) 852-3291
http://cmgm.stanford.edu/~pwlepp  [EMAIL PROTECTED]

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


Re: [R] heatmap2

2004-09-15 Thread Sean Davis
Paul,
It is called heatmap.2 and it is in the gregmisc package.
Sean
On Sep 15, 2004, at 1:03 PM, Paul Lepp wrote:
Can anybody tell me where to find a copy of heatmap2?  I've seen it in 
my
travels across the web but didn't bookmark it and can't find it again.
Thanks in advance.

Paul
 `-:-.   ,-;`-:-.   ,-;`-:-.   ,-;`-:-.   ,-;`-:-.   ,-;`-:-.
   `=`,'=/ `=`,'=/ `=`,'=/ `=`,'=/ `=`,'=/ `=`
==/==/==/==/==/
   ,=,-=`.,=,-=`.,=,-=`.,=,-=`.,=,-=`.,=,
,-'-'   `-=_,-'-'   `-=_,-'-'   `-=_,-'-'   `-=_,-'-'   `-=_,-'-'
Paul Lepp, Ph.D.   Stanford School of Medicine
VAPAHCS, 154T   Dept. of Microbiology  Immunology
3801 Miranda Ave   Stanford University
Palo Alto, CA 94304   Stanford, CA
(650) 493-5000 x66762  fax: (650) 852-3291
http://cmgm.stanford.edu/~pwlepp  [EMAIL PROTECTED]
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! 
http://www.R-project.org/posting-guide.html
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] lomb periodogram package

2004-09-15 Thread Richard Bonneau
Hi,
Does anyone know the name of the package that
includes a function for computing the lomb periodogram on irregular
spaced ts data? I saw the package once ~ 1 month ago but cannot
find it now ...
,
Rich
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Problem installing source packages on OS X

2004-09-15 Thread Aric Gregson
I am attempting to install the Hmisc, rreport and Design packages, but
am not able to do so. I am running R v1.9.1 on Mac OS 10.3.5. 

I have tried installation using both the GUI version of R and also
running R from sudo in a terminal session. In the terminal I receive the
following error:

* Installing *source* package 'Design' ...
** libs
g77   -fno-common  -g -O2 -c lrmfit.f -o lrmfit.o
make: g77: Command not found
make: *** [lrmfit.o] Error 127
ERROR: compilation failed for package 'Design'
** Removing
'/Library/Frameworks/R.framework/Versions/1.9.1/Resources/library/Design
'

I get the same error for Hmisc (rreport is not on CRAN). It looks like
it is trying to use g77 to compile the source package. How can I change
the default compiler? Will this solve the problem? I cannot find a
binary version of either package.

thanks,

aric

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


Re: [R] Problem installing source packages on OS X

2004-09-15 Thread Jari Oksanen
On 15 Sep 2004, at 20:29, Aric Gregson wrote:
I am attempting to install the Hmisc, rreport and Design packages, but
am not able to do so. I am running R v1.9.1 on Mac OS 10.3.5.
I get the same error for Hmisc (rreport is not on CRAN). It looks like
it is trying to use g77 to compile the source package. How can I change
the default compiler? Will this solve the problem? I cannot find a
binary version of either package.
R is trying to build a Fortran program, and it needs a Fortran 
compiler. Fortran compiler does not ship with MacOS X, but you got to 
get one. See the MacOS FAQ for R. If I remember correctly, it tells you 
to go http://hpc.sourceforge.net/ for the compiler.

Normally I wouldn't remember addresses like this, but just today I had 
to make a visit there: I had installed g77 using fink, and that puts 
its stuff into /sw instead of /usr/local. Some R routines had hardcoded 
the g77 path to /usr/local/bin/g77 and so building a package failed in 
the false claim of missing g77 (yeah, it was in the path).

cheers, jari oksanen
--
Jari Oksanen, Oulu, Finland
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] BUGS and OS X

2004-09-15 Thread Tamas K Papp
Is there a way of running BUGS on OS X (from R)?  I only see Windows
versions on their website.

If not, what are the alternatives for Bayesian analysis?

Thanks,

Tamas

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


Re: [R] Bessel function

2004-09-15 Thread cstrato
You can find C++ source code for Bessel function and similar
functions for example here:
http://root.cern.ch/root/htmldoc/src/TMath.cxx.html
Hope this helps
Best regards
Christian
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
C.h.r.i.s.t.i.a.n. .S.t.r.a.t.o.w.a
V.i.e.n.n.a. .A.u.s.t.r.i.a
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
Ravi Varadhan wrote:
Dear Beat,
 
You might also want to look at the book by Zang and Jin -  Computation of Special Functions, John Wiley.  They have Fortran sources for all the special functions covered in there.
 
Ravi.
 
Ravi Varadhan, Ph.D.
Assistant Professor,  The Center on Aging and Health
Division of Geriatric Medicine and Gerontology,
Johns Hopkins University,
2024 E. Monument Street, Suite 2-700
Baltimore, MD 21205.
(410) 502 - 7806.


From: [EMAIL PROTECTED] on behalf of [EMAIL PROTECTED]
Sent: Wed 9/15/2004 4:24 AM
To: [EMAIL PROTECTED]
Subject: [R] Bessel function


 Dear all
Currently, I'm implementing the generalized hyperbolic distribution into
Splus. Unfortunately the Bessel function is not implemented in Splus. In
R the Bessel function does exist but it is an internal function and I'm
not able to look at the code.
Is there any possibility to see the code of the Bessel function in R or
does anybody has an implementation of the Bessel function in Splus?
Thanks a lot for your help.
Beat
[[alternative HTML version deleted]]
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[[alternative HTML version deleted]]
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

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


[R] Slightly off-topic --- distribution name.

2004-09-15 Thread Rolf Turner

I've built R functions to ``effect'' a particular distribution, and
would like to find out if that distribution is already ``known'' by
an existing name.  (I.e. suppose it were called the ``Melvin''
distribution --- I've built dmelvin, pmelvin, qmelvin, and rmelvin as
it were, but I need a real name to substitute for melvin.)

The distribution is really just a toy --- but it provides a nice (and
``non-obviouse'') example of a two parameter distribution where both
the moment and maximum likelihood equations for the parameter
estimators are readily solvable, but at the same time are
``interesting''.  So it's good for exercises in an intro math-stats
course.

The distribution is simply that of the ***difference*** of two
independent exponential variates, with different parameters.

I.e.  X = U - V  where U ~ exp(beta) and V ~ exp(alpha) (where
E(U) = beta, E(V) = alpha).

This makes the distribution of X something like an asymetric Laplace
distribution, with its mode at 0.  (One could shift the mode too, but
that would add a third parameter, which would be de trop.)

Anyhow:  Is this a ``known'' distribution?  Does it have a name?
(I've never seen it mentioned in any of the intro math-stat books
that I've looked into.) If not, can anyone suggest a good name for
it?  (Don't be rude now!)

cheers,

Rolf Turner
[EMAIL PROTECTED]

P. S.  To save you putting pen to paper and working it out,
   the density function is

   { exp(x/alpha)/(alpha + beta) for x = 0
f(x) = {
   { exp(-x/beta)/(alpha + beta) for x = 0

   The mean and variance are mu = beta - alpha and
   sigma^2 = alpha^2 + beta^2 respectfully. :-)

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


Re: [R] BUGS and OS X

2004-09-15 Thread Tamas K Papp
On Wed, Sep 15, 2004 at 02:21:18PM -0400, Liaw, Andy wrote:
 That's more of a question for the BUGS developers.  BUGS is not open source,
 so whatever binary is provided, that's all you can use.  If I'm not
 mistaken, WinBUGS is the only version under development.

I found something called JAGS, and I am still exploring it.  It
appears to be an open-source BUGS replacement, thought with
limitations.

I was asking what software people would recommend for the same
functionality, not a drop-in replacement.  I am just baffled by the
bewildering array of R packages, and would be so happy if somebody
told me what THEY use for Bayesian analysis, so I could read the docs
and get started.  MCMC? Boa? etc.  Suggestions on how experienced
users do bayesian analysis in R would be welcome.

Thanks,

Tamas

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


Re: [R] Slightly off-topic --- distribution name.

2004-09-15 Thread A.J. Rossini

Have you checked Johnson and Kotz?  That's the obvious place to start
looking for distributions beyond the usual.

Rolf Turner [EMAIL PROTECTED] writes:

 I've built R functions to ``effect'' a particular distribution, and
 would like to find out if that distribution is already ``known'' by
 an existing name.  (I.e. suppose it were called the ``Melvin''
 distribution --- I've built dmelvin, pmelvin, qmelvin, and rmelvin as
 it were, but I need a real name to substitute for melvin.)

 The distribution is really just a toy --- but it provides a nice (and
 ``non-obviouse'') example of a two parameter distribution where both
 the moment and maximum likelihood equations for the parameter
 estimators are readily solvable, but at the same time are
 ``interesting''.  So it's good for exercises in an intro math-stats
 course.

 The distribution is simply that of the ***difference*** of two
 independent exponential variates, with different parameters.

 I.e.  X = U - V  where U ~ exp(beta) and V ~ exp(alpha) (where
 E(U) = beta, E(V) = alpha).

 This makes the distribution of X something like an asymetric Laplace
 distribution, with its mode at 0.  (One could shift the mode too, but
 that would add a third parameter, which would be de trop.)

 Anyhow:  Is this a ``known'' distribution?  Does it have a name?
 (I've never seen it mentioned in any of the intro math-stat books
 that I've looked into.) If not, can anyone suggest a good name for
 it?  (Don't be rude now!)

   cheers,

   Rolf Turner
   [EMAIL PROTECTED]

 P. S.  To save you putting pen to paper and working it out,
the density function is

{ exp(x/alpha)/(alpha + beta) for x = 0
   f(x) = {
{ exp(-x/beta)/(alpha + beta) for x = 0

The mean and variance are mu = beta - alpha and
sigma^2 = alpha^2 + beta^2 respectfully. :-)

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


-- 
Anthony Rossini Research Associate Professor
[EMAIL PROTECTED]http://www.analytics.washington.edu/ 
Biomedical and Health Informatics   University of Washington
Biostatistics, SCHARP/HVTN  Fred Hutchinson Cancer Research Center
UW (Tu/Th/F): 206-616-7630 FAX=206-543-3461 | Voicemail is unreliable
FHCRC  (M/W): 206-667-7025 FAX=206-667-4812 | use Email

CONFIDENTIALITY NOTICE: This e-mail message and any attachme...{{dropped}}

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


[R] replacing NA's with 0 in a dataframe for specified columns

2004-09-15 Thread David Kane
I know that there must be a cool way of doing this, but I can't think
of it. Let's say I have an dataframe with NA's.

 x - data.frame(a = c(0,1,2,NA), b = c(0,NA,1,2), c = c(NA, 0, 1, 2))
 x
   a  b  c
1  0  0 NA
2  1 NA  0
3  2  1  1
4 NA  2  2
 

I know it is easy to replace all the NA's with zeroes.

 x[is.na(x)] - 0
 x
  a b c
1 0 0 0
2 1 0 0
3 2 1 1
4 0 2 2
 

But how do I do this for just columns a and c, leaving the NA in
column b alone?

Thanks,

Dave Kane

 R.version
 _
platform i686-pc-linux-gnu
arch i686 
os   linux-gnu
system   i686, linux-gnu  
status
major1
minor9.1  
year 2004 
month06   
day  21   
language R


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


Re: [R] lomb periodogram package

2004-09-15 Thread Earl F. Glynn
 Does anyone know the name of the package that
 includes a function for computing the lomb periodogram on irregular
 spaced ts data? I saw the package once ~ 1 month ago but cannot
 find it now ...

I have a LombScargleLibrary.R file that I will be talking about at CAMDA '04
in November:

Searching for Periodic Microarray Expression Patterns Using Lomb-Scargle
Periodograms
http://research.stowers-institute.org/efg/2004/CAMDA/

I'll get the file online by the time of the conference on Nov 10:
http://www.camda.duke.edu/camda04

I could E-mail you the preliminary version before then if you'd like.  Just
drop me an E-mail.

efg
Earl F. Glynn
Scientific Programmer
Bioinformatics Department
Stowers Institute for Medical Research

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


[R] [R-pkgs] Announcing snowFT 0.1

2004-09-15 Thread A.J. Rossini

Parallel programming with snowFT

Our package snowFT is now available at CRAN. It is an extention of the
package snow, which adds fault tolerance (in the sense of recomputing
computational units when hardware/network failures occur on compute
nodes) and a tighter notion of reproducibility for computations
running on clusters.  It additionally provides tools for flexible
management of cluster size as well as computation transparency.

snowFT is written by Hana Sevcikova and Tony Rossini.

(this tool currently requires rpvm for the SNOW backend, though we are
exploring the possibility of extensions for the Rmpi backend.  It is
unlikely that these extensions will be implemented for the socket
backend.

best,
-tony



-- 
Anthony Rossini Research Associate Professor
[EMAIL PROTECTED]http://www.analytics.washington.edu/ 
Biomedical and Health Informatics   University of Washington
Biostatistics, SCHARP/HVTN  Fred Hutchinson Cancer Research Center
UW (Tu/Th/F): 206-616-7630 FAX=206-543-3461 | Voicemail is unreliable
FHCRC  (M/W): 206-667-7025 FAX=206-667-4812 | use Email

CONFIDENTIALITY NOTICE: This e-mail message and any attachme...{{dropped}}

___
R-packages mailing list
[EMAIL PROTECTED]
https://stat.ethz.ch/mailman/listinfo/r-packages

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


Re: [R] replacing NA's with 0 in a dataframe for specified columns

2004-09-15 Thread Spencer Graves
Have you considered the following: 

 x - data.frame(a = c(0,1,2,NA), b = c(0,NA,1,2), c = c(NA, 0, 1, 2))
 x$a[is.na(x$a)] - 0
 x$c[is.na(x$c)] - 0
 x
 a  b c
1 0  0 0
2 1 NA 0
3 2  1 1
4 0  2 2
 hope this helps.  spencer graves
David Kane wrote:
I know that there must be a cool way of doing this, but I can't think
of it. Let's say I have an dataframe with NA's.
 

x - data.frame(a = c(0,1,2,NA), b = c(0,NA,1,2), c = c(NA, 0, 1, 2))
x
   

  a  b  c
1  0  0 NA
2  1 NA  0
3  2  1  1
4 NA  2  2
 

I know it is easy to replace all the NA's with zeroes.
 

x[is.na(x)] - 0
x
   

 a b c
1 0 0 0
2 1 0 0
3 2 1 1
4 0 2 2
 

But how do I do this for just columns a and c, leaving the NA in
column b alone?
Thanks,
Dave Kane
 

R.version
   

_
platform i686-pc-linux-gnu
arch i686 
os   linux-gnu
system   i686, linux-gnu  
status
major1
minor9.1  
year 2004 
month06   
day  21   
language R
 

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

--
Spencer Graves, PhD, Senior Development Engineer
O:  (408)938-4420;  mobile:  (408)655-4567
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] replacing NA's with 0 in a dataframe for specified columns

2004-09-15 Thread Corey Moffet
try:

x[is.na(x$a) | is.na(x$c),] - 0

At 02:44 PM 9/15/2004 -0400, David Kane wrote:
I know that there must be a cool way of doing this, but I can't think
of it. Let's say I have an dataframe with NA's.

 x - data.frame(a = c(0,1,2,NA), b = c(0,NA,1,2), c = c(NA, 0, 1, 2))
 x
   a  b  c
1  0  0 NA
2  1 NA  0
3  2  1  1
4 NA  2  2
 

I know it is easy to replace all the NA's with zeroes.

 x[is.na(x)] - 0
 x
  a b c
1 0 0 0
2 1 0 0
3 2 1 1
4 0 2 2
 

But how do I do this for just columns a and c, leaving the NA in
column b alone?

Thanks,

Dave Kane

 R.version
 _
platform i686-pc-linux-gnu
arch i686 
os   linux-gnu
system   i686, linux-gnu  
status
major1
minor9.1  
year 2004 
month06   
day  21   
language R


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

With best wishes and kind regards I am

Sincerely,

Corey A. Moffet
Rangeland Scientist

##
 
USDA-ARS#
Northwest Watershed Research Center #
800 Park Blvd, Plaza IV, Suite 105  ###   
Boise, ID 83712-7716   ##  # #
Voice: (208) 422-0718  ##     
FAX:   (208) 334-1502  ## #   #   
   ###
##

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


RE: [R] replacing NA's with 0 in a dataframe for specified columns

2004-09-15 Thread Berton Gunter

But Spencer's solution would require looping to generalize to a large
numbers of columns.

In fact, you've given the answer already in your post:

xsub-x[,somecols]
xsub[is.na(xsub)]-0
x[,somecols]-xsub

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

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of Spencer Graves
 Sent: Wednesday, September 15, 2004 12:01 PM
 To: David Kane
 Cc: [EMAIL PROTECTED]
 Subject: Re: [R] replacing NA's with 0 in a dataframe for 
 specified columns
 
 Have you considered the following: 
 
   x - data.frame(a = c(0,1,2,NA), b = c(0,NA,1,2), c = 
 c(NA, 0, 1, 2))
   x$a[is.na(x$a)] - 0
   x$c[is.na(x$c)] - 0
   x
   a  b c
 1 0  0 0
 2 1 NA 0
 3 2  1 1
 4 0  2 2
 
   hope this helps.  spencer graves
 
 David Kane wrote:
 
 I know that there must be a cool way of doing this, but I can't think
 of it. Let's say I have an dataframe with NA's.
 
   
 
 x - data.frame(a = c(0,1,2,NA), b = c(0,NA,1,2), c = c(NA, 
 0, 1, 2))
 x
 
 
a  b  c
 1  0  0 NA
 2  1 NA  0
 3  2  1  1
 4 NA  2  2
   
 
 
 I know it is easy to replace all the NA's with zeroes.
 
   
 
 x[is.na(x)] - 0
 x
 
 
   a b c
 1 0 0 0
 2 1 0 0
 3 2 1 1
 4 0 2 2
   
 
 
 But how do I do this for just columns a and c, leaving the NA in
 column b alone?
 
 Thanks,
 
 Dave Kane
 
   
 
 R.version
 
 
  _
 platform i686-pc-linux-gnu
 arch i686 
 os   linux-gnu
 system   i686, linux-gnu  
 status
 major1
 minor9.1  
 year 2004 
 month06   
 day  21   
 language R
   
 
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
   
 
 
 -- 
 Spencer Graves, PhD, Senior Development Engineer
 O:  (408)938-4420;  mobile:  (408)655-4567
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

Re: [R] replacing NA's with 0 in a dataframe for specified columns

2004-09-15 Thread John Fox
Dear David,

How about the following?

cols - c(1,3)
x[,cols][is.na(x[,cols])] - 0

I hope that this helps,
 John

On Wed, 15 Sep 2004 14:44:53 -0400
 David Kane [EMAIL PROTECTED] wrote:
 I know that there must be a cool way of doing this, but I can't think
 of it. Let's say I have an dataframe with NA's.
 
  x - data.frame(a = c(0,1,2,NA), b = c(0,NA,1,2), c = c(NA, 0, 1,
 2))
  x
a  b  c
 1  0  0 NA
 2  1 NA  0
 3  2  1  1
 4 NA  2  2
  
 
 I know it is easy to replace all the NA's with zeroes.
 
  x[is.na(x)] - 0
  x
   a b c
 1 0 0 0
 2 1 0 0
 3 2 1 1
 4 0 2 2
  
 
 But how do I do this for just columns a and c, leaving the NA in
 column b alone?
 
 Thanks,
 
 Dave Kane
 
  R.version
  _
 platform i686-pc-linux-gnu
 arch i686 
 os   linux-gnu
 system   i686, linux-gnu  
 status
 major1
 minor9.1  
 year 2004 
 month06   
 day  21   
 language R
 
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide!
 http://www.R-project.org/posting-guide.html


John Fox
Department of Sociology
McMaster University
Hamilton, Ontario, Canada
http://socserv.mcmaster.ca/jfox/

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


Re: [R] replacing NA's with 0 in a dataframe for specified columns

2004-09-15 Thread Chuck Cleland
mydata - data.frame(a = c(0,1,2,NA), b = c(0,NA,1,2), c = c(NA, 0, 1, 2))
mydata
   a  b  c
1  0  0 NA
2  1 NA  0
3  2  1  1
4 NA  2  2
mydata[,c(a, c)] -
apply(mydata[,c(a,c)], 2, function(x){replace(x, is.na(x), 0)})
mydata
  a  b c
1 0  0 0
2 1 NA 0
3 2  1 1
4 0  2 2
David Kane wrote:
I know that there must be a cool way of doing this, but I can't think
of it. Let's say I have an dataframe with NA's.

x - data.frame(a = c(0,1,2,NA), b = c(0,NA,1,2), c = c(NA, 0, 1, 2))
x
   a  b  c
1  0  0 NA
2  1 NA  0
3  2  1  1
4 NA  2  2
I know it is easy to replace all the NA's with zeroes.

x[is.na(x)] - 0
x
  a b c
1 0 0 0
2 1 0 0
3 2 1 1
4 0 2 2
But how do I do this for just columns a and c, leaving the NA in
column b alone?
Thanks,
Dave Kane

R.version
 _
platform i686-pc-linux-gnu
arch i686 
os   linux-gnu
system   i686, linux-gnu  
status
major1
minor9.1  
year 2004 
month06   
day  21   
language R

__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
--
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 452-1424 (M, W, F)
fax: (917) 438-0894
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] replacing NA's with 0 in a dataframe for specified columns

2004-09-15 Thread Gabor Grothendieck

The col funtion can be helpful here.  We want to satisfy two conditions: 

1. the element is an NA
2. the element lies in one of the specified columns

The first two lines below calculate logical vectors for these two,
respectively, and the last line assigns 0 to those elements.

isna - is.na(x)
iscol - col(isna) %in% c(1,3)
x[isna  iscol] - 0

To specify the columns by name rather than number first
calculate their column numbers and then proceed as before.
That is, replace the second line by:

cols - match(c(a, c), colnames(x))
iscol - col(isna) %in% cols





David Kane dave at kanecap.com writes:

: 
: I know that there must be a cool way of doing this, but I can't think
: of it. Let's say I have an dataframe with NA's.
: 
:  x - data.frame(a = c(0,1,2,NA), b = c(0,NA,1,2), c = c(NA, 0, 1, 2))
:  x
:a  b  c
: 1  0  0 NA
: 2  1 NA  0
: 3  2  1  1
: 4 NA  2  2
:  
: 
: I know it is easy to replace all the NA's with zeroes.
: 
:  x[is.na(x)] - 0
:  x
:   a b c
: 1 0 0 0
: 2 1 0 0
: 3 2 1 1
: 4 0 2 2
:  
: 
: But how do I do this for just columns a and c, leaving the NA in
: column b alone?
: 
: Thanks,
: 
: Dave Kane
: 
:  R.version
:  _
: platform i686-pc-linux-gnu
: arch i686 
: os   linux-gnu
: system   i686, linux-gnu  
: status
: major1
: minor9.1  
: year 2004 
: month06   
: day  21   
: language R
: 
: 
: __
: R-help at stat.math.ethz.ch mailing list
: https://stat.ethz.ch/mailman/listinfo/r-help
: PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
: 
:

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


RE: [R] Splitting vector into individual elements

2004-09-15 Thread Liaw, Andy
do.call() is good for this, I believe:

 offred.rgb - c(1, 0, 0) * 0.60
 offred.col - do.call(rgb, c(as.list(offred.rgb), names=offred))
 offred.col
[1] #99

HTH,
Andy

 From: Paul Roebuck
 
 Is there a means to split a vector into its individual
 elements without going the brute-force route for arguments
 to a predefined function call?
 
 offred.rgb - c(1, 0, 0) * 0.60;
 
 ## Brute force style
 offred.col - rgb(offred.rgb[1],
   offred.rgb[2],
   offred.rgb[3],
   names = offred)
 ## Desired style
 offred.col - rgb(silver.bullet(offred.rgb),
   names = offred)
 
 Neither of my attempts gets it right.
 
 silver.bullet.try1 - function(x) {
 expr - cat(x, sep = ,)
 return(parse(text = expr))
 }
 
 silver.bullet.try2 - function(x) {
 expr - expression(cat(x, sep = ,))
 return(eval(expr))
 }
 
 --
 SIGSIG -- signature too long (core dumped)
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 


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


Re: [R] Splitting vector into individual elements

2004-09-15 Thread Spencer Graves
 Have you considered do.call: 

 do.call(rgb, as.list((1:3)/10))
[1] #1A334C
same as: 
 rgb(.1, .2, .3)
[1] #1A334C

 Hope this helps.  spencer graves
Paul Roebuck wrote:
Is there a means to split a vector into its individual
elements without going the brute-force route for arguments
to a predefined function call?
   offred.rgb - c(1, 0, 0) * 0.60;
   ## Brute force style
   offred.col - rgb(offred.rgb[1],
 offred.rgb[2],
 offred.rgb[3],
 names = offred)
   ## Desired style
   offred.col - rgb(silver.bullet(offred.rgb),
 names = offred)
Neither of my attempts gets it right.
   silver.bullet.try1 - function(x) {
   expr - cat(x, sep = ,)
   return(parse(text = expr))
   }
   silver.bullet.try2 - function(x) {
   expr - expression(cat(x, sep = ,))
   return(eval(expr))
   }
--
SIGSIG -- signature too long (core dumped)
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
 

--
Spencer Graves, PhD, Senior Development Engineer
O:  (408)938-4420;  mobile:  (408)655-4567
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] RWAVE axis notation

2004-09-15 Thread Jonathan M. Lees
I am using Rwave wavelets and I need better axis notation.
Does anyone have code similar to
matlab's
centfreq
or
scale2frq
functions that turn a scale for a wavelet transform into
a good looking scale to plot on my wavelet transfoms?
I am using the rwave package to investigate
seismic signals from an exploding volcano.
Jonathan Lees
--
==
Prof. Jonathan M. Lees
Department of Geological Sciences
CB #3315, Mitchell Hall
University of North Carolina
Chapel Hill, NC  27599-3315
(919) 962-0695
FAX (919) 966-4519
[EMAIL PROTECTED]
http://www.unc.edu/~leesj
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Splitting vector into individual elements

2004-09-15 Thread Gabor Grothendieck
Paul Roebuck roebuck at odin.mdacc.tmc.edu writes:
 
: Is there a means to split a vector into its individual
: elements without going the brute-force route for arguments
: to a predefined function call?
: 
: offred.rgb - c(1, 0, 0) * 0.60;
: 
: ## Brute force style
: offred.col - rgb(offred.rgb[1],
:   offred.rgb[2],
:   offred.rgb[3],
:   names = offred)
: ## Desired style
: offred.col - rgb(silver.bullet(offred.rgb),
:   names = offred)


See:

http://maths.newcastle.edu.au/~rking/R/help/03a/7417.html

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


Re: [R] Splitting vector into individual elements

2004-09-15 Thread John Fox
Dear Paul,

How about do.call(rgb, as.list(offred.rgb)) ?

I hope that this helps,
 John

On Wed, 15 Sep 2004 15:20:24 -0500 (CDT)
 Paul Roebuck [EMAIL PROTECTED] wrote:
 Is there a means to split a vector into its individual
 elements without going the brute-force route for arguments
 to a predefined function call?
 
 offred.rgb - c(1, 0, 0) * 0.60;
 
 ## Brute force style
 offred.col - rgb(offred.rgb[1],
   offred.rgb[2],
   offred.rgb[3],
   names = offred)
 ## Desired style
 offred.col - rgb(silver.bullet(offred.rgb),
   names = offred)
 
 Neither of my attempts gets it right.
 
 silver.bullet.try1 - function(x) {
 expr - cat(x, sep = ,)
 return(parse(text = expr))
 }
 
 silver.bullet.try2 - function(x) {
 expr - expression(cat(x, sep = ,))
 return(eval(expr))
 }
 
 --
 SIGSIG -- signature too long (core dumped)
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide!
 http://www.R-project.org/posting-guide.html


John Fox
Department of Sociology
McMaster University
Hamilton, Ontario, Canada
http://socserv.mcmaster.ca/jfox/

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


Re: [R] Splitting vector into individual elements

2004-09-15 Thread Peter Dalgaard
Paul Roebuck [EMAIL PROTECTED] writes:

 Is there a means to split a vector into its individual
 elements without going the brute-force route for arguments
 to a predefined function call?
 
 offred.rgb - c(1, 0, 0) * 0.60;
 
 ## Brute force style
 offred.col - rgb(offred.rgb[1],
   offred.rgb[2],
   offred.rgb[3],
   names = offred)
 ## Desired style
 offred.col - rgb(silver.bullet(offred.rgb),
   names = offred)

The closest is probably this:

offred.col - do.call(rgb, c(as.list(offred.rgb), 
   list(names=offred)))

(ever read/seen The Handmaid's Tale, btw?)

-- 
   O__   Peter Dalgaard Blegdamsvej 3  
  c/ /'_ --- Dept. of Biostatistics 2200 Cph. N   
 (*) \(*) -- University of Copenhagen   Denmark  Ph: (+45) 35327918
~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907

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


[R] efficient submatrix extraction

2004-09-15 Thread Rajarshi Guha
Hi,
  I have a matrix of say 1024x1024 and I want to look at it in chunks.
That is I'd like to divide into a series of submatrices of order 2x2.

| 1 2 3 4 5 6 7 8 ... |
| 1 2 3 4 5 6 7 8 ... |
| 1 2 3 4 5 6 7 8 ... |
| 1 2 3 4 5 6 7 8 ... |
...

So the first submatrix would be

| 1 2 |
| 1 2 |

the second one would be

| 3 4 |
| 3 4 |

and so on. That is I want the matrix to be evenly divided into 2x2
submatrices. Now I'm also doing this subdivision into 4x4, 8x8 ...
256x256 submatrices.

Currently I'm using loops and I'm sure there is a mroe efficient way to
do it:

m - matrix(runif(1024*1024), nrow=1024)
boxsize - 2^(1:8)

for (b in boxsize) {
bcount - 0
bstart - seq(1,1024, by=b)
for (x in bstart) {
for (y in bstart) {
xend - x + b - 1
yend - y + b - 1
if (length(which( m[ x:xend, y:yend ]  0.7))  0) {
bcount - bcount + 1
}
}
}
}

Is there any way to vectorize the two inner loops?

Thanks,

---
Rajarshi Guha [EMAIL PROTECTED] http://jijo.cjb.net
GPG Fingerprint: 0CCA 8EE2 2EEB 25E2 AB04 06F7 1BB9 E634 9B87 56EE
---
The way to love anything is to realize that it might be lost.

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


Re: [R] Splitting vector into individual elements

2004-09-15 Thread Paul Roebuck
On Wed, 15 Sep 2004, Peter Dalgaard wrote:

 Paul Roebuck [EMAIL PROTECTED] writes:

  Is there a means to split a vector into its individual
  elements without going the brute-force route for arguments
  to a predefined function call?
 
  offred.rgb - c(1, 0, 0) * 0.60;
 
  ## Brute force style
  offred.col - rgb(offred.rgb[1],
offred.rgb[2],
offred.rgb[3],
names = offred)
  ## Desired style
  offred.col - rgb(silver.bullet(offred.rgb),
names = offred)

 The closest is probably this:

 offred.col - do.call(rgb, c(as.list(offred.rgb),
list(names=offred)))


Everyone offered 'do.call' as the solution. While that
works, is it to say that there is no means of expanding
the expression as an argument to the original function?

 (ever read/seen The Handmaid's Tale, btw?)


Not yet. Though renaming my sample variable 'off.red.col'
would avoid future confusion with oppressed handmaids.

--
SIGSIG -- signature too long (core dumped)

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


Re: [R] efficient submatrix extraction

2004-09-15 Thread Tony Plate
I think you should be able to do something with reassigning the dim 
attribute, and then using apply(), something along the lines of the 
following (which doesn't do your computation on the data in the subarrays, 
but merely illustrates how to create and access them):

 x - matrix(1:64,ncol=8)
 x
 [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
[1,]19   17   25   33   41   49   57
[2,]2   10   18   26   34   42   50   58
[3,]3   11   19   27   35   43   51   59
[4,]4   12   20   28   36   44   52   60
[5,]5   13   21   29   37   45   53   61
[6,]6   14   22   30   38   46   54   62
[7,]7   15   23   31   39   47   55   63
[8,]8   16   24   32   40   48   56   64
 x2 - x
 dim(x2) - c(2,4,2,4)
 x2[,1,,1]
 [,1] [,2]
[1,]19
[2,]2   10
 x2[,2,,1]
 [,1] [,2]
[1,]3   11
[2,]4   12
 x2[,1,,2]
 [,1] [,2]
[1,]   17   25
[2,]   18   26
 x4 - x
 dim(x4) - c(4,2,4,2)
 x4[,1,,1]
 [,1] [,2] [,3] [,4]
[1,]19   17   25
[2,]2   10   18   26
[3,]3   11   19   27
[4,]4   12   20   28
 invisible(apply(x4, c(2,4), print))
 [,1] [,2] [,3] [,4]
[1,]19   17   25
[2,]2   10   18   26
[3,]3   11   19   27
[4,]4   12   20   28
 [,1] [,2] [,3] [,4]
[1,]5   13   21   29
[2,]6   14   22   30
[3,]7   15   23   31
[4,]8   16   24   32
 [,1] [,2] [,3] [,4]
[1,]   33   41   49   57
[2,]   34   42   50   58
[3,]   35   43   51   59
[4,]   36   44   52   60
 [,1] [,2] [,3] [,4]
[1,]   37   45   53   61
[2,]   38   46   54   62
[3,]   39   47   55   63
[4,]   40   48   56   64

hope this helps,
Tony Plate
At Wednesday 03:10 PM 9/15/2004, Rajarshi Guha wrote:
Hi,
  I have a matrix of say 1024x1024 and I want to look at it in chunks.
That is I'd like to divide into a series of submatrices of order 2x2.
| 1 2 3 4 5 6 7 8 ... |
| 1 2 3 4 5 6 7 8 ... |
| 1 2 3 4 5 6 7 8 ... |
| 1 2 3 4 5 6 7 8 ... |
...
So the first submatrix would be
| 1 2 |
| 1 2 |
the second one would be
| 3 4 |
| 3 4 |
and so on. That is I want the matrix to be evenly divided into 2x2
submatrices. Now I'm also doing this subdivision into 4x4, 8x8 ...
256x256 submatrices.
Currently I'm using loops and I'm sure there is a mroe efficient way to
do it:
m - matrix(runif(1024*1024), nrow=1024)
boxsize - 2^(1:8)
for (b in boxsize) {
bcount - 0
bstart - seq(1,1024, by=b)
for (x in bstart) {
for (y in bstart) {
xend - x + b - 1
yend - y + b - 1
if (length(which( m[ x:xend, y:yend ]  0.7))  0) {
bcount - bcount + 1
}
}
}
}
Is there any way to vectorize the two inner loops?
Thanks,
---
Rajarshi Guha [EMAIL PROTECTED] http://jijo.cjb.net
GPG Fingerprint: 0CCA 8EE2 2EEB 25E2 AB04 06F7 1BB9 E634 9B87 56EE
---
The way to love anything is to realize that it might be lost.
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


RE: [R] Splitting vector into individual elements

2004-09-15 Thread Liaw, Andy
 From: Paul Roebuck
 
 On Wed, 15 Sep 2004, Peter Dalgaard wrote:
 
  Paul Roebuck [EMAIL PROTECTED] writes:
 
   Is there a means to split a vector into its individual
   elements without going the brute-force route for arguments
   to a predefined function call?
  
   offred.rgb - c(1, 0, 0) * 0.60;
  
   ## Brute force style
   offred.col - rgb(offred.rgb[1],
 offred.rgb[2],
 offred.rgb[3],
 names = offred)
   ## Desired style
   offred.col - rgb(silver.bullet(offred.rgb),
 names = offred)
 
  The closest is probably this:
 
  offred.col - do.call(rgb, c(as.list(offred.rgb),
 list(names=offred)))
 
 Everyone offered 'do.call' as the solution. While that
 works, is it to say that there is no means of expanding
 the expression as an argument to the original function?

What would be the point?  That's what do.call() does for you internally.

Andy
 
  (ever read/seen The Handmaid's Tale, btw?)
 
 
 Not yet. Though renaming my sample variable 'off.red.col'
 would avoid future confusion with oppressed handmaids.
 
 --
 SIGSIG -- signature too long (core dumped)
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 


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


Re: [R] Splitting vector into individual elements

2004-09-15 Thread Peter Dalgaard
Paul Roebuck [EMAIL PROTECTED] writes:

 Everyone offered 'do.call' as the solution. While that
 works, is it to say that there is no means of expanding
 the expression as an argument to the original function?

Not really. You need an explicit expansion of the argument to a list
somehow, and there's no silver bullet that can convert one function
argument to several. There are solutions without do.call, like

 offred.rgb - c(1, 0, 0) * 0.60
 x - quote(rgb(.,.,.,names=offred))
 x[2:4] - as.list(offred.rgb)
 eval(x)
   offred
#99

but you might find it difficult to explain how it works a year later


-- 
   O__   Peter Dalgaard Blegdamsvej 3  
  c/ /'_ --- Dept. of Biostatistics 2200 Cph. N   
 (*) \(*) -- University of Copenhagen   Denmark  Ph: (+45) 35327918
~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907

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


RE: [R] Splitting vector into individual elements

2004-09-15 Thread Liaw, Andy
 From: Liaw, Andy
 
  From: Paul Roebuck
  
  On Wed, 15 Sep 2004, Peter Dalgaard wrote:
  
   Paul Roebuck [EMAIL PROTECTED] writes:
  
Is there a means to split a vector into its individual
elements without going the brute-force route for arguments
to a predefined function call?
   
offred.rgb - c(1, 0, 0) * 0.60;
   
## Brute force style
offred.col - rgb(offred.rgb[1],
  offred.rgb[2],
  offred.rgb[3],
  names = offred)
## Desired style
offred.col - rgb(silver.bullet(offred.rgb),
  names = offred)
  
   The closest is probably this:
  
   offred.col - do.call(rgb, c(as.list(offred.rgb),
  list(names=offred)))
  
  Everyone offered 'do.call' as the solution. While that
  works, is it to say that there is no means of expanding
  the expression as an argument to the original function?
 
 What would be the point?  That's what do.call() does for you 
 internally.

Is this what you're after?

 toCall - c(as.name(rgb), as.list(offred.rgb), names=offred)
 toCall
[[1]]
rgb

[[2]]
[1] 0.6

[[3]]
[1] 0

[[4]]
[1] 0

$names
[1] offred

 toCall - as.call(toCall)
 toCall
rgb(0.6, 0, 0, names = offred)
 eval(toCall)
[1] #99

Andy
 
 Andy
  
   (ever read/seen The Handmaid's Tale, btw?)
  
  
  Not yet. Though renaming my sample variable 'off.red.col'
  would avoid future confusion with oppressed handmaids.
  
  --
  SIGSIG -- signature too long (core dumped)
  
  __
  [EMAIL PROTECTED] mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide! 
  http://www.R-project.org/posting-guide.html
  
 
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 
 
 --
 
 Notice:  This e-mail message, together with any attachments, 
 contains information of Merck  Co., Inc. (One Merck Drive, 
 Whitehouse Station, New Jersey, USA 08889), and/or its 
 affiliates (which may be known outside the United States as 
 Merck Frosst, Merck Sharp  Dohme or MSD and in Japan, as 
 Banyu) that may be confidential, proprietary copyrighted 
 and/or legally privileged. It is intended solely for the use 
 of the individual or entity named on this message.  If you 
 are not the intended recipient, and have received this 
 message in error, please notify us immediately by reply 
 e-mail and then delete it from your system.
 --
 


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


Re: [R] Slightly off-topic --- distribution name.

2004-09-15 Thread David Scott
I believe this is the skew-Laplace distribution, although the skew-Laplace 
does allow for the location of the mode of the distribution to vary.

Have a look at the function dskewlap in HyperbolicDist. The help on that 
function gives a reference to a paper by Feiller et al which describes the 
distribution.

David Scott

On Wed, 15 Sep 2004, Rolf Turner wrote:
I've built R functions to ``effect'' a particular distribution, and
would like to find out if that distribution is already ``known'' by
an existing name.  (I.e. suppose it were called the ``Melvin''
distribution --- I've built dmelvin, pmelvin, qmelvin, and rmelvin as
it were, but I need a real name to substitute for melvin.)
The distribution is really just a toy --- but it provides a nice (and
``non-obviouse'') example of a two parameter distribution where both
the moment and maximum likelihood equations for the parameter
estimators are readily solvable, but at the same time are
``interesting''.  So it's good for exercises in an intro math-stats
course.
The distribution is simply that of the ***difference*** of two
independent exponential variates, with different parameters.
I.e.  X = U - V  where U ~ exp(beta) and V ~ exp(alpha) (where
E(U) = beta, E(V) = alpha).
This makes the distribution of X something like an asymetric Laplace
distribution, with its mode at 0.  (One could shift the mode too, but
that would add a third parameter, which would be de trop.)
Anyhow:  Is this a ``known'' distribution?  Does it have a name?
(I've never seen it mentioned in any of the intro math-stat books
that I've looked into.) If not, can anyone suggest a good name for
it?  (Don't be rude now!)
cheers,
Rolf Turner
[EMAIL PROTECTED]
P. S.  To save you putting pen to paper and working it out,
  the density function is
  { exp(x/alpha)/(alpha + beta) for x = 0
f(x) = {
  { exp(-x/beta)/(alpha + beta) for x = 0
  The mean and variance are mu = beta - alpha and
  sigma^2 = alpha^2 + beta^2 respectfully. :-)
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
_
David Scott Department of Statistics, Tamaki Campus
The University of Auckland, PB 92019
AucklandNEW ZEALAND
Phone: +64 9 373 7599 ext 86830 Fax: +64 9 373 7000
Email:  [EMAIL PROTECTED]
Graduate Officer, Department of Statistics
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Splitting vector into individual elements

2004-09-15 Thread Spencer Graves
Slightly more transparent but arguably uglier: 

 offred.rgb - c(1, 0, 0) * 0.60
 ofr - paste(offred.rgb, collapse=,)
 ofr. - paste(rgb(, ofr, ',names=offred)')
 ofr.
[1] rgb( 0.6,0,0 ,names=\offred\)
 eval(parse(text=ofr.))
  offred
#99

 As long as I can remember eval(parse(text=, this is for me the 
most transparent and works for constructing virtually any R command. 

 hope this helps.  spencer graves
Peter Dalgaard wrote:
Paul Roebuck [EMAIL PROTECTED] writes:
 

Everyone offered 'do.call' as the solution. While that
works, is it to say that there is no means of expanding
the expression as an argument to the original function?
   

Not really. You need an explicit expansion of the argument to a list
somehow, and there's no silver bullet that can convert one function
argument to several. There are solutions without do.call, like
 

offred.rgb - c(1, 0, 0) * 0.60
x - quote(rgb(.,.,.,names=offred))
x[2:4] - as.list(offred.rgb)
eval(x)
   

  offred
#99
but you might find it difficult to explain how it works a year later
 

--
Spencer Graves, PhD, Senior Development Engineer
O:  (408)938-4420;  mobile:  (408)655-4567
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Signs of loadings from princomp on Windows

2004-09-15 Thread Francisco Chamu
I am sorry to insist, but we have three other people that were able to
reproduce the behavior I mentioned.  I have also installed R 1.9.1
from the CRAN binaries on a different Windows machine and again I see
the differents signs as mentioned before.  What would be causing the
difference?

-Francisco


On Tue, 14 Sep 2004 11:04:29 -0600, Tony Plate
[EMAIL PROTECTED] wrote:
 FWIW, I see the same behavior as Francisco on my Windows machine (also an
 installation of the windows binary without trying to install any special
 BLAS libraries):
 
   library(MASS)
   data(painters)
   pca.painters - princomp(painters[ ,1:4])
   loadings(pca.painters)
 
 Loadings:
  Comp.1 Comp.2 Comp.3 Comp.4
 Composition  0.484 -0.376  0.784 -0.101
 Drawing  0.424  0.187 -0.280 -0.841
 Colour  -0.381 -0.845 -0.211 -0.310
 Expression   0.664 -0.330 -0.513  0.432
 
 Comp.1 Comp.2 Comp.3 Comp.4
 SS loadings  1.00   1.00   1.00   1.00
 Proportion Var   0.25   0.25   0.25   0.25
 Cumulative Var   0.25   0.50   0.75   1.00
   pca.painters - princomp(painters[ ,1:4])
   loadings(pca.painters)
 
 Loadings:
  Comp.1 Comp.2 Comp.3 Comp.4
 Composition -0.484 -0.376  0.784 -0.101
 Drawing -0.424  0.187 -0.280 -0.841
 Colour   0.381 -0.845 -0.211 -0.310
 Expression  -0.664 -0.330 -0.513  0.432
 
 Comp.1 Comp.2 Comp.3 Comp.4
 SS loadings  1.00   1.00   1.00   1.00
 Proportion Var   0.25   0.25   0.25   0.25
 Cumulative Var   0.25   0.50   0.75   1.00
   R.version
   _
 platform i386-pc-mingw32
 arch i386
 os   mingw32
 system   i386, mingw32
 status
 major1
 minor9.1
 year 2004
 month06
 day  21
 language R
  
 
 My machine is a dual-processor hp xw8000.
 
 I also get the same results with R 2.0.0 dev as in
   R.version
   _
 platform i386-pc-mingw32
 arch i386
 os   mingw32
 system   i386, mingw32
 status   Under development (unstable)
 major2
 minor0.0
 year 2004
 month09
 day  13
 language R
  
 
 -- Tony Plate
 
 
 
 At Tuesday 10:25 AM 9/14/2004, Prof Brian Ripley wrote:
 On Tue, 14 Sep 2004, Francisco Chamu wrote:
 
   I have run this on both Windows 2000 and XP.  All I did was install
   the binaries from CRAN so I think I am using the standard Rblas.dll.
  
   To reproduce what I see you must run the code at the beginning of the
   R session.
 
 We did, as you said `start a clean session'.
 
 I think to reproduce what you see we have to be using your account on your
 computer.
 
   After the second run, all subsequent runs give the same
   result as the second set.
  
   Thanks,
   Francisco
  
  
   On Tue, 14 Sep 2004 08:29:25 +0200, Uwe Ligges
   [EMAIL PROTECTED] wrote:
Prof Brian Ripley wrote:
 I get the second set each time, on Windows, using the build from CRAN.
 Which BLAS are you using?
   
   
Works also well for me with a self compiled R-1.9.1 (both with standard
Rblas as well as with the Rblas.dll for Athlon CPU from CRAN).
Is this a NT-based version of Windows (NT, 2k, XP)?
   
Uwe
   
   
   
   
 On Tue, 14 Sep 2004, Francisco Chamu wrote:


I start a clean session of R 1.9.1 on Windows and I run the
  following code:


library(MASS)
data(painters)
pca.painters - princomp(painters[ ,1:4])
loadings(pca.painters)

Loadings:
Comp.1 Comp.2 Comp.3 Comp.4
Composition  0.484 -0.376  0.784 -0.101
Drawing  0.424  0.187 -0.280 -0.841
Colour  -0.381 -0.845 -0.211 -0.310
Expression   0.664 -0.330 -0.513  0.432

   Comp.1 Comp.2 Comp.3 Comp.4
SS loadings  1.00   1.00   1.00   1.00
Proportion Var   0.25   0.25   0.25   0.25
Cumulative Var   0.25   0.50   0.75   1.00

However, if I rerun the same analysis, the loadings of the first
component have the opposite sign (see below), why is that?  I have
read the note
in the princomp help that says

The signs of the columns of the loadings and scores are arbitrary,
 and so may differ between different programs for PCA, and even
 between different builds of R.

However, I still would expect the same signs for two runs in the
  same session.


pca.painters - princomp(painters[ ,1:4])
loadings(pca.painters)

Loadings:
Comp.1 Comp.2 Comp.3 Comp.4
Composition -0.484 -0.376  0.784 -0.101
Drawing -0.424  0.187 -0.280 -0.841
Colour   0.381 -0.845 -0.211 -0.310
Expression  -0.664 -0.330 -0.513  0.432

   Comp.1 Comp.2 Comp.3 Comp.4
SS loadings  1.00   1.00   1.00   1.00
Proportion Var   0.25   0.25   0.25   0.25
Cumulative Var   0.25   0.50   0.75   1.00

R.version

 _
platform i386-pc-mingw32
arch i386
os   mingw32
system   i386, mingw32
status
major1
minor9.1
year 2004

[R] Cross-validation for Linear Discrimitant Analysis

2004-09-15 Thread Yu Shao
Hello:

I am new to R and statistics and I have two questions.

First I need help to interpret the cross-validation result from the R
linear discriminant analysis function lda. I did the following:

lda (group ~ Var1 + Var2, CV=T)

where CV=T tells the lda to do cross-validation. The output of lda are
the posterior probabilities among other things, but I can't find an error
term (like delta returned by cv.glm). My question is how to get such an
error term from the output? Can I just simply calculate the prediction
accuracy using the posterior probabilities from the cross-validation, and
use that to measure the quality of the model?

Another question is more basic: how to determine if a lda model is
significant? (There is no p-value.) Thanks,

Yu Shao

Wadsworth Research Center
Department of Health of New York State
Albany, NY 12208

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


Re: [R] Signs of loadings from princomp on Windows

2004-09-15 Thread Tony Plate
You could investigate this yourself by looking at the code of princomp (try 
getAnywhere(princomp.default)).  I'd suggest making a file that in-lines 
the body of princomp.default into the commands you had below.  See if you 
still get the difference.  (I'd be surprised if you didn't).  Then try 
commenting out lines the second pass through the commands produces the same 
results as the first.  The very last thing you commented out might help to 
answer your question What would be causing the
difference?  (The fact that various people chimed in to say they could 
reproduce the behavior that bothered you, but didn't bother dig deeper 
suggests it didn't bother them that much, which further suggests that you 
are the person most motivated by this and thus the best candidate for 
investigating it further...)

-- Tony Plate
At Wednesday 05:07 PM 9/15/2004, Francisco Chamu wrote:
I am sorry to insist, but we have three other people that were able to
reproduce the behavior I mentioned.  I have also installed R 1.9.1
from the CRAN binaries on a different Windows machine and again I see
the differents signs as mentioned before.  What would be causing the
difference?
-Francisco
On Tue, 14 Sep 2004 11:04:29 -0600, Tony Plate
[EMAIL PROTECTED] wrote:
 FWIW, I see the same behavior as Francisco on my Windows machine (also an
 installation of the windows binary without trying to install any special
 BLAS libraries):

   library(MASS)
   data(painters)
   pca.painters - princomp(painters[ ,1:4])
   loadings(pca.painters)

 Loadings:
  Comp.1 Comp.2 Comp.3 Comp.4
 Composition  0.484 -0.376  0.784 -0.101
 Drawing  0.424  0.187 -0.280 -0.841
 Colour  -0.381 -0.845 -0.211 -0.310
 Expression   0.664 -0.330 -0.513  0.432

 Comp.1 Comp.2 Comp.3 Comp.4
 SS loadings  1.00   1.00   1.00   1.00
 Proportion Var   0.25   0.25   0.25   0.25
 Cumulative Var   0.25   0.50   0.75   1.00
   pca.painters - princomp(painters[ ,1:4])
   loadings(pca.painters)

 Loadings:
  Comp.1 Comp.2 Comp.3 Comp.4
 Composition -0.484 -0.376  0.784 -0.101
 Drawing -0.424  0.187 -0.280 -0.841
 Colour   0.381 -0.845 -0.211 -0.310
 Expression  -0.664 -0.330 -0.513  0.432

 Comp.1 Comp.2 Comp.3 Comp.4
 SS loadings  1.00   1.00   1.00   1.00
 Proportion Var   0.25   0.25   0.25   0.25
 Cumulative Var   0.25   0.50   0.75   1.00
   R.version
   _
 platform i386-pc-mingw32
 arch i386
 os   mingw32
 system   i386, mingw32
 status
 major1
 minor9.1
 year 2004
 month06
 day  21
 language R
  

 My machine is a dual-processor hp xw8000.

 I also get the same results with R 2.0.0 dev as in
   R.version
   _
 platform i386-pc-mingw32
 arch i386
 os   mingw32
 system   i386, mingw32
 status   Under development (unstable)
 major2
 minor0.0
 year 2004
 month09
 day  13
 language R
  

 -- Tony Plate



 At Tuesday 10:25 AM 9/14/2004, Prof Brian Ripley wrote:
 On Tue, 14 Sep 2004, Francisco Chamu wrote:
 
   I have run this on both Windows 2000 and XP.  All I did was install
   the binaries from CRAN so I think I am using the standard Rblas.dll.
  
   To reproduce what I see you must run the code at the beginning of the
   R session.
 
 We did, as you said `start a clean session'.
 
 I think to reproduce what you see we have to be using your account on your
 computer.
 
   After the second run, all subsequent runs give the same
   result as the second set.
  
   Thanks,
   Francisco
  
  
   On Tue, 14 Sep 2004 08:29:25 +0200, Uwe Ligges
   [EMAIL PROTECTED] wrote:
Prof Brian Ripley wrote:
 I get the second set each time, on Windows, using the build 
from CRAN.
 Which BLAS are you using?
   
   
Works also well for me with a self compiled R-1.9.1 (both with 
standard
Rblas as well as with the Rblas.dll for Athlon CPU from CRAN).
Is this a NT-based version of Windows (NT, 2k, XP)?
   
Uwe
   
   
   
   
 On Tue, 14 Sep 2004, Francisco Chamu wrote:


I start a clean session of R 1.9.1 on Windows and I run the
  following code:


library(MASS)
data(painters)
pca.painters - princomp(painters[ ,1:4])
loadings(pca.painters)

Loadings:
Comp.1 Comp.2 Comp.3 Comp.4
Composition  0.484 -0.376  0.784 -0.101
Drawing  0.424  0.187 -0.280 -0.841
Colour  -0.381 -0.845 -0.211 -0.310
Expression   0.664 -0.330 -0.513  0.432

   Comp.1 Comp.2 Comp.3 Comp.4
SS loadings  1.00   1.00   1.00   1.00
Proportion Var   0.25   0.25   0.25   0.25
Cumulative Var   0.25   0.50   0.75   1.00

However, if I rerun the same analysis, the loadings of the first
component have the opposite sign (see below), why is that?  I have
read the note
in the princomp help that says

The signs of the columns of the loadings and scores are 
arbitrary,
 and so may differ between 

[R] Newbie q. need some help understanding this code.

2004-09-15 Thread sean kim
dear all. 

Would someone be kind and willing to explain the code
below for a person who has never used R?  ( that is if
one has enough time and inclination) 

It implements gillepsie's stochastic algorithm for
Lotka Volterra model. 

What would help me tremendously is to see the
breakdown of the line by line code into plain english.


thanks for any insights or other comments. 

sean 



library(stepfun)

lv - function(N=1000,cvec=c(1,0.005,0.6),x=c(50,100))
{
m-length(cvec)
n-length(x)
xmat-matrix(nrow=N+1,ncol=n)
tvec-vector(numeric,N)
h-vector(numeric,m)
t-0
xmat[1,]-x
for (i in 1:N) {
h[1]-cvec[1]*x[1]
h[2]-cvec[2]*x[1]*x[2]
h[3]-cvec[3]*x[2]
h0=sum(h)
tp-rexp(1,h0)
t-t+tp
u-runif(1,0,1)
if ( u  h[1]/h0 ) {
x[1] - x[1]+1
} else if ( u  (h[1]+h[2])/h0 ) {
x[1] - x[1]-1
x[2] - x[2]+1
} else {
x[2] - x[2]-1
}
xmat[i+1,]-x
tvec[i]-t
}
   
list(stepfun(tvec,xmat[,1]),stepfun(tvec,xmat[,2]))
}



results - lv(N=1)


op-par(mfrow=c(2,1))
plot(results[[1]],do.points=True)
plot(results[[2]],do.points=False)
par(op)

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


Re: [R] Newbie q. need some help understanding this code.

2004-09-15 Thread Kevin Wang
Hi,
sean kim wrote:
thanks for any insights or other comments. 
I would suggest that you run the code line by line, to see what it does 
yourself.  It is the best way to learn!

library(stepfun)
This just loads the package stepfun.
lv - function(N=1000,cvec=c(1,0.005,0.6),x=c(50,100))
{
Declaration of a function named lv, where N, cvec and x are the 
parameters/arguments with default values.

m-length(cvec)
m is a new variable, with in this case is just a number, the length of 
the vector cvec.  In this case, 3.

n-length(x)
Ditto.
xmat-matrix(nrow=N+1,ncol=n)
Generate a matrix with N + 1 (1001) rows and  n columns.
tvec-vector(numeric,N)
h-vector(numeric,m)
Generating two vectors.
t-0
xmat[1,]-x
Assign x to the first row of the xmat matrix.
for (i in 1:N) {
h[1]-cvec[1]*x[1]
The first element in cvec and x, multiply them together and the result 
becomes the first element of h.

h[2]-cvec[2]*x[1]*x[2]
h[3]-cvec[3]*x[2]
Ditto.
h0=sum(h)
Get the sum of h.
tp-rexp(1,h0)
Generating an exponential random number with a rate of h0.
t-t+tp
u-runif(1,0,1)
Generating a uniform random number, [0, 1]
if ( u  h[1]/h0 ) {
x[1] - x[1]+1
If the uniform random number, u, is smaller than h[1]/h0 then do...
otherwise do the following...
} else if ( u  (h[1]+h[2])/h0 ) {
x[1] - x[1]-1
x[2] - x[2]+1
} else {
x[2] - x[2]-1
}
xmat[i+1,]-x
Increment to the next row.
tvec[i]-t
}
   
list(stepfun(tvec,xmat[,1]),stepfun(tvec,xmat[,2]))
}
Put things together into a list.
Again, check these line by line and you'll have a better understanding!
Kev
--
Ko-Kang Kevin Wang
PhD Student
Centre for Mathematics and its Applications
Building 27, Room 1004
Mathematical Sciences Institute (MSI)
Australian National University
Canberra, ACT 0200
Australia
Homepage: http://wwwmaths.anu.edu.au/~wangk/
Ph (W): +61-2-6125-2431
Ph (H): +61-2-6125-7407
Ph (M): +61-40-451-8301
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Splitting vector into individual elements

2004-09-15 Thread Gabor Grothendieck
Gabor Grothendieck ggrothendieck at myway.com writes:

: 
: Paul Roebuck roebuck at odin.mdacc.tmc.edu writes:
: 
: : 
: : On Wed, 15 Sep 2004, Peter Dalgaard wrote:
: : 
: :  Paul Roebuck roebuck at odin.mdacc.tmc.edu writes:
: : 
: :   Is there a means to split a vector into its individual
: :   elements without going the brute-force route for arguments
: :   to a predefined function call?
: :  
: :   offred.rgb - c(1, 0, 0) * 0.60;
: :  
: :   ## Brute force style
: :   offred.col - rgb(offred.rgb[1],
: : offred.rgb[2],
: : offred.rgb[3],
: : names = offred)
: :   ## Desired style
: :   offred.col - rgb(silver.bullet(offred.rgb),
: : names = offred)
: : 
: :  The closest is probably this:
: : 
: :  offred.col - do.call(rgb, c(as.list(offred.rgb),
: : list(names=offred)))
: : 
: : 
: : Everyone offered 'do.call' as the solution. While that
: : works, is it to say that there is no means of expanding
: : the expression as an argument to the original function?
: 
: This is not a true answer to the question of expanding a list
: into arguments without using do.call but it does allow you to 
: carry out either syntax in this particular case using S3
: dispatch:
: 
: R silver.bullet - as.list
: R rgb - function(x, ...) UseMethod(rgb)
: R rgb.list - function(x, ...) rgb(x[[1]],x[[2]],x[[3]],...)
: R rgb.default - graphics::rgb
: 
: R offred.rgb - c(1, 0, 0) * 0.60;
: R # original syntax
: R rgb(offred.rgb[1], offred.rgb[2], offred.rgb[3], names = offred)
:offred 
: #99 
: R # list syntax
: R rgb(silver.bullet(offred.rgb), names = offred)
:offred 
: #99

Here is a second, different approach.  

Again, it is not exactly what you are asking for since silver.bullet,
here called flatten.args, is applied to the function rather than the
argument in question and operates on all arguments, not just one but
I think its closer in spirit to your query than my previous solution:

flatten.args - function(f)
function(...) {
L - list()
for(i in list(...)) L - c(L, unlist(i))
do.call(as.character(substitute(f)), L)
}

flatten.args(rgb)(0.6 * c(1,0,0), name = offred)

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


Re: [R] Cross-validation for Linear Discrimitant Analysis

2004-09-15 Thread Prof Brian Ripley
On Wed, 15 Sep 2004, Yu Shao wrote:

 I am new to R and statistics and I have two questions.

Perhaps then you need to start by explaining why you are using LDA.
Please take a good look at the posting guide.

 First I need help to interpret the cross-validation result from the R
 linear discriminant analysis function lda. 

You mean Professor Ripley's function lda in package MASS, I guess.

 I did the following:
 
 lda (group ~ Var1 + Var2, CV=T)

R allows you to use meaningful names, so please do so.

 where CV=T tells the lda to do cross-validation. The output of lda are
 the posterior probabilities among other things, but I can't find an error
 term (like delta returned by cv.glm). My question is how to get such an
 error term from the output? Can I just simply calculate the prediction
 accuracy using the posterior probabilities from the cross-validation, and
 use that to measure the quality of the model?

cv.glm as in Dr Canty's package boot?  If you are trying to predict
classifications, LDA is not the right tool, and LOO CV probably is not
either.  There is no unique definition of `error term' (true for cv.glm as
well), and people have written whole books about how to assess
classifiers.  LDA is about `discrimination' not `allocation' in the jargon 
used ca 1960.

 Another question is more basic: how to determine if a lda model is
 significant? (There is no p-value.) Thanks,

Please do read the references on the ?lda page.  It's not a useful
question, as LDA is about discriminating between populations and makes the
unrealistic assumption of multivariate normality.  (Analogously for linear
regression, there are ways to test if that is (statistically)
`significant', but knowledgable users almost never do so.)

Perhaps more realistic advice is to suggest you seek some statistical 
consultancy.

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

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


[R] sapply() with cat()

2004-09-15 Thread Kevin Wang
Hi,
Suppose I've got a data frame:
   inter.df
  V1  V2  V3  V4  V5  V6  V7  V8  V9 V10 V11 V12
  1  3.3  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
  2  0.0 0.1  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
  3  1.0 0.9 0.2  NA  NA  NA  NA  NA  NA  NA  NA  NA
  4  1.6 0.0 2.9 0.7  NA  NA  NA  NA  NA  NA  NA  NA
  5  0.0 0.1 2.9 0.1 0.1  NA  NA  NA  NA  NA  NA  NA
  6  2.4 1.0 0.6 0.4 1.9 0.1  NA  NA  NA  NA  NA  NA
  7  2.8 1.4 1.2 7.5 0.0 0.0 4.2  NA  NA  NA  NA  NA
  8  0.3 3.1 0.8 3.7 5.7 0.0 0.8 0.0  NA  NA  NA  NA
  9  0.1 2.9 0.3 1.3 0.2 0.2 0.5 1.4 0.9  NA  NA  NA
  10 0.8 2.6 0.0 0.0 0.1 4.1 0.8 4.3 0.6 2.2  NA  NA
  11 0.0 4.0 0.0 0.3 0.5 0.9 0.0 1.5 0.2 0.7 0.8  NA
  12 0.6 0.8 0.3 0.0 0.2 1.2 0.0 0.8 1.5 0.9 0.4   0
  
   foo - function(x) {
  +   ifelse(x  3.8, NA, x)
  + }
   sapply(inter.df, foo)
 V1  V2  V3  V4  V5  V6  V7  V8  V9 V10 V11 V12
   [1,] 3.3  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
   [2,] 0.0 0.1  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
   [3,] 1.0 0.9 0.2  NA  NA  NA  NA  NA  NA  NA  NA  NA
   [4,] 1.6 0.0 2.9 0.7  NA  NA  NA  NA  NA  NA  NA  NA
   [5,] 0.0 0.1 2.9 0.1 0.1  NA  NA  NA  NA  NA  NA  NA
   [6,] 2.4 1.0 0.6 0.4 1.9 0.1  NA  NA  NA  NA  NA  NA
   [7,] 2.8 1.4 1.2  NA 0.0 0.0  NA  NA  NA  NA  NA  NA
   [8,] 0.3 3.1 0.8 3.7  NA 0.0 0.8 0.0  NA  NA  NA  NA
   [9,] 0.1 2.9 0.3 1.3 0.2 0.2 0.5 1.4 0.9  NA  NA  NA
  [10,] 0.8 2.6 0.0 0.0 0.1  NA 0.8  NA 0.6 2.2  NA  NA
  [11,] 0.0  NA 0.0 0.3 0.5 0.9 0.0 1.5 0.2 0.7 0.8  NA
  [12,] 0.6 0.8 0.3 0.0 0.2 1.2 0.0 0.8 1.5 0.9 0.4   0
Up to here, sapply() does what I want, replacing values that are greater 
than 3.8 to NA, as per my foo() function.  But...

   goo - function(x) {
  +   ifelse(x  3.8, cat(\\textbf{, x, }, sep = ), x)
  + }
   sapply(inter.df, goo)
  \textbf{NA0.10.900.111.43.12.92.640.8}Error in [-(`*tmp*`, test, 
value = rep(yes, length.out = length(ans))[test]) :
  incompatible types

If instead whenever I get a value that's greater than 3.8 I want to 
change it what goo() does, e.g. the value 4.0 should become:
  \textbf{4.0}
but I got the above error.

What I'm intending is to then pass the inter.df (with the \textbf{} 
markups) into xtable() so the values will be bolded in the resulting 
LaTeX table (as so far, I cannot see how I can achieve bolding only 
certain numbers with xtable() alone).

Any suggestions will be welcome!
Kevin
--
Ko-Kang Kevin Wang
PhD Student
Centre for Mathematics and its Applications
Building 27, Room 1004
Mathematical Sciences Institute (MSI)
Australian National University
Canberra, ACT 0200
Australia
Homepage: http://wwwmaths.anu.edu.au/~wangk/
Ph (W): +61-2-6125-2431
Ph (H): +61-2-6125-7407
Ph (M): +61-40-451-8301
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] sapply() with cat()

2004-09-15 Thread Prof Brian Ripley
On Thu, 16 Sep 2004, Kevin Wang wrote:

 Hi,
 
 Suppose I've got a data frame:
 inter.df
V1  V2  V3  V4  V5  V6  V7  V8  V9 V10 V11 V12
1  3.3  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
2  0.0 0.1  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
3  1.0 0.9 0.2  NA  NA  NA  NA  NA  NA  NA  NA  NA
4  1.6 0.0 2.9 0.7  NA  NA  NA  NA  NA  NA  NA  NA
5  0.0 0.1 2.9 0.1 0.1  NA  NA  NA  NA  NA  NA  NA
6  2.4 1.0 0.6 0.4 1.9 0.1  NA  NA  NA  NA  NA  NA
7  2.8 1.4 1.2 7.5 0.0 0.0 4.2  NA  NA  NA  NA  NA
8  0.3 3.1 0.8 3.7 5.7 0.0 0.8 0.0  NA  NA  NA  NA
9  0.1 2.9 0.3 1.3 0.2 0.2 0.5 1.4 0.9  NA  NA  NA
10 0.8 2.6 0.0 0.0 0.1 4.1 0.8 4.3 0.6 2.2  NA  NA
11 0.0 4.0 0.0 0.3 0.5 0.9 0.0 1.5 0.2 0.7 0.8  NA
12 0.6 0.8 0.3 0.0 0.2 1.2 0.0 0.8 1.5 0.9 0.4   0

 foo - function(x) {
+   ifelse(x  3.8, NA, x)
+ }
 sapply(inter.df, foo)
   V1  V2  V3  V4  V5  V6  V7  V8  V9 V10 V11 V12
 [1,] 3.3  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
 [2,] 0.0 0.1  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
 [3,] 1.0 0.9 0.2  NA  NA  NA  NA  NA  NA  NA  NA  NA
 [4,] 1.6 0.0 2.9 0.7  NA  NA  NA  NA  NA  NA  NA  NA
 [5,] 0.0 0.1 2.9 0.1 0.1  NA  NA  NA  NA  NA  NA  NA
 [6,] 2.4 1.0 0.6 0.4 1.9 0.1  NA  NA  NA  NA  NA  NA
 [7,] 2.8 1.4 1.2  NA 0.0 0.0  NA  NA  NA  NA  NA  NA
 [8,] 0.3 3.1 0.8 3.7  NA 0.0 0.8 0.0  NA  NA  NA  NA
 [9,] 0.1 2.9 0.3 1.3 0.2 0.2 0.5 1.4 0.9  NA  NA  NA
[10,] 0.8 2.6 0.0 0.0 0.1  NA 0.8  NA 0.6 2.2  NA  NA
[11,] 0.0  NA 0.0 0.3 0.5 0.9 0.0 1.5 0.2 0.7 0.8  NA
[12,] 0.6 0.8 0.3 0.0 0.2 1.2 0.0 0.8 1.5 0.9 0.4   0
 
 Up to here, sapply() does what I want, replacing values that are greater 
 than 3.8 to NA, as per my foo() function.  But...
 
 goo - function(x) {
+   ifelse(x  3.8, cat(\\textbf{, x, }, sep = ), x)
+ }
 sapply(inter.df, goo)
\textbf{NA0.10.900.111.43.12.92.640.8}Error in [-(`*tmp*`, test, 
 value = rep(yes, length.out = length(ans))[test]) :
incompatible types
 
 If instead whenever I get a value that's greater than 3.8 I want to 
 change it what goo() does, e.g. the value 4.0 should become:
\textbf{4.0}
 but I got the above error.

cat writes to connection (here stdout) and then returns NULL.  I think you 
meant to use paste():

goo - function(x)ifelse(x  3.8, paste(\\textbf{, x, }, sep = ), x)

works.

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

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


Re: [R] sapply() with cat()

2004-09-15 Thread Gabor Grothendieck
Prof Brian Ripley ripley at stats.ox.ac.uk writes:

: 
: On Thu, 16 Sep 2004, Kevin Wang wrote:
: 
:  Hi,
:  
:  Suppose I've got a data frame:
:  inter.df
: V1  V2  V3  V4  V5  V6  V7  V8  V9 V10 V11 V12
: 1  3.3  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
: 2  0.0 0.1  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
: 3  1.0 0.9 0.2  NA  NA  NA  NA  NA  NA  NA  NA  NA
: 4  1.6 0.0 2.9 0.7  NA  NA  NA  NA  NA  NA  NA  NA
: 5  0.0 0.1 2.9 0.1 0.1  NA  NA  NA  NA  NA  NA  NA
: 6  2.4 1.0 0.6 0.4 1.9 0.1  NA  NA  NA  NA  NA  NA
: 7  2.8 1.4 1.2 7.5 0.0 0.0 4.2  NA  NA  NA  NA  NA
: 8  0.3 3.1 0.8 3.7 5.7 0.0 0.8 0.0  NA  NA  NA  NA
: 9  0.1 2.9 0.3 1.3 0.2 0.2 0.5 1.4 0.9  NA  NA  NA
: 10 0.8 2.6 0.0 0.0 0.1 4.1 0.8 4.3 0.6 2.2  NA  NA
: 11 0.0 4.0 0.0 0.3 0.5 0.9 0.0 1.5 0.2 0.7 0.8  NA
: 12 0.6 0.8 0.3 0.0 0.2 1.2 0.0 0.8 1.5 0.9 0.4   0
: 
:  foo - function(x) {
: +   ifelse(x  3.8, NA, x)
: + }
:  sapply(inter.df, foo)
:V1  V2  V3  V4  V5  V6  V7  V8  V9 V10 V11 V12
:  [1,] 3.3  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
:  [2,] 0.0 0.1  NA  NA  NA  NA  NA  NA  NA  NA  NA  NA
:  [3,] 1.0 0.9 0.2  NA  NA  NA  NA  NA  NA  NA  NA  NA
:  [4,] 1.6 0.0 2.9 0.7  NA  NA  NA  NA  NA  NA  NA  NA
:  [5,] 0.0 0.1 2.9 0.1 0.1  NA  NA  NA  NA  NA  NA  NA
:  [6,] 2.4 1.0 0.6 0.4 1.9 0.1  NA  NA  NA  NA  NA  NA
:  [7,] 2.8 1.4 1.2  NA 0.0 0.0  NA  NA  NA  NA  NA  NA
:  [8,] 0.3 3.1 0.8 3.7  NA 0.0 0.8 0.0  NA  NA  NA  NA
:  [9,] 0.1 2.9 0.3 1.3 0.2 0.2 0.5 1.4 0.9  NA  NA  NA
: [10,] 0.8 2.6 0.0 0.0 0.1  NA 0.8  NA 0.6 2.2  NA  NA
: [11,] 0.0  NA 0.0 0.3 0.5 0.9 0.0 1.5 0.2 0.7 0.8  NA
: [12,] 0.6 0.8 0.3 0.0 0.2 1.2 0.0 0.8 1.5 0.9 0.4   0
:  
:  Up to here, sapply() does what I want, replacing values that are greater 
:  than 3.8 to NA, as per my foo() function.  But...
:  
:  goo - function(x) {
: +   ifelse(x  3.8, cat(\\textbf{, x, }, sep = ), x)
: + }
:  sapply(inter.df, goo)
: \textbf{NA0.10.900.111.43.12.92.640.8}Error in [-(`*tmp*`, test, 
:  value = rep(yes, length.out = length(ans))[test]) :
: incompatible types
:  
:  If instead whenever I get a value that's greater than 3.8 I want to 
:  change it what goo() does, e.g. the value 4.0 should become:
: \textbf{4.0}
:  but I got the above error.
: 
: cat writes to connection (here stdout) and then returns NULL.  I think you 
: meant to use paste():
: 
: goo - function(x)ifelse(x  3.8, paste(\\textbf{, x, }, sep = ), x)
: 
: works.
: 


Also, if you convert your data frame to a matrix you won't need sapply:

   inter.mat - as.matrix(inter.df)
   ifelse(inter.mat  3.8, 
  paste(\\textbf{, inter.mat, }, sep = ), 
  inter.mat)

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