[R] Odp: Getting mean in a dataframe

2007-11-14 Thread Petr PIKAL
Hi

[EMAIL PROTECTED] napsal dne 13.11.2007 23:29:25:

> Sorry to hijack this thread. I have a similar but slightly different 
> situation. Using the original poster's example, how to elegantly get 
> the mean of column V2 when column V1 is either A or C and F1 is 0?

I do not fully understand. You want to select subset of a data frame and 
then to compute column mean?

If yes just use ?subset or construct a logical vector according to your 
condition

lll <- (V1 == "A" | V1 == "B") & F1 == 0
mean(data[lll,"V2"])

Regards
Petr


> 
> Thanks,
> Gang
> 
> 
> On Nov 13, 2007, at 5:30 AM, Petr PIKAL wrote:
> > Hi
> >
> > [EMAIL PROTECTED] napsal dne 13.11.2007 10:59:09:
> >
> >> Dear list,
> >>
> >> I have this dataframe
> >>
> >>V1 V2 F1
> >> 1  A  2  0
> >> 2  A  3  0
> >> 3  A  4  1
> >> 4  B  3  0
> >> 5  B  2  1
> >> 6  C  6  0
> >> 7  C  2  0
> >> 8  C  6  0
> >>
> >> and would like to calculate a new column
> >> with mean-values, following this rule
> >>
> >> 1. If F1 = 0 calculate the mean from V2
> >> for each factor in V1.
> >>
> >> 2. If F1 = 1, then F1_mean = 0
> >>
> >> So, the new DF should look like this
> >>
> >>V1 V2 F1 F1_mean
> >> 1  A  2  0 2.5
> >> 2  A  3  0 2.5
> >> 3  A  4  1 0.0
> >> 4  B  3  0 3.0
> >> 5  B  2  1 0.0
> >> 6  C  6  0 7.0
> >> 7  C  2  0 7.0
> >> 8  C  6  0 7.0
> >
> > I would use ave for computing mean for combination of V1 and F1 and 
> > then I
> > put all values of F1_mean for which F1 is 1 to 0
> >
> > test$F1_mean <- ave(test$V2, test$V1, factor(test$F1), FUN = mean)
> > test$F1_mean[test$F1==1] <- 0
> >
> > Regards
> > Petr
> >
> >>
> >> Thank you for any help!
> >>
> >> Patrick Hausmann
> 
> 
> 
>[[alternative HTML version deleted]]
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

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


[R] geotiff calculations

2007-11-14 Thread Laura Poggio
Dear list,
I have to compare two digital elevation models in raster format (geotiff).
I then have to calculate the differences in altitude for each cell and make
some statistics (basic as mean, median, std, range but also more "advanced"
as RMSE) on that.
I do not know very much how to proceed:
1) is it possible to import the geotiff in R? If so with which package? if
not which is the best way to import such files?
2) is it better to perform the calculations of the differences in a GIS
software and then to use R only for statistical analysis? or it is better to
do everything in R?
3) is there any specific package for doing this kind of analysis?

Thank you very much in advance

Laura

[[alternative HTML version deleted]]

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


[R] 50% Off Sale, Wednesday Only

2007-11-14 Thread Kohls.com
http://e.kohls.com/a/tBHOp0BBBZVhBBiZzO2BVHiEYU4/[EMAIL PROTECTED]



50% Off Sale
Wednesday Only

Give a warm welcome to cooler temperatures! Stock up and save on
sweaters for everyone in all the latest fashions! Shop your local
Kohl's Department Store for great prices on everything you need to
greet the new season in style!



View Today's Ad online to see what's on sale storewide!



Kohls.com

FEATURING
50% Off Sweaters



It's time to Deck-orate the Halls!
Browse our Home for the Holidays gift shop.



$25,000 GIVEAWAY
Enter for your chance to WIN BIG during the Kohls.com $25,000 Giveaway!
Now through December 29 only.



60-80% Clearance--NEW MARKDOWNS TAKEN
Shop Kohls.com Clearance for up to 60-80% off* original prices on items
from every department. But hurry--quantities are limited and deals
like these won't last.



*Clearance prices represent savings off original prices. Interim
markdowns may have been taken. Sorry, no price adjustments on
prior purchases.

This mailbox is unattended, so please do not reply to this message.
Instead, e-mail us at [EMAIL PROTECTED], or write to us at Kohl's
Department Stores, Attention: Customer Service, N54 W13600 Woodale Drive,
Menomonee Falls, WI 53051. If you no longer wish to receive e-mails from
Kohls.com, unsubscribe by pasting this link into the Address field of
your Internet browser: 

http://e.kohls.com/a/tBHOp0BBBZVhBBiZzO2BVHiEYU4/[EMAIL PROTECTED]

50% Off Sale prices good November 14, 2007.


[[alternative HTML version deleted]]

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


[R] Hottelings T2-test for multivariate lingitudinal data

2007-11-14 Thread Mauricio Malfert
Dear R-users
I've simulated a longitudinal multivariate normal data set from which
I've simulated missing-patterns such as MCAR MAR and a simple kind of
non-MAR. I've imputated the values so I now have 'complete' data sets. I'm
trying to perform a T2-test as done in the multivariate case under th
enormal assumption. Is there something I've to think about when performing
this test on a longitudinal multivariate data set when dependencies between
the time-dependent variables are present. Any help will much appreciated.

Thanks in advance

Mauricio Malfert

[[alternative HTML version deleted]]

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


Re: [R] cronbach's alpha

2007-11-14 Thread Dimitris Rizopoulos
look also at functions cronbach.alpha() and descript(), in the 'ltm' 
package.

Best,
Dimitris


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

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


- Original Message - 
From: "raymond chiruka" <[EMAIL PROTECTED]>
To: "r" <[EMAIL PROTECTED]>
Sent: Wednesday, November 14, 2007 12:33 AM
Subject: [R] cronbach's alpha


> hie
>1...i'm trying to carryout a relibility testusing cronbach's 
> alpha what fuctin do i use.
>
>2.. this is more of a statistical question.if the alpha value 
> for all the variables  is negative what does it mean. and if the 
> alpha value is negative for all tyha variables but is greater than 
> 0.7  for some sections of the variables what does that mean
>
>  thanks in advance
>
>
>
> -
> Never miss a thing.   Make Yahoo your homepage.
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 


Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm

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


[R] Newbie question about data preparation

2007-11-14 Thread Heidemeier Dr, Joachim
I'm new to the R language and still struggling with compactness of R. I
haven't got the right compass into the Documentation, too. So, please
apologize the possibly stupidiy of my question.
I have the following problem:
I have two data sets combined in a data.frame
 x  y
1.3 2.2
2.5 3.4
3.1 3.7
8.2 9.5 
7.5 8.3
For the analyses of the data I want to group one column (like the
classes in a histogram).
So I'd like to add one column with the center of each group with width=2
for an x value in the interval of the class.
So the output should look like 
  x y   x-grouped
1.3 2.2 1
2.5 3.4 3
3.1 3.7 3   
8.2 9.5 9
7.5 8.3 7
What would be the R'ish idiom to do this operation?

TIA
joachim
-- 
Dr. Joachim Heidemeier
c/o Umweltbundesamt FG II 2.2
Tel.: +49340 2103-2780 
eMail: [EMAIL PROTECTED]
 

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


[R] Help about nnet library

2007-11-14 Thread Pierre-Henri BUSSIERE
Hello,

First of all I am french, so please forgive me, if there are some big 
language mistakes in my sentences.

I think, it is the good mail address to send my question, if not please 
tell me and forgive me.

I am working on a project, and I use the nnet library. Our customers do 
not want us to install R on their machine, so we just use R for making 
the training of our neurons network. This part is good and we have the 
weights for each neurons and each links.

The neuron network we used is a very simple one : 4 inputs for *one 
linear output* unit with just one hidden layer (with 5 neurons) between 
the inputs and the ouputs. So in R it is translated by the following 
network : 4-5-1. We have also a decay of 0.3, but I do not think that 
fact interfers in our computation.

And our problem is to recompute the output with those weights.

When we program it or making the computation by hand, we do not find the 
same result than R. We check our program, we remake the computation by 
hand and we find always the same result.

So one of our hypothesis is that, we do not use the good transfert 
function or the good activation function.

For example for the neuron 1 of the hidden layer we have the weights as 
follows :
b(which is I suppose the input set to one)->h1=w0
i1->h1=w1
i2->h1=w2
i3->h1=w3
i4->h1=w4

So first we make this computation :
sum=w0 + i1*w1 + i2*w2 + i3*w3 + i4*w4 (the transfert function)

And for computing h1 we use the classic sigmoid function :

h1=1/(1+exp(-sum)).

And for calculating the output we have those weights :
b->o=wo0
i1->o=wo1
i2->o=wo2
i3->o=wo3
i4->o=wo4
h1->o=wo5
h2->o=wo6
h3->o=wo7
h4->o=wo8
h5->o=wo9

And the ouput is computed as follows :

o=wo0+wo1*i1+wo2*i2+wo3*i3+wo4*i4+wo5*h1+wo6*h2+wo7*h3+wo8*h4+wo9*h5.

So I would like to know, where is our error. Is it the sigmoid function 
that we used ? The transfert function ? Please answer me, I really need 
to know it.

Thanks in advance.

Pierre-Henri BUSSIERE

-- 
Pierre-Henri BUSSIERE
Ingénieur Informatique

NUMTECH
6, allée Alan Turing
Parc Technologique de la Pardieu
BP 30242
63175 AUBIERE Cedex
FRANCE

Tél. 04.73.28.75.95
Fax. 04.73.28.75.99 


[[alternative HTML version deleted]]

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


[R] R Crashes on certain calls of Adapt

2007-11-14 Thread Omkar Muralidharan
I'm having trouble with adapt. I'm trying to use it in a Bayesian setting,
to integrate the posterior distribution, and to find posterior means. I
tried using the following script, and things went ok:

data = rnorm(100,0.2,1.1)
data = c(data,rnorm(10,3,1))
data = data[abs(data)<2*sd(data)]

prior = function(x){
dgamma(x[2],shape=2,scale=1)*dnorm(x[1],0,.5)
}

liklihood = function(x,val){
prod(dnorm(val,m=x[1],sd=sqrt(.5/x[2])))#/(pnorm(2,x[1],x[2])-pnorm(-2,x[1],x[2]))^length(val)
}

unsc.post = function(x,val){
prior(x)*liklihood(x,val)
}

cons = adapt(2,c(-5,0),c(5,10),f=unsc.post,min=1e+04,max=5e+05,val=data)

However, if I try to add a skewness parameter using sn as follows, I get an
error.

data = rnorm(20,0.2,1.1)
data = c(data,rnorm(1,3,1))
data = data[abs(data)<2*sd(data)]

prior = function(x){
dgamma(x[2],shape=2,scale=1)*dnorm(x[1],0,.5)*dnorm(x[3],0,.1)
}

liklihood = function(x,val){
prod(dsn(val,loc=x[1],sc=sqrt(.5/x[2]),sh=x[3]))
}

unsc.post = function(x,val){
prior(x)*liklihood(x,val)
}

cons = adapt(2,c(-5,0,-.5),c(5,10,.5),f=unsc.post,max=5e+05,val=data)

This produces the error:

 *** caught segfault ***
address 0x7eaa8580, cause 'memory not mapped'

Traceback:
 1: .C("cadapt", as.integer(ndim), as.double(lower), as.double(upper),
minpts = as.integer(minpts), maxpts = as.integer(maxpts), ff, rho =
environment(), as.double(eps), relerr = double(1), lenwrk =
as.integer(lenwrk),
value = double(1), ifail = integer(1), PACKAGE = "adapt")
 2: adapt(2, c(-5, 0, -0.5), c(5, 10, 0.5), f = unsc.post, max = 5e+05,
val = data)

On a perhaps related note, if I increase the length of data to 1000+10 in
the original example, I get a garbage answer and:

Warning message:
Ifail=2, lenwrk was too small. -- fix adapt() !
 Check the returned relerr! in: adapt(2, c(-5, 0), c(5, 10), f = unsc.post,
min = 1, max = 5e+05.

And If I uncomment the line there with data size 100+10, adapt takes forever
and gives me an error saying maxpts is too low (I don't think this is
related to the previous problem, but I thought I'd provide the information
in case it is useful).

Any help would be greatly appreciated.

Thanks,
Omkar Muralidharan

[[alternative HTML version deleted]]

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


Re: [R] time plotting problem

2007-11-14 Thread Jim Lemon
John Kane wrote:
> I clearly spoke too soon.  
> 
> With the actual data I am not getting sensible x-axis
> units.  The program with the actual data below.  Graph
> output is here:
> http://ca.geocities.com/jrkrideau/R/hd.png .  
> 
> I seem to be getting only a single entry for the
> x-axis of "2007". However dates range from
>  First DateLast Date 
> "2006-09-26" "2007-11-10" 
> 
> I must be missing something blindingly obvious but I
> don't see it. 
> 
> Thanks for any suggestions.
> 
> 
> Actual data and test program
> 
> mydata <-
> read.table("http://ca.geocities.com/jrkrideau/R/heartdata.txt";,
> sep="\t",
>  header=FALSE)
> mydata[,1] <- as.Date(mydata[,1],"%m/%d/%y")
> names(mydata) <- Cs(dates, sy,dys,pulse, weight)
> minmax  <-  c(max(mydata[,2]),min(mydata[,2]),
> max(mydata[,3]),min(mydata[,3]))
> names(minmax) <- c("maxsys", "minsys", "maxdsy",
> "mindys") ; minmax 
> 
> min.max.dates <- c(min(mydata[,1]), max(mydata[,1]))
> names(min.max.dates) <- c("First Date", "Last Date");
> min.max.dates
> 
> ss <- lm(mydata[,2]~mydata[,1])
> plot(mydata[,1],mydata[,2], xlab="Dates", ylab="Blood
> pressure",  
>ylim=c(minmax[4], minmax[1]), col= "red",
> type="l")
> abline (ss, col ="yellow")
> dd <- lm(mydata[,3]~mydata[,1])
> points(mydata[,1], mydata[,3], type="l", col="blue" )
> abline(dd, col= "yellow", lwd=2)
> 
> 
> 
Hi John,

It can be done like this:

plot(mydata[,1],mydata[,2], xlab="Dates", ylab="Blood
  pressure",ylim=c(minmax[4], minmax[1]), col= "red",
  type="l",xaxt="n")
...
axisdates<-c("2006-11-01","2007-01-01","2007-03-01",
  "2007-05-01","2007-07-01","2007-09-01","2007-11-01")
staxlab(at=as.Date(axisdates),labels=axisdates)

Jim

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


Re: [R] combine two dataframe

2007-11-14 Thread sun
Thanks all for the answers. Both Merge and sqldf works perfectly for me. 
Well, I feel sqldf run a littile bit slower. And I failed to install this 
package (sqldf ) on my linux box.

Denver, your approach also works, but in my case, data frame A has much more 
rows then B, so B has to be duplicated many many times.

kind regards,
Sun

- Original Message - 
From: "Gabor Grothendieck" <[EMAIL PROTECTED]>
To: "sun" <[EMAIL PROTECTED]>
Cc: <[EMAIL PROTECTED]>
Sent: Tuesday, November 13, 2007 6:07 PM
Subject: Re: [R] combine two dataframe


> Try this:
>
>> A <- data.frame(a1 = c(1, 2, 1), a2 = c(2, 3, 3), a3 = c(3, 1, 2))
>> B <- data.frame(b1 = 1:2, b2 = 2:1)
>>
>> library(sqldf)
>> sqldf("select * from A, B")
>  a1 a2 a3 b1 b2
> 1  1  2  3  1  2
> 2  1  2  3  2  1
> 3  2  3  1  1  2
> 4  2  3  1  2  1
> 5  1  3  2  1  2
> 6  1  3  2  2  1
>
>
> On Nov 13, 2007 6:49 AM, sun <[EMAIL PROTECTED]> wrote:
>> I have two data frame A and B adn want to cross them.
>> A has format as:
>>
>> a1  a2 a3
>> 1   23
>> 2   31
>> 1   32
>> ...
>>
>> B:
>>
>> b1 b2
>> 1   2
>> 2   1
>> ...
>>
>> the combine result shall be something like
>>
>> a1 a2 a3 b1 b2
>> 1   2   3   1  2
>> 1   2   3   2  1
>> 2   3   1   1  2
>> 2   3   1   2  1
>> 1   3   2   1  2
>> 1   3   2   2  1
>> 
>>
>>
>> is there a function able of doing this instead of  loops?
>>
>> Thanks,
>> Sun
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide 
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>

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


[R] problems with splinefun()

2007-11-14 Thread david csongor
I am working with the function: splinefun() ...
When plugging in the variables, I get the function program as if
though having only entered   '"splinefun". only way to get the values
is by
spline(xxx,yyy, n=length(xxx)/10, ties = mean)$x  and  spline(xxx,yyy,
n=length(xxx)/10, ties = mean)$y.
I'm just wondering if there is something wrong with the package or if
I'm doing something wrong... Is there a way to get the values right
away? Has this happened to anyone?

Thanx in advance for the the ever great help one gets asking stuff this way!!!

/David

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


Re: [R] Newbie question about data preparation

2007-11-14 Thread Dieter Menne
Heidemeier Dr, Joachim  uba.de> writes:

> 
...
> For the analyses of the data I want to group one column (like the
> classes in a histogram).
> So I'd like to add one column with the center of each group with width=2
> for an x value in the interval of the class.
> So the output should look like 
>   x   y   x-grouped
> 1.3   2.2 1
> 2.5   3.4 3
> 3.1   3.7 3   
> 8.2   9.5 9
> 7.5   8.3 7
> What would be the R'ish idiom to do this operation?

df = data.frame(x=rnorm(10,10),y=rnorm(10,20))
df$xgrouped = (df$x+df$y) %/%2

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


[R] polygon overlap in multidimensional space

2007-11-14 Thread newRuser
Hi,

I was wondering if there is any function in R that can be used to
check for overlapping polygons in multidimensional space? I've used
overlap.xypolygon for 2D, but what if I want to check, for example, in
5D?

Any help would be greatly appreciated!
Thanks,

newRuser

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


[R] Generating log transformed random numbers

2007-11-14 Thread Tobias Schlottmann
Dear R users,   
  My question is that how it is possible to generate some random numbers using 
rnorm( ) function but in log transformed values. 
   
  Thank you,
   
  Tobias

   
-

[[alternative HTML version deleted]]

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


[R] test for existance of a method for given class

2007-11-14 Thread William Valdar
Dear All,

I want to test whether a method exists for given object. For example, 
whether a function "deviance" is defined for an object of the "lm" class.

My imperfect understanding leads me to think something like

  hasMethod("deviance", object)
  hasMethod("deviance", "lm")
  existsMethod("deviance", signature(class="lm"))

or similar might work (I don't fully understand how to manipulate 
signatures), but all the variations on this I have tried return FALSE. 
(Except, interestingly, when I first load library lme4, after which all 
return TRUE even for non-existant classes and functions).

I realize there are several ways in which R implements function 
polymorphism and that this is all documented somewhere but a hint would 
save me considerable time. I would also prefer not to resort to the hack 
solution of try()ing the function with the object and then catching the 
error to determine whether it was defined.

Thanks,

Will

=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
Dr William Valdar   ++44 (0)1865 287 589
Wellcome Trust Centre   [EMAIL PROTECTED]
for Human Genetics, Oxford  www.well.ox.ac.uk/~valdar

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


Re: [R] ESRI Shapefile for EU-25

2007-11-14 Thread Millo Giovanni
Dear Albrecht,

on ESRI's site you should be able to download a shapefile with details
down to NUTS2 for the area you're interested in.
Check http://www.esri.com/data/download/basemap/how_to_download.html
out.

Please also consider posting such questions on the R-sig-Geo mailing
list (submissions to [EMAIL PROTECTED]).

HTH 
Giovanni

Giovanni Millo
Research Dept.,
Assicurazioni Generali SpA
Via Machiavelli 4, 
34131 Trieste (Italy)
tel. +39 040 671184 
fax  +39 040 671160 

## Original message: 

Date: Tue, 13 Nov 2007 08:35:39 -0600
From: "hadley wickham" <[EMAIL PROTECTED]>
Subject: Re: [R] ESRI Shapefile for EU-25
To: "Albrecht Kauffmann" <[EMAIL PROTECTED]>
Cc: r-help@r-project.org
Message-ID:
<[EMAIL PROTECTED]>
Content-Type: text/plain; charset=ISO-8859-1

http://www.openstreetmap.org/ might be one place to start.

Hadley

On Nov 13, 2007 6:30 AM, Albrecht Kauffmann <[EMAIL PROTECTED]>
wrote:
> Hi all,
>
> who knows how to get an ESRI Shapefile for the NUTS-2 Regions of the 
> enlarged European Union? Particularly I want to draw maps of Germany, 
> Poland, Czech Republik, Hungary and Austria. I've found Shapefiles for

> the US, Russia and other countries elsewhere in the web, but for 
> Europe it seems really difficult.
>
> With many thanks for any hint
> Albrecht
>
> __
> R-help@r-project.org mailing list 
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
 
Ai sensi del D.Lgs. 196/2003 si precisa che le informazi...{{dropped:13}}

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


Re: [R] test for existance of a method for given class

2007-11-14 Thread Prof Brian Ripley
On Wed, 14 Nov 2007, William Valdar wrote:

> Dear All,
>
> I want to test whether a method exists for given object. For example,
> whether a function "deviance" is defined for an object of the "lm" class.

For an S3 generic 'f' and with an S3 object or an S3 class 'x', try

hasS3method <- function(f, x)
{
 if(is.object(x)) x <- oldClass(x)
 m <- methods(f)
 cl <- sub(paste("^", f, ".", sep=""), "", m)
 any(c("default", x) %in% cl)
}

(You can break this, e.g. by f="resid" or using implicit classes: it needs 
inside knowledge to know if the latter would be invoked.  Also, the set 
of available methods is in principle scope-specific.)

For S4 generics and classes, look at selectMethod(optional=TRUE): this is 
documented to return NULL if and only if there is no applicable method.


> My imperfect understanding leads me to think something like
>
>  hasMethod("deviance", object)
>  hasMethod("deviance", "lm")
>  existsMethod("deviance", signature(class="lm"))
>
> or similar might work (I don't fully understand how to manipulate
> signatures), but all the variations on this I have tried return FALSE.
> (Except, interestingly, when I first load library lme4, after which all
> return TRUE even for non-existant classes and functions).
>
> I realize there are several ways in which R implements function
> polymorphism and that this is all documented somewhere but a hint would
> save me considerable time. I would also prefer not to resort to the hack
> solution of try()ing the function with the object and then catching the
> error to determine whether it was defined.
>
> Thanks,
>
> Will
>
> =-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
> Dr William Valdar   ++44 (0)1865 287 589
> Wellcome Trust Centre   [EMAIL PROTECTED]
> for Human Genetics, Oxford  www.well.ox.ac.uk/~valdar

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

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


[R] Weights argument in loess (stats)

2007-11-14 Thread Felix Lamp
Hello,

I have a question concerning the weights argument of the loess
function in the stats package.

Do the weights correspond to multiplying the local regression
equation by the weights or by the square root of the weights (like in lm).

In "http://www.netlib.org/a/cloess.pdf"; page 7 it appears that the former
is true.
A look at the source code loess.f lets me conclude that the latter is   true.

Since I do not know Fortran, I would be really delighted if someone could
clarify this.

Best,

Felix

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


Re: [R] combine two dataframe

2007-11-14 Thread Gabor Grothendieck
On Nov 14, 2007 4:59 AM, sun <[EMAIL PROTECTED]> wrote:
> Thanks all for the answers. Both Merge and sqldf works perfectly for me.
> Well, I feel sqldf run a littile bit slower. And I failed to install this
> package (sqldf ) on my linux box.

Have never heard of anyone not being able to install sqldf on Linux
before.  sqldf
is written in 100% R so it should run on all platforms R runs on and for
which its dependencies work, mainly RSQLite (or RMySQL).   As mentioned on the
home page, http://sqldf.googlecode.com, sqldf is optimized for convenience,
not speed, so I would not think it would be the fastest.

>
> Denver, your approach also works, but in my case, data frame A has much more
> rows then B, so B has to be duplicated many many times.
>
> kind regards,
> Sun
>
>
> - Original Message -
> From: "Gabor Grothendieck" <[EMAIL PROTECTED]>
> To: "sun" <[EMAIL PROTECTED]>
> Cc: <[EMAIL PROTECTED]>
> Sent: Tuesday, November 13, 2007 6:07 PM
> Subject: Re: [R] combine two dataframe
>
>
> > Try this:
> >
> >> A <- data.frame(a1 = c(1, 2, 1), a2 = c(2, 3, 3), a3 = c(3, 1, 2))
> >> B <- data.frame(b1 = 1:2, b2 = 2:1)
> >>
> >> library(sqldf)
> >> sqldf("select * from A, B")
> >  a1 a2 a3 b1 b2
> > 1  1  2  3  1  2
> > 2  1  2  3  2  1
> > 3  2  3  1  1  2
> > 4  2  3  1  2  1
> > 5  1  3  2  1  2
> > 6  1  3  2  2  1
> >
> >
> > On Nov 13, 2007 6:49 AM, sun <[EMAIL PROTECTED]> wrote:
> >> I have two data frame A and B adn want to cross them.
> >> A has format as:
> >>
> >> a1  a2 a3
> >> 1   23
> >> 2   31
> >> 1   32
> >> ...
> >>
> >> B:
> >>
> >> b1 b2
> >> 1   2
> >> 2   1
> >> ...
> >>
> >> the combine result shall be something like
> >>
> >> a1 a2 a3 b1 b2
> >> 1   2   3   1  2
> >> 1   2   3   2  1
> >> 2   3   1   1  2
> >> 2   3   1   2  1
> >> 1   3   2   1  2
> >> 1   3   2   2  1
> >> 
> >>
> >>
> >> is there a function able of doing this instead of  loops?
> >>
> >> Thanks,
> >> Sun
> >>
> >> __
> >> R-help@r-project.org mailing list
> >> https://stat.ethz.ch/mailman/listinfo/r-help
> >> PLEASE do read the posting guide
> >> http://www.R-project.org/posting-guide.html
> >> and provide commented, minimal, self-contained, reproducible code.
> >>
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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


Re: [R] Generating log transformed random numbers

2007-11-14 Thread Scionforbai
If what you want is a lognormal distribution of n values you can use
the following transformations:

lognorm1 <- M*exp((rnorm(n)*sigma)-sigma^2/2.)

which gives a lognormal distribution such that:
mean(lognorm1)=M ;
var(lognorm1)=M^2*(exp(sigma^2)-1);
Changing the sigma (standard deviation) you always obtain the same
arithmetic mean.

Or, alternatively,

lognorm2 <- exp(m + sigma * rnorm(n))
such that:
exp(mean(log(lognorm2))=exp(m) [geometric mean]
mean(lognorm2)=exp(m + sigma^2/2);
var(lognorm2)=exp(2*m + sigma^2)*(exp(sigma^2/2)-1)
In this case, for different sigma values is the geometric mean to stay
constant, not the arithmetic.

Did it answer your question?

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


Re: [R] Generating log transformed random numbers

2007-11-14 Thread David Scott

If you are after the lognormal distribution, I would recommend you use the 
lognormal functions from R itself. Check them out with

?dlnorm

David Scott


On Wed, 14 Nov 2007, Scionforbai wrote:

> If what you want is a lognormal distribution of n values you can use
> the following transformations:
>
> lognorm1 <- M*exp((rnorm(n)*sigma)-sigma^2/2.)
>
> which gives a lognormal distribution such that:
> mean(lognorm1)=M ;
> var(lognorm1)=M^2*(exp(sigma^2)-1);
> Changing the sigma (standard deviation) you always obtain the same
> arithmetic mean.
>
> Or, alternatively,
>
> lognorm2 <- exp(m + sigma * rnorm(n))
> such that:
> exp(mean(log(lognorm2))=exp(m) [geometric mean]
> mean(lognorm2)=exp(m + sigma^2/2);
> var(lognorm2)=exp(2*m + sigma^2)*(exp(sigma^2/2)-1)
> In this case, for different sigma values is the geometric mean to stay
> constant, not the arithmetic.
>
> Did it answer your question?
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

_
David Scott Department of Statistics, Tamaki Campus
The University of Auckland, PB 92019
Auckland 1142,NEW ZEALAND
Phone: +64 9 373 7599 ext 86830 Fax: +64 9 373 7000
Email:  [EMAIL PROTECTED]

Graduate Officer, Department of Statistics
Director of Consulting, Department of Statistics

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


[R] Projection of Shapefiles using Maptools

2007-11-14 Thread Jens Oldeland
Dear R-Users

I have successfully imported my Shapefile using maptools with the 
following command lines, and I hoped that my object (nc) will be 
recognized as projected data, so I have included the CRS command:

library(maptools)
p4s = CRS("+proj=longlat +datum=WGS84")
nc = 
readShapePoints("D:/Projekte_INPROGRESS/Termiten/DATA/ERF_Termite.shp", 
proj4string=p4s)
plot(nc, axes = TRUE)
summary(nc)

but in the summary it tells me:

Object of class SpatialPointsDataFrame
Coordinates:
min max
coords.x1 16.93515 16.94475
coords.x2 -21.60583 -21.59686
Is projected: FALSE
proj4string : [+proj=longlat +datum=WGS84]
Number of points: 388

Is projected: FALSE ?? does this mean that it is a geographic data set 
and not a projected or does it mean that it won´t accept my proj4string ??

does anyone has an idea how to write the correct proj4string to import 
data from UTM WGS 84 S33 ? would be nice

thank you
JEns

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


Re: [R] logistic regression model specification

2007-11-14 Thread Peter Dalgaard
Dylan Beaudette wrote:
> On Tuesday 13 November 2007, Peter Dalgaard wrote:
>   
>> Prof Brian Ripley wrote:
>> 
>>> On Tue, 13 Nov 2007, Dylan Beaudette wrote:
>>>   
 Hi,

 I have setup a simple logistic regression model with the glm() function,
 with the follow formula:

 y ~ a + b

 where:
 'a' is a continuous variable stratified by
 the levels of 'b'


 Looking over the manual for model specification, it seems that
 coefficients for unordered factors are given 'against' the first level
 of that factor.
 
>>> Only for the default coding.
>>>
>>>   
 This makes for difficult interpretation when using factor 'b' as a
 stratifying model term.
 
>>> Really?  You realize that you have not 'stratified' on 'b', which would
>>> need the model to be a*b?  What you have is a model with parallel linear
>>> predictors for different levels of 'b', and if the coefficients are not
>>> telling you what you want you should change the coding.
>>>   
>
> Thanks for the comments Peter. Note that use/abuse of jargon on my part is 
> augmented by my limited understanding.
>
>   
>> I have to differ slightly here. "Stratification", at least in the fields
>> that I connect with, usually means that you combine information from
>> several groups via an assumption that they have a common value of a
>> parameter, which in the present case is essentially the same as assuming
>> an additive model y~a+b.
>> 
>
> This was the definition of 'stratification' that I was using when formulating 
> this model.
>
>   
>> I share your confusion as to why the parametrization of the effects of
>> factor b should matter, though. Surely, the original poster has already
>> noticed that the estimated effect of a is the same whether or not the
>> intercept is included? The only difference I see is that the running
>> anova() or drop1() would not give interesting results for the effect of
>> b in the no-intercept variation.
>> 
>
> Yes, I have noticed that the estimated effect is the same. It looks like I am 
> having trouble interpreting the model formula and coefficients. An example 
> from MASS, using the default R contrasts:
>
> library(MASS)
> data(whiteside)
>
> # simple additive model, relative to first level of 'Insul'
> # 'parallel regressions'
> summary(lm(Gas ~ Temp + Insul, data=whiteside))
>
> Estimate Std. Error t value Pr(>|t|)
> (Intercept)  6.551330.11809   55.48   <2e-16 ***
> Temp-0.336700.01776  -18.95   <2e-16 ***
> InsulAfter  -1.565200.09705  -16.13   <2e-16 ***
>
>
> # same model without the intercept
> summary(lm(Gas ~ Temp + Insul -1 , data=whiteside))
>
> Estimate Std. Error t value Pr(>|t|)
> Temp-0.336700.01776  -18.95   <2e-16 ***
> InsulBefore  6.551330.11809   55.48   <2e-16 ***
> InsulAfter   4.986120.10268   48.56   <2e-16 ***
>
> In presenting the terms of this model, I would like to be able to talk about 
> the physical meaning of the coefficients. Minus the intercept term, the 
> second example presents the slopes before and after. I now understand how to 
> interpret the first example, as:
>
> Intercept - 1.56520 = InsulAfter
> 6.55133 - 1.56520 = 4.98612
>
> ... the meaning of the coeficient for InsulAfter is relative to that of the 
> Intercept (treatment style contrasts).
>
> With the above example the significance terms do no really change when the 
> intercept is removed. Within the context of my data, removing the intercept 
> has the following effect:
>
> # with intercept
>Estimate Std. Error z value Pr(>|z|)
> (Intercept)   7.296e+00  1.404e+00   5.196 2.03e-07 ***
> bd_sum   -8.331e-04  1.533e-04  -5.433 5.55e-08 ***
> clastic_volcanic   -1.396e-01  8.209e-01  -0.170  0.86495
> coarse_sedimentary -2.720e+00  8.634e-01  -3.150  0.00163 ** 
> felsic_intrusive   -3.862e-01  8.404e-01  -0.460  0.64583
> fine_sedimentary   -1.010e+00  8.795e-01  -1.149  0.25069
> rhyolite   -8.420e-01  8.531e-01  -0.987  0.32365
> tuff1.569e+01  1.008e+03   0.016  0.98758
>
> # without intercept:
>Estimate Std. Error z value Pr(>|z|)
> bd_sum   -8.331e-04  1.533e-04  -5.433 5.55e-08 ***
> andesite7.296e+00  1.404e+00   5.196 2.03e-07 ***
> clastic_volcanic7.156e+00  1.276e+00   5.608 2.05e-08 ***
> coarse_sedimentary  4.576e+00  1.158e+00   3.951 7.77e-05 ***
> felsic_intrusive6.910e+00  1.252e+00   5.520 3.39e-08 ***
> fine_sedimentary6.286e+00  1.279e+00   4.916 8.83e-07 ***
> rhyolite6.454e+00  1.289e+00   5.005 5.57e-07 ***
> tuff2.299e+01  1.008e+03   0.0230.982
>
>
> ... the meaning of the coefficients now makes sense (I think!), however the 
> significance of each level of the 'stratifying' factor has changed 
> c

[R] geotiff calculations

2007-11-14 Thread Monica Pisica

Laura,
 
As far as i know RGDAL and maybe proj4 will help you with geotiff, but it may 
be tricky if you want to preserve the geotiff projection. Your safest bet is to 
use a GIS software and ESRI ArcGIS can do everything you want without using R. 
To calculate differences, mean, and so on is a simple task in ArcGIS, no need 
for R. Evem more advanced statistics can be done in ArcGIS in their "expanded" 
calculator. You also can do some geostatistics, analysis of variances, 
correlations, "hot spots" and local outliers.
 
I hope this helps,
 
Monica
 
-- Message: 110Date: Wed, 14 Nov 2007 09:22:52 
+0100From: "Laura Poggio" <[EMAIL PROTECTED]>Subject: [R] geotiff 
calculationsTo: [EMAIL PROTECTED]:<[EMAIL PROTECTED]>Content-Type: text/plain 
Dear list,I have to compare two digital elevation models in raster format 
(geotiff).I then have to calculate the differences in altitude for each cell 
and makesome statistics (basic as mean, median, std, range but also more 
"advanced"as RMSE) on that.I do not know very much how to proceed:1) is it 
possible to import the geotiff in R? If so with which package? ifnot which is 
the best way to import such files?2) is it better to perform the calculations 
of the differences in a GISsoftware and then to use R only for statistical 
analysis? or it is better todo everything in R?3) is there any specific package 
for doing this kind of analysis? Thank you very much in advance Laura 
[[alternative HTML version del!
 eted]]   --
_
[[replacing trailing spam]]

[[alternative HTML version deleted]]

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


Re: [R] test for existance of a method for given class

2007-11-14 Thread William Valdar
Thanks Brian, that helps a lot. For others interested, a few variants for 
testing existance of methods are below:

# modified sub/grep of BDR's example
hasS3method.1 <- function(f, x)
{
 if(is.object(x)) x <- oldClass(x)
 m <- methods(f)
 pattern <- paste("^", f, ".", sep="")
 cl <- sub(pattern, "", grep(pattern, m, value=TRUE))
 any(c("default", x) %in% cl)
}

# almost equivalently...
hasS3method.2 <- function(f, x, include.default=TRUE)
{
 if(is.object(x)) x <- oldClass(x)
 !is.null(getS3method(f, x, optional=TRUE))
}

hasS4method <- function(f, x)
{
 if (is.object(x)) x <- class(x)
 for (cl in x)
 {
 m <- selectMethod(f, signature(object=cl), optional=TRUE)
 if (!is.null(m)) return (TRUE)
 }
 FALSE
}

Will

On Wed, 14 Nov 2007, Prof Brian Ripley wrote:
> On Wed, 14 Nov 2007, William Valdar wrote:
>
>> Dear All,
>> 
>> I want to test whether a method exists for given object. For example,
>> whether a function "deviance" is defined for an object of the "lm" class.
>
> For an S3 generic 'f' and with an S3 object or an S3 class 'x', try
>
> hasS3method <- function(f, x)
> {
>if(is.object(x)) x <- oldClass(x)
>m <- methods(f)
>cl <- sub(paste("^", f, ".", sep=""), "", m)
>any(c("default", x) %in% cl)
> }
>
> (You can break this, e.g. by f="resid" or using implicit classes: it needs 
> inside knowledge to know if the latter would be invoked.  Also, the set of 
> available methods is in principle scope-specific.)
>
> For S4 generics and classes, look at selectMethod(optional=TRUE): this is 
> documented to return NULL if and only if there is no applicable method.
>
>
>> My imperfect understanding leads me to think something like
>>
>>  hasMethod("deviance", object)
>>  hasMethod("deviance", "lm")
>>  existsMethod("deviance", signature(class="lm"))
>> 
>> or similar might work (I don't fully understand how to manipulate
>> signatures), but all the variations on this I have tried return FALSE.
>> (Except, interestingly, when I first load library lme4, after which all
>> return TRUE even for non-existant classes and functions).
>> 
>> I realize there are several ways in which R implements function
>> polymorphism and that this is all documented somewhere but a hint would
>> save me considerable time. I would also prefer not to resort to the hack
>> solution of try()ing the function with the object and then catching the
>> error to determine whether it was defined.
>> 
>> Thanks,
>> 
>> Will
>>

=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
Dr William Valdar   ++44 (0)1865 287 589
Wellcome Trust Centre   [EMAIL PROTECTED]
for Human Genetics, Oxford  www.well.ox.ac.uk/~valdar

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


[R] label plotting on nmds diagram

2007-11-14 Thread Antonio Olinto
Hi,

I'm using nmds command (library vegan) to analyze some fishing data.

I'd like to plot not only points, but also the names of species and stations in
a specified position.

I used the command
text(nmds$points[,1], nmds $points[,2],labels=row.names(nmds
$points),pos=3,cex=0.5)

But the labels are sometimes overlapped.

Is there any way to use identify, or a similar command, to plot the row names in
a given position? Identify would be perfect but it indicates the row number and
I'd like to have the row name given in nmds$points.

Many thanks,

Antonio Olinto
Sao Paulo Fisheries Institute
Brasil






-
BCMG Internet Webmail
www.bcmg.com.br

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


Re: [R] problems with splinefun()

2007-11-14 Thread jim holtman
Exactly what values do you want "right away"?

You can do:

result <- spline(.)

and then reference 'result$x' and 'result$y'.  Can you be more
specific on your request and provide an example of what you are
currently doing (with data) and what you expect the results to be.

On Nov 14, 2007 5:31 AM, david csongor <[EMAIL PROTECTED]> wrote:
> I am working with the function: splinefun() ...
> When plugging in the variables, I get the function program as if
> though having only entered   '"splinefun". only way to get the values
> is by
> spline(xxx,yyy, n=length(xxx)/10, ties = mean)$x  and  spline(xxx,yyy,
> n=length(xxx)/10, ties = mean)$y.
> I'm just wondering if there is something wrong with the package or if
> I'm doing something wrong... Is there a way to get the values right
> away? Has this happened to anyone?
>
> Thanx in advance for the the ever great help one gets asking stuff this way!!!
>
> /David
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



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

What is the problem you are trying to solve?

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


Re: [R] label plotting on nmds diagram

2007-11-14 Thread Gavin Simpson
On Wed, 2007-11-14 at 11:22 -0300, Antonio Olinto wrote:
> Hi,
> 
> I'm using nmds command (library vegan) to analyze some fishing data.

I doubt it - there is no nmds command. Do you mean metaMDS() ? This uses
isoMDS() from package MASS with some extra features.

> 
> I'd like to plot not only points, but also the names of species and stations 
> in
> a specified position.
> 
> I used the command
> text(nmds$points[,1], nmds $points[,2],labels=row.names(nmds
> $points),pos=3,cex=0.5)
> 
> But the labels are sometimes overlapped.
> 
> Is there any way to use identify, or a similar command, to plot the row names 
> in
> a given position? Identify would be perfect but it indicates the row number 
> and
> I'd like to have the row name given in nmds$points.

Have a look at ?orditorp in vegan for an alternative approach, but yes
you can use identify, but you have to capture the output of the plot
call and use that in your identify call:

require(vegan)
data(varespec)
mod <- cca(varespec)
mod.plt <- plot(mod, display = "species", type = "p")
identify(mod.plt)

See ?identify.ordiplot for more info.

HTH

> 
> Many thanks,
> 
> Antonio Olinto
> Sao Paulo Fisheries Institute
> Brasil
> 
> 
> 
> 
> 
> 
> -
> BCMG Internet Webmail
> www.bcmg.com.br
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
-- 
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
 Dr. Gavin Simpson [t] +44 (0)20 7679 0522
 ECRC, UCL Geography,  [f] +44 (0)20 7679 0565
 Pearson Building, [e] gavin.simpsonATNOSPAMucl.ac.uk
 Gower Street, London  [w] http://www.ucl.ac.uk/~ucfagls/
 UK. WC1E 6BT. [w] http://www.freshwaters.org.uk
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%

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


Re: [R] time plotting problem

2007-11-14 Thread John Kane
Excellent, that looks very nice!

 I had not realised that plotrix would do that and I
was definately hesitating to try and figure out how to
do it myself. 

Any idea why all I was getting was a single date on
the x-axis when just doing plot?  


--- Jim Lemon <[EMAIL PROTECTED]> wrote:

> John Kane wrote:
> > I clearly spoke too soon.  
> > 
> > With the actual data I am not getting sensible
> x-axis
> > units.  The program with the actual data below. 
> Graph
> > output is here:
> > http://ca.geocities.com/jrkrideau/R/hd.png .  
> > 
> > I seem to be getting only a single entry for the
> > x-axis of "2007". However dates range from
> >  First DateLast Date 
> > "2006-09-26" "2007-11-10" 
> > 
> > I must be missing something blindingly obvious but
> I
> > don't see it. 
> > 
> > Thanks for any suggestions.
> > 
> > 
> > Actual data and test program
> > 
> > mydata <-
> >
>
read.table("http://ca.geocities.com/jrkrideau/R/heartdata.txt";,
> > sep="\t",
> >  header=FALSE)
> > mydata[,1] <- as.Date(mydata[,1],"%m/%d/%y")
> > names(mydata) <- Cs(dates, sy,dys,pulse, weight)
> > minmax  <-  c(max(mydata[,2]),min(mydata[,2]),
> > max(mydata[,3]),min(mydata[,3]))
> > names(minmax) <- c("maxsys", "minsys", "maxdsy",
> > "mindys") ; minmax 
> > 
> > min.max.dates <- c(min(mydata[,1]),
> max(mydata[,1]))
> > names(min.max.dates) <- c("First Date", "Last
> Date");
> > min.max.dates
> > 
> > ss <- lm(mydata[,2]~mydata[,1])
> > plot(mydata[,1],mydata[,2], xlab="Dates",
> ylab="Blood
> > pressure",  
> >ylim=c(minmax[4], minmax[1]), col= "red",
> > type="l")
> > abline (ss, col ="yellow")
> > dd <- lm(mydata[,3]~mydata[,1])
> > points(mydata[,1], mydata[,3], type="l",
> col="blue" )
> > abline(dd, col= "yellow", lwd=2)
> > 
> > 
> > 
> Hi John,
> 
> It can be done like this:
> 
> plot(mydata[,1],mydata[,2], xlab="Dates",
> ylab="Blood
>   pressure",ylim=c(minmax[4], minmax[1]), col=
> "red",
>   type="l",xaxt="n")
> ...
> axisdates<-c("2006-11-01","2007-01-01","2007-03-01",
>  
> "2007-05-01","2007-07-01","2007-09-01","2007-11-01")
> staxlab(at=as.Date(axisdates),labels=axisdates)
> 
> Jim
>

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


[R] convex optimization package for R, specifically semidefinite programming

2007-11-14 Thread Galkowski, Jan
Recently, a package for convex optimization was announced for Python,
based upon the LP solver GLPK, the SDP solver
in DSDP5, and the LP and QP solvers in MOSEK.  I'm aware GLPK is
available for R, but wondered if anyone had good
packages for convex optimization along these lines for R.

TIA.

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


[R] as.data.frame.array

2007-11-14 Thread Wollkind, Steven
I've noticed what seems like inconsistent behavior from this function,
but I can't figure out if this is intentional or not.

If an array has only one element, the names attached to that array are
not transferred to the data frame that is returned, while If the array
has two or more elements the names are preserved.  The documentation for
as.data.frame indicates that if the row.names argument is NULL then the
names attached to the object in question will be used for row names of
the resulting data frame.

Some sample code to demonstrate this behavior:

a <- array(c(1))
names(a) <- c("A")

b <- array(c(1,2))
names(b) <- c("A","B")

as.data.frame(a)
as.data.frame(b)


The code for as.data.frame.array reveals that drop is being called on
the array and then it is passed on to as.data.frame.vector.  It looks
like drop strips the names off of the array that is being passed, so it
seems like we should be checking if the row.names argument is null and
fetching those names before casting them aside.

Here's the code for as.data.frame.array:

function (x, row.names = NULL, optional = FALSE, ...) 
{
d <- dim(x)
if (length(d) == 1L) {
value <- as.data.frame.vector(drop(x), row.names, optional, 
...)
if (!optional) 
names(value) <- deparse(substitute(x))[[1L]]
value
}
else if (length(d) == 2L) {
as.data.frame.matrix(x, row.names, optional, ...)
}
else {
dn <- dimnames(x)
dim(x) <- c(d[1L], prod(d[-1]))
if (!is.null(dn)) {
if (length(dn[[1L]])) 
rownames(x) <- dn[[1L]]
for (i in 2L:length(d)) if (is.null(dn[[i]])) 
dn[[i]] <- seq_len(d[i])
colnames(x) <- interaction(expand.grid(dn[-1]))
}
as.data.frame.matrix(x, row.names, optional, ...)
}
}


Can someone tell me if there's a good reason for this behavior, or if
this is a bug in as.data.frame.array?  Also, is this the sort of
question that belongs on R devel rather than R help?  I'm only signed on
this list at present.

Thanks
Steve


Steve Wollkind
Associate Analyst
Geode Capital Management, LLC
1 Post Office Square / 28th Floor / Boston, MA 02109
[EMAIL PROTECTED]
Tel:   (617) 392-8991
Fax:  (617) 476-6389

This e-mail, and any attachments hereto, are intended fo...{{dropped:11}}

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


Re: [R] cronbach's alpha

2007-11-14 Thread John Kane
alpha(psychometric) ?  
--- raymond chiruka <[EMAIL PROTECTED]> wrote:

> hie 
> 1...i'm trying to carryout a relibility
> testusing cronbach's alpha what fuctin do i use.
> 
> 2.. this is more of a statistical question.if
> the alpha value  for all the variables  is negative
> what does it mean. and if the  alpha value is
> negative for all tyha variables but is greater than
> 0.7  for some sections of the variables what does
> that mean
>   
>   thanks in advance

Could be several things but you may just have a coding
problem.  Is there any chance you have reversed your
coding of some items?

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


[R] no applicable method for "names"

2007-11-14 Thread Schiller Judith 1541 EB
hi,

after installing R-2.6.0 the function "names" doesn't work anymore on my
windows xp machine. 
for example for a simple vector i get

> z <- 1:3
> names(x)
Error in UseMethod("name"): no applicable method for "names"

... instead of NULL. the same is true for lists and dataframes. attr(z,
"names") is a workaround, but i don't want to change all my functions.

is this a know bug? thanks for any comments, 
judith schiller 

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


Re: [R] convex optimization package for R, specifically semidefinite programming

2007-11-14 Thread roger koenker
Convex optimization is a large problem domain with many interesting
new developments,  and MOSEK is a particularly good implementation
of many of these developments.  It is quite easy to get R to speak to
MOSEK via R.matlab -- i've been doing this recently for some work on  
density
estimation, but given the commercial nature of MOSEK and most of the
other convex optimization development efforts closer coordination seems
unlikely.  (I would be delighted to be proved wrong about this...)


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 Nov 14, 2007, at 9:31 AM, Galkowski, Jan wrote:

> Recently, a package for convex optimization was announced for Python,
> based upon the LP solver GLPK, the SDP solver
> in DSDP5, and the LP and QP solvers in MOSEK.  I'm aware GLPK is
> available for R, but wondered if anyone had good
> packages for convex optimization along these lines for R.
>
> TIA.
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting- 
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

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


[R] forecasting package installation errors

2007-11-14 Thread Dan Robb
R gurus,

I've exhausted my search of online help.  This is my last resort.

We're running R-2.1.1.  I've been told by one of the R users here that 
we cannot upgrade to the lastest version because of some python rpy 
wrapper dependency that hasn't caught up to the latest version of R.  
The software is running on Solaris 10 x86.  Gcc version is 3.4.1.  R is 
installed in /opt/app/R

 > Sys.getlocale()
[1] "/iso_8859_1/C/C/C/C/C"

 > install.packages("forecasting",dependencies=TRUE)
Warning in install.packages("forecasting", dependencies = TRUE) :
argument 'lib' is missing: using /opt/app/R/lib/R/library
trying URL
'http://www.biometrics.mtu.edu/CRAN/src/contrib/forecasting_1.07.tar.gz'
Content type 'application/x-gzip' length 1059921 bytes
opened URL
==
downloaded 1035Kb

* Installing *source* package 'forecast' ...
** libs
gcc -I/opt/app/R-2.1.1/lib/R/include -I/opt/app/include
-D__NO_MATH_INLINES -fPIC -I/opt/app/include -c etscalc.c -o etscalc.o
etscalc.c:8:1: warning: "HUGE" redefined
In file included from etscalc.c:1:
/opt/app/g++lib6/gcc-3.4/lib/gcc/i386-pc-solaris2.10/3.4.1/include/math.h:118:1:
warning: this is the location of the previous definition
gcc -G -L/opt/app/lib -R/opt/app/lib -o forecast.so etscalc.o
-L/opt/app/R-2.1.1/lib/R/lib -lR
** R
** data
** moving datasets to lazyload DB
** preparing package for lazy loading
Error in importIntoEnv(impenv, impnames, ns, impvars) :
object 'simulate' is not exported by 'namespace:stats'
Execution halted
ERROR: lazy loading failed for package 'forecast'
** Removing '/opt/app/R-2.1.1/lib/R/library/forecast'
** Removing '/opt/app/R-2.1.1/lib/R/library/fma'
** Removing '/opt/app/R-2.1.1/lib/R/library/Mcomp'

The downloaded packages are in
/tmp/RtmplNa4PN/downloaded_packages
Warning message:
installation of package 'forecasting' had non-zero exit status in:
install.packages("forecasting", dependencies = TRUE)

Any help would be greatly appreciated.

-Dan

-- 
---
Daniel P. Robb
71 S. Wacker
Suite 1900
Chicago, IL 60606
[EMAIL PROTECTED]
direct: (312) 264-2080
fax: (312) 264-2001

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


Re: [R] time plotting problem

2007-11-14 Thread John Kane
Thank you. Almost too fancy for the target user but
very nice indeed.  For some reason I had not thought
of using ggplot2.  


--- hadley wickham <[EMAIL PROTECTED]> wrote:

> Here's a version using reshape and ggplot:
> 
> mydata <-
>
read.table("http://ca.geocities.com/jrkrideau/R/heartdata.txt";,
> sep="\t", header=FALSE)
> 
> mydata[,1] <- as.Date(mydata[,1],"%m/%d/%y")
> names(mydata) <- c("dates", "sy","dys","pulse",
> "weight")
> 
> molten <- melt(mydata, m = c("sy", "dys"))
> qplot(dates, value, data=molten, geom="line",
> colour=variable) +
> geom_smooth(method=lm)
> 
> 
> ggplot has somewhat better defaults for dates in
> this case.  See
> http://had.co.nz/ggplot2/scale_date.html for how to
> tweak the
> defaults.
> 
> Hadley
> 
> -- 
> http://had.co.nz/
>

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


[R] How to get row numbers of a subset of rows

2007-11-14 Thread affy snp
Hello list,

I read in a txt file using

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


Re: [R] time plotting problem

2007-11-14 Thread John Kane
Thank you Jim.  That does give better results.  I had
not realised how complicated a question I was asking. 


Both yours and Jim Lemon's solutions work very nicely.
  I am still messing up the syntax with Gabour's
approach. Thanks to all for the fast and valuable
help.

--- jim holtman <[EMAIL PROTECTED]> wrote:

> I have had some problems with dates on axis in
> certain ranges.  Here
> is a version of axis.POSIXct that I have been using
> that seems to do a
> better job:
> 
>  my.axis.POSIXct <-
> function (side, x, format, inc)
> {
> .diff <- diff(range(unclass(x), na.rm = TRUE))
> if (.diff < 30) {
> .by <- "sec"
> .for <- "%M:%S"
> }
> else if (.diff < 5 * 60) {
> .by <- "15 secs"
> .for <- "%M:%S"
> }
> else if (.diff < 30 * 60) {
> .by <- "2 min"
> .for <- "%H:%M"
> }
> else if (.diff < 2 * 3600) {
> .by <- "15 mins"
> .for <- "%H:%M"
> }
> else if (.diff < 12 * 3600) {
> .by <- "hour"
> .for <- "%H:00"
> }
> else if (.diff < 48 * 3600) {
> .by <- "6 hours"
> .for <- "%d-%H"
> }
> else if (.diff < 7 * 24 * 3600) {
> .by <- "12 hours"
> .for <- "%d-%H"
> }
> else if (.diff < 14 * 24 * 3600) {
> .by <- "day"
> .for <- "%m/%d"
> }
> else {
> .inc <- diff(pretty(x))[1]
> if (.inc < 7 * 84600) {
> .by <- paste(as.integer(.inc/86400),
> "days")
> }
> else {
> .by <- paste(as.integer(.inc/(86400 *
> 7)) * 7, "days")
> }
> .for <- "%m/%d/%y"
> }
> .range <- range(x, na.rm = TRUE)
> axis.POSIXct(side, at = seq(trunc(.range[1],
> "day"), trunc(.range[2] +
> 86400, "day"), by = ifelse(missing(inc),
> .by, inc)),
> format = ifelse(missing(format), .for,
> format))
> }
> 
> With your code I did the following which may give
> you better results:
> 
> plot(as.POSIXct(mydata[,1]),mydata[,2],
> xlab="Dates", ylab="Blood
> pressure",
>   ylim=c(minmax[4], minmax[1]), col= "red",
> type="l", xaxt='n')
> my.axis.POSIXct(1, as.POSIXct(mydata[,1]))
> 
> 
> On Nov 13, 2007 9:08 AM, John Kane
> <[EMAIL PROTECTED]> wrote:
> > I clearly spoke too soon.
> >
> > With the actual data I am not getting sensible
> x-axis
> > units.  The program with the actual data below. 
> Graph
> > output is here:
> > http://ca.geocities.com/jrkrideau/R/hd.png .
> >
> > I seem to be getting only a single entry for the
> > x-axis of "2007". However dates range from
> >  First DateLast Date
> > "2006-09-26" "2007-11-10"
> >
> > I must be missing something blindingly obvious but
> I
> > don't see it.
> >
> > Thanks for any suggestions.
> >
> >
> > Actual data and test program
> > 
> > mydata <-
> >
>
read.table("http://ca.geocities.com/jrkrideau/R/heartdata.txt";,
> > sep="\t",
> > header=FALSE)
> > mydata[,1] <- as.Date(mydata[,1],"%m/%d/%y")
> > names(mydata) <- Cs(dates, sy,dys,pulse, weight)
> > minmax  <-  c(max(mydata[,2]),min(mydata[,2]),
> > max(mydata[,3]),min(mydata[,3]))
> > names(minmax) <- c("maxsys", "minsys", "maxdsy",
> > "mindys") ; minmax
> >
> > min.max.dates <- c(min(mydata[,1]),
> max(mydata[,1]))
> > names(min.max.dates) <- c("First Date", "Last
> Date");
> > min.max.dates
> >
> > ss <- lm(mydata[,2]~mydata[,1])
> > plot(mydata[,1],mydata[,2], xlab="Dates",
> ylab="Blood
> > pressure",
> >   ylim=c(minmax[4], minmax[1]), col= "red",
> > type="l")
> > abline (ss, col ="yellow")
> > dd <- lm(mydata[,3]~mydata[,1])
> > points(mydata[,1], mydata[,3], type="l",
> col="blue" )
> > abline(dd, col= "yellow", lwd=2)
> >
> >
> >
> >
> >
> >
> > > Thanks to Gabor and Jim. I am not sure if the
> first
> > > entry year = 2009 is all the problem I'm getting
> but
> > > it is certainly seems like the worst of it.
> > >
> > > My stupidity:  Someone sent me the data set in
> Excel
> > >
> > > and I didn't do the basic data checks on. I
> _KNEW_
> > > the
> > > data went from 2006 to 2007.
> > >
> > > --- Gabor Grothendieck <[EMAIL PROTECTED]>
> > > wrote:
> > >
> > > > In your examples the first line of your data
> > > refers
> > > > to the
> > > > year 2009 and Oct 1st is repeated.  Is that
> really
> > > > what
> > > > you meant?
> > > >
> > > > I can't tell what your problem is from your
> > > > description
> > > > other than the data problems cited but there
> are
> > > > lots of
> > > > examples of plotting with zoo in the following
> > > which
> > > > may
> > > > help you:
> > > > vignette("zoo")
> > > > vignette("zoo-quickref")
> > > > ?plot.zoo
> > > > ?xyplot.zoo
> > > >
> > > > Note that zoo series must be time series, i.e.
> > > they
> > > > must
> > > > have unique times.
> > > >
> > > > On Nov 12, 2007 1:47 PM, John Kane
> >
> > > > > I am completely misunderstanding how to
> handle
> > > > dates.
> > > > > I want to plot a couple of data series
> against
> > > > s

Re: [R] How to get row numbers of a subset of rows

2007-11-14 Thread jim holtman
Here is a way of doing it using 'rle':

> x <- read.table(textConnection(" SNPChromosome  
> PhysicalPosition
+ 1 SNP_A-1909444  1   7924293
+ 2 SNP_A-2237149  1   8173763
+ 3 SNP_A-4303947  1   8191853
+ 4 SNP_A-2236359  1   8323433
+ 5 SNP_A-2205441  1   8393263
+ 6 SNP_A-1909445  1   7924293
+ 7 SNP_A-2237146  2   8173763
+ 8 SNP_A-4303946  2   8191853
+ 9 SNP_A-2236357  2   8323433
+ 10 SNP_A-2205442 2   8393263"), header=TRUE)
> # use rle to get the 'runs'
> y <- rle(x$Chromosome)
> # create dataframe with start/ends and values
> start <- head(cumsum(c(1, y$lengths)), -1)
> index <- data.frame(values=y$values, start=start, end=start + y$lengths - 1)
>
> index
  values start end
1  1 1   6
2  2 7  10
>


On Nov 14, 2007 10:56 AM, affy snp <[EMAIL PROTECTED]> wrote:
> Hello list,
>
> I read in a txt file using
>
> 
> by specifying the row.names=NULL so that the rows are numbered.
> Below is an example after how the table looks like using
> 
>
>  SNPChromosome  PhysicalPosition
> 1 SNP_A-1909444  1   7924293
> 2 SNP_A-2237149  1   8173763
> 3 SNP_A-4303947  1   8191853
> 4 SNP_A-2236359  1   8323433
> 5 SNP_A-2205441  1   8393263
> 6 SNP_A-1909445  1   7924293
> 7 SNP_A-2237146  2   8173763
> 8 SNP_A-4303946  2   8191853
> 9 SNP_A-2236357  2   8323433
> 10 SNP_A-2205442 2   8393263
>
> I am wondering if there is a way to return the start and end row numbers
> for a subset of rows.
>
> For example, If I specify B[,2]=1, I would like to get
> start=1 and end=6
>
> if B[,2]=2, then start=7 and end=10
>
> Is there any way in R to quickly do this?
>
> Thanks a bunch!
>
> Allen
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



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

What is the problem you are trying to solve?

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


Re: [R] How to get row numbers of a subset of rows

2007-11-14 Thread Bert Gunter
Am I missing something? ...

Why not: range(seq(nrow(B))[B[,2]==1] ) ?? ## note: "==" not "=" 

Alternatively, and easily generalized (to start with a frame which is a
subset of the original and any subset of rows, contiguous or not)

range(as.numeric(row.names(B)[B[,2]==1]))

Again, am I missing something that makes this "obvious" solution impossible?
(Wouldn't be the first time.)

Bert Gunter
Genentech Nonclinical Statistics


-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On
Behalf Of jim holtman
Sent: Wednesday, November 14, 2007 8:39 AM
To: affy snp
Cc: r-help@r-project.org
Subject: Re: [R] How to get row numbers of a subset of rows

Here is a way of doing it using 'rle':

> x <- read.table(textConnection(" SNPChromosome
PhysicalPosition
+ 1 SNP_A-1909444  1   7924293
+ 2 SNP_A-2237149  1   8173763
+ 3 SNP_A-4303947  1   8191853
+ 4 SNP_A-2236359  1   8323433
+ 5 SNP_A-2205441  1   8393263
+ 6 SNP_A-1909445  1   7924293
+ 7 SNP_A-2237146  2   8173763
+ 8 SNP_A-4303946  2   8191853
+ 9 SNP_A-2236357  2   8323433
+ 10 SNP_A-2205442 2   8393263"), header=TRUE)
> # use rle to get the 'runs'
> y <- rle(x$Chromosome)
> # create dataframe with start/ends and values
> start <- head(cumsum(c(1, y$lengths)), -1)
> index <- data.frame(values=y$values, start=start, end=start + y$lengths -
1)
>
> index
  values start end
1  1 1   6
2  2 7  10
>


On Nov 14, 2007 10:56 AM, affy snp <[EMAIL PROTECTED]> wrote:
> Hello list,
>
> I read in a txt file using
>
> 
> by specifying the row.names=NULL so that the rows are numbered.
> Below is an example after how the table looks like using
> 
>
>  SNPChromosome  PhysicalPosition
> 1 SNP_A-1909444  1   7924293
> 2 SNP_A-2237149  1   8173763
> 3 SNP_A-4303947  1   8191853
> 4 SNP_A-2236359  1   8323433
> 5 SNP_A-2205441  1   8393263
> 6 SNP_A-1909445  1   7924293
> 7 SNP_A-2237146  2   8173763
> 8 SNP_A-4303946  2   8191853
> 9 SNP_A-2236357  2   8323433
> 10 SNP_A-2205442 2   8393263
>
> I am wondering if there is a way to return the start and end row numbers
> for a subset of rows.
>
> For example, If I specify B[,2]=1, I would like to get
> start=1 and end=6
>
> if B[,2]=2, then start=7 and end=10
>
> Is there any way in R to quickly do this?
>
> Thanks a bunch!
>
> Allen
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



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

What is the problem you are trying to solve?

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

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


Re: [R] logistic regression model specification

2007-11-14 Thread Dylan Beaudette
On Nov 14, 2007 5:17 AM, Peter Dalgaard <[EMAIL PROTECTED]> wrote:
>
> Dylan Beaudette wrote:
> > On Tuesday 13 November 2007, Peter Dalgaard wrote:
> >
> >> Prof Brian Ripley wrote:
> >>
> >>> On Tue, 13 Nov 2007, Dylan Beaudette wrote:
> >>>
>  Hi,
> 
>  I have setup a simple logistic regression model with the glm() function,
>  with the follow formula:
> 
>  y ~ a + b
> 
>  where:
>  'a' is a continuous variable stratified by
>  the levels of 'b'
> 
> 
>  Looking over the manual for model specification, it seems that
>  coefficients for unordered factors are given 'against' the first level
>  of that factor.
> 
> >>> Only for the default coding.
> >>>
> >>>
>  This makes for difficult interpretation when using factor 'b' as a
>  stratifying model term.
> 
> >>> Really?  You realize that you have not 'stratified' on 'b', which would
> >>> need the model to be a*b?  What you have is a model with parallel linear
> >>> predictors for different levels of 'b', and if the coefficients are not
> >>> telling you what you want you should change the coding.
> >>>
> >
> > Thanks for the comments Peter. Note that use/abuse of jargon on my part is
> > augmented by my limited understanding.
> >
> >
> >> I have to differ slightly here. "Stratification", at least in the fields
> >> that I connect with, usually means that you combine information from
> >> several groups via an assumption that they have a common value of a
> >> parameter, which in the present case is essentially the same as assuming
> >> an additive model y~a+b.
> >>
> >
> > This was the definition of 'stratification' that I was using when 
> > formulating
> > this model.
> >
> >
> >> I share your confusion as to why the parametrization of the effects of
> >> factor b should matter, though. Surely, the original poster has already
> >> noticed that the estimated effect of a is the same whether or not the
> >> intercept is included? The only difference I see is that the running
> >> anova() or drop1() would not give interesting results for the effect of
> >> b in the no-intercept variation.
> >>
> >
> > Yes, I have noticed that the estimated effect is the same. It looks like I 
> > am
> > having trouble interpreting the model formula and coefficients. An example
> > from MASS, using the default R contrasts:
> >
> > library(MASS)
> > data(whiteside)
> >
> > # simple additive model, relative to first level of 'Insul'
> > # 'parallel regressions'
> > summary(lm(Gas ~ Temp + Insul, data=whiteside))
> >
> > Estimate Std. Error t value Pr(>|t|)
> > (Intercept)  6.551330.11809   55.48   <2e-16 ***
> > Temp-0.336700.01776  -18.95   <2e-16 ***
> > InsulAfter  -1.565200.09705  -16.13   <2e-16 ***
> >
> >
> > # same model without the intercept
> > summary(lm(Gas ~ Temp + Insul -1 , data=whiteside))
> >
> > Estimate Std. Error t value Pr(>|t|)
> > Temp-0.336700.01776  -18.95   <2e-16 ***
> > InsulBefore  6.551330.11809   55.48   <2e-16 ***
> > InsulAfter   4.986120.10268   48.56   <2e-16 ***
> >
> > In presenting the terms of this model, I would like to be able to talk about
> > the physical meaning of the coefficients. Minus the intercept term, the
> > second example presents the slopes before and after. I now understand how to
> > interpret the first example, as:
> >
> > Intercept - 1.56520 = InsulAfter
> > 6.55133 - 1.56520 = 4.98612
> >
> > ... the meaning of the coeficient for InsulAfter is relative to that of the
> > Intercept (treatment style contrasts).
> >
> > With the above example the significance terms do no really change when the
> > intercept is removed. Within the context of my data, removing the intercept
> > has the following effect:
> >
> > # with intercept
> >Estimate Std. Error z value Pr(>|z|)
> > (Intercept)   7.296e+00  1.404e+00   5.196 2.03e-07 ***
> > bd_sum   -8.331e-04  1.533e-04  -5.433 5.55e-08 ***
> > clastic_volcanic   -1.396e-01  8.209e-01  -0.170  0.86495
> > coarse_sedimentary -2.720e+00  8.634e-01  -3.150  0.00163 **
> > felsic_intrusive   -3.862e-01  8.404e-01  -0.460  0.64583
> > fine_sedimentary   -1.010e+00  8.795e-01  -1.149  0.25069
> > rhyolite   -8.420e-01  8.531e-01  -0.987  0.32365
> > tuff1.569e+01  1.008e+03   0.016  0.98758
> >
> > # without intercept:
> >Estimate Std. Error z value Pr(>|z|)
> > bd_sum   -8.331e-04  1.533e-04  -5.433 5.55e-08 ***
> > andesite7.296e+00  1.404e+00   5.196 2.03e-07 ***
> > clastic_volcanic7.156e+00  1.276e+00   5.608 2.05e-08 ***
> > coarse_sedimentary  4.576e+00  1.158e+00   3.951 7.77e-05 ***
> > felsic_intrusive6.910e+00  1.252e+00   5.520 3.39e-08 ***
> > fine_sedimentary6.286e+00  1.279e+00   4.916 8.83e-07 ***
> > rhyolite6.454e+00  1.289e+00   5.005 5.57e-07 ***
> > tuf

Re: [R] How to get row numbers of a subset of rows

2007-11-14 Thread jim holtman
That works for the specific value of '1', but you would have to repeat
it for other values in the column.  If you had 100 different ranges in
that column, what would you do?  Here is another solution using
'range' on the same data:

> tapply(seq_len(nrow(x)), x$Chromosome, range)
$`1`
[1] 1 6

$`2`
[1]  7 10


On Nov 14, 2007 12:04 PM, Bert Gunter <[EMAIL PROTECTED]> wrote:
> Am I missing something? ...
>
> Why not: range(seq(nrow(B))[B[,2]==1] ) ?? ## note: "==" not "="
>
> Alternatively, and easily generalized (to start with a frame which is a
> subset of the original and any subset of rows, contiguous or not)
>
> range(as.numeric(row.names(B)[B[,2]==1]))
>
> Again, am I missing something that makes this "obvious" solution impossible?
> (Wouldn't be the first time.)
>
> Bert Gunter
> Genentech Nonclinical Statistics
>
>
>
> -Original Message-
> From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On
> Behalf Of jim holtman
> Sent: Wednesday, November 14, 2007 8:39 AM
> To: affy snp
> Cc: r-help@r-project.org
> Subject: Re: [R] How to get row numbers of a subset of rows
>
> Here is a way of doing it using 'rle':
>
> > x <- read.table(textConnection(" SNPChromosome
> PhysicalPosition
> + 1 SNP_A-1909444  1   7924293
> + 2 SNP_A-2237149  1   8173763
> + 3 SNP_A-4303947  1   8191853
> + 4 SNP_A-2236359  1   8323433
> + 5 SNP_A-2205441  1   8393263
> + 6 SNP_A-1909445  1   7924293
> + 7 SNP_A-2237146  2   8173763
> + 8 SNP_A-4303946  2   8191853
> + 9 SNP_A-2236357  2   8323433
> + 10 SNP_A-2205442 2   8393263"), header=TRUE)
> > # use rle to get the 'runs'
> > y <- rle(x$Chromosome)
> > # create dataframe with start/ends and values
> > start <- head(cumsum(c(1, y$lengths)), -1)
> > index <- data.frame(values=y$values, start=start, end=start + y$lengths -
> 1)
> >
> > index
>  values start end
> 1  1 1   6
> 2  2 7  10
> >
>
>
> On Nov 14, 2007 10:56 AM, affy snp <[EMAIL PROTECTED]> wrote:
> > Hello list,
> >
> > I read in a txt file using
> >
> >  >
> > by specifying the row.names=NULL so that the rows are numbered.
> > Below is an example after how the table looks like using
> >  >
> >
> >  SNPChromosome  PhysicalPosition
> > 1 SNP_A-1909444  1   7924293
> > 2 SNP_A-2237149  1   8173763
> > 3 SNP_A-4303947  1   8191853
> > 4 SNP_A-2236359  1   8323433
> > 5 SNP_A-2205441  1   8393263
> > 6 SNP_A-1909445  1   7924293
> > 7 SNP_A-2237146  2   8173763
> > 8 SNP_A-4303946  2   8191853
> > 9 SNP_A-2236357  2   8323433
> > 10 SNP_A-2205442 2   8393263
> >
> > I am wondering if there is a way to return the start and end row numbers
> > for a subset of rows.
> >
> > For example, If I specify B[,2]=1, I would like to get
> > start=1 and end=6
> >
> > if B[,2]=2, then start=7 and end=10
> >
> > Is there any way in R to quickly do this?
> >
> > Thanks a bunch!
> >
> > Allen
> >
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
> >
>
>
>
> --
> Jim Holtman
> Cincinnati, OH
> +1 513 646 9390
>
> What is the problem you are trying to solve?
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>



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

What is the problem you are trying to solve?

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


[R] executable script

2007-11-14 Thread mysimbaa

Dear All,

Apologies for this simple question and thanks in advance for any help given.

I want to make from my .R script an .exe file.

Is there any way to transfort my script to an autolaunch file?
It means it runs the script by double clicking on it.

p.s.: I'm using windows 
-- 
View this message in context: 
http://www.nabble.com/executable-script-tf4806651.html#a13751752
Sent from the R help mailing list archive at Nabble.com.

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


[R] intersect character

2007-11-14 Thread Kevin J. Thompson
hello,

  can anyone help me concatenate a string containing the intersect character?  
i need to use it in a figure legend.

thanks,
kevin

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


Re: [R] logistic regression model specification

2007-11-14 Thread Dylan Beaudette
On Nov 13, 2007 4:58 PM, Emmanuel Charpentier
<[EMAIL PROTECTED]> wrote:
> Pardon me for intruding, but I had recently to explain something
> analogous to this to a distinguished physician who had hurt himself
> playing with survival models and was bleeding graphs and analyses all
> over the place...

Thanks for taking the time to reply. I would rather be embarrassed now
than later!

> Dylan Beaudette a écrit :
>
> > On Tuesday 13 November 2007, Prof Brian Ripley wrote:
> >> On Tue, 13 Nov 2007, Dylan Beaudette wrote:
> >>> Hi,
> >>>
> >>> I have setup a simple logistic regression model with the glm() function,
> >>> with the follow formula:
> >>>
> >>> y ~ a + b
> >>>
> >>> where:
> >>> 'a' is a continuous variable stratified by
> >>> the levels of 'b'
> >>>
> >>>
> >>> Looking over the manual for model specification, it seems that
> >>> coefficients for unordered factors are given 'against' the first level of
> >>> that factor.
> >
> > Thanks for the quick reply.
> >
> >> Only for the default coding.
> >
> > Indeed, I should have added that to my initial message.
> >
> >>> This makes for difficult interpretation when using factor 'b' as a
> >>> stratifying model term.
> >> Really?  You realize that you have not 'stratified' on 'b', which would
> >> need the model to be a*b?  What you have is a model with parallel linear
> >> predictors for different levels of 'b', and if the coefficients are not
> >> telling you what you want you should change the coding.
> >
> > I should have specified that interpretation was difficult, not because of 
> > the
> > default behaviour, rather my limitations and the nature of the data. Perhaps
> > an example would help.
> >
> > y ~ a + b
> >
> > 'a' is a continuous predictor (i.e. temperature)
> > observed on the levels of 'b' (geology)
>
> Hmmm... are they fixed (as in "reproductible" or "having well known
> properties which could be expected to remain stable if the
> observation/experiment is repeated") or random (i. e. label for whatever
> variation is thought to be attribuable to local conditions). In other
> words, are you interested in the effect of a, b being a nuisance (source
> of variability), or in the effects of b by themselves ?

This gets at the real heart of the hypothesis: the logit of my
response (a soil taxonomic feature) is directly related to a energetic
parameter (temperature, solar radiance, etc.), but the relationship is
not necessarily the same across geologic type. Since we do not have an
exhaustive set of data which would encompass all of the parameters
that geology would influence, we are using it as a stratifying
variable. However, interpretation of the coefficients associated with
the levels of geology might be give some other insight.


> What glm does supposes that you are considering b as a fixed effects.
> Treating it as a random effect would involve something like Douglas
> Bates' lme4.
>
> Your reference to geology leaves the interpretation open... more on this
> later.

ok -- i do not really want to use geology as a 'treatment' per se ...
however I see where this may be leading.

>
> > The form of the model (or at least what I was hoping for) would account for
> > the variation in 'y' as predicted by 'a', within each level of 'b' . Am I
> > specifying this model incorrectly?
> >
> >> Much of what I am trying to get across is that you have a lot of choice as
> >> to how you specify a model to R. There has to be a default, which is
> >> chosen because it is often a good choice.  It does rely on factors being
> >> coded well: the 'base level' (to quote ?contr.treatment) needs to be
> >> interpretable.  And also bear in mind that the default choices of
> >> statistical software in this area almost all differ (and R's differs from
> >> S, GLIM, some ways to do this in SAS ...), so people's ideas of a 'good
> >> choice' do differ.
> >
> > Understood. I was not implying a level of 'goodness', rather hoping to gain
> > some insight into a (possibly) mis-coded model specification.
> >
> >>> Setting up the model, minus the intercept term, gives me what appear to
> >>> be more meaningful coefficients. However, I am not sure if I am
> >>> interpreting the results from a linear model without an intercept term.
> >>> Model predictions from both specifications (with and without an intercept
> >>> term) are nearly identical (different by about 1E-16 in probability
> >>> space).
> >>>
> >>> Are there any gotchas to look out for when removing the intercept term
> >>> from such a model?
> >> It is just a different parametrization of the linear predictor.
> >> Anything interpretable in terms of the predictions of the model will be
> >> unchanged.  That is the crux: the default coefficients of 'b' will be
> >> log odds-ratios that are directly interpretable, and those in the
> >> per-group coding will be log-odds for a zero value of 'a'. Does a zero
> >> value of 'a' make sense?
> >
> > In the case of this experiment, a zero-level for 'a' does not make sense.

[R] enumeration variable by groups

2007-11-14 Thread lamack lamack





Dear all, How can I create an enumeration variable by groups?
 
I have: 
gender score  1 48  1 45  2 50  2 42  1 41  2 51  1 52  1 43  2 52
 
and Y would like to get:
 
genderscoreindex
 148   1
 145   2
 141   3
 152   4
 143   5
 250   1
 242   2
 251   3
 252   4
best regards
_

[[replacing trailing spam]]

[[alternative HTML version deleted]]

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


[R] Name collisions

2007-11-14 Thread Tim Fennell
I am receiving a list of name collisions when I launch R as seen below.  
I'm new to R and any suggestion or help with how I can go about getting 
rid of these collisions would be greatly appreciated.

% R

R version 2.5.1 (2007-06-27)
Copyright (C) 2007 The R Foundation for Statistical Computing
ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

Name collision between M0 M0
Name collision between M1 M1
Name collision between M5 M5
Name collision between MI MI
Name collision between Mx Mx
Name collision between My My
Name collision between Mz Mz
Name collision between Mu Mu Mu
Name collision between Mu Mu
Name collision between M2 M2
Name collision between M3 M3


-- 
Tim Fennell
Network Administrator
Channing Laboratory


The information transmitted in this electronic communica...{{dropped:8}}

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


Re: [R] How to get row numbers of a subset of rows

2007-11-14 Thread affy snp
Thanks a lot, Jim and Bert. It worked pretty well.

Best,
 Allen

On Nov 14, 2007 12:11 PM, jim holtman <[EMAIL PROTECTED]> wrote:
> That works for the specific value of '1', but you would have to repeat
> it for other values in the column.  If you had 100 different ranges in
> that column, what would you do?  Here is another solution using
> 'range' on the same data:
>
> > tapply(seq_len(nrow(x)), x$Chromosome, range)
> $`1`
> [1] 1 6
>
> $`2`
> [1]  7 10
>
>
>
> On Nov 14, 2007 12:04 PM, Bert Gunter <[EMAIL PROTECTED]> wrote:
> > Am I missing something? ...
> >
> > Why not: range(seq(nrow(B))[B[,2]==1] ) ?? ## note: "==" not "="
> >
> > Alternatively, and easily generalized (to start with a frame which is a
> > subset of the original and any subset of rows, contiguous or not)
> >
> > range(as.numeric(row.names(B)[B[,2]==1]))
> >
> > Again, am I missing something that makes this "obvious" solution impossible?
> > (Wouldn't be the first time.)
> >
> > Bert Gunter
> > Genentech Nonclinical Statistics
> >
> >
> >
> > -Original Message-
> > From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On
> > Behalf Of jim holtman
> > Sent: Wednesday, November 14, 2007 8:39 AM
> > To: affy snp
> > Cc: r-help@r-project.org
> > Subject: Re: [R] How to get row numbers of a subset of rows
> >
> > Here is a way of doing it using 'rle':
> >
> > > x <- read.table(textConnection(" SNPChromosome
> > PhysicalPosition
> > + 1 SNP_A-1909444  1   7924293
> > + 2 SNP_A-2237149  1   8173763
> > + 3 SNP_A-4303947  1   8191853
> > + 4 SNP_A-2236359  1   8323433
> > + 5 SNP_A-2205441  1   8393263
> > + 6 SNP_A-1909445  1   7924293
> > + 7 SNP_A-2237146  2   8173763
> > + 8 SNP_A-4303946  2   8191853
> > + 9 SNP_A-2236357  2   8323433
> > + 10 SNP_A-2205442 2   8393263"), header=TRUE)
> > > # use rle to get the 'runs'
> > > y <- rle(x$Chromosome)
> > > # create dataframe with start/ends and values
> > > start <- head(cumsum(c(1, y$lengths)), -1)
> > > index <- data.frame(values=y$values, start=start, end=start + y$lengths -
> > 1)
> > >
> > > index
> >  values start end
> > 1  1 1   6
> > 2  2 7  10
> > >
> >
> >
> > On Nov 14, 2007 10:56 AM, affy snp <[EMAIL PROTECTED]> wrote:
> > > Hello list,
> > >
> > > I read in a txt file using
> > >
> > >  > >
> > > by specifying the row.names=NULL so that the rows are numbered.
> > > Below is an example after how the table looks like using
> > >  > >
> > >
> > >  SNPChromosome  PhysicalPosition
> > > 1 SNP_A-1909444  1   7924293
> > > 2 SNP_A-2237149  1   8173763
> > > 3 SNP_A-4303947  1   8191853
> > > 4 SNP_A-2236359  1   8323433
> > > 5 SNP_A-2205441  1   8393263
> > > 6 SNP_A-1909445  1   7924293
> > > 7 SNP_A-2237146  2   8173763
> > > 8 SNP_A-4303946  2   8191853
> > > 9 SNP_A-2236357  2   8323433
> > > 10 SNP_A-2205442 2   8393263
> > >
> > > I am wondering if there is a way to return the start and end row numbers
> > > for a subset of rows.
> > >
> > > For example, If I specify B[,2]=1, I would like to get
> > > start=1 and end=6
> > >
> > > if B[,2]=2, then start=7 and end=10
> > >
> > > Is there any way in R to quickly do this?
> > >
> > > Thanks a bunch!
> > >
> > > Allen
> > >
> > > __
> > > R-help@r-project.org mailing list
> > > https://stat.ethz.ch/mailman/listinfo/r-help
> > > PLEASE do read the posting guide
> > http://www.R-project.org/posting-guide.html
> > > and provide commented, minimal, self-contained, reproducible code.
> > >
> >
> >
> >
> > --
> > Jim Holtman
> > Cincinnati, OH
> > +1 513 646 9390
> >
> > What is the problem you are trying to solve?
> >
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
> >
> >
>
>
>
> --
> Jim Holtman
> Cincinnati, OH
> +1 513 646 9390
>
> What is the problem you are trying to solve?
>

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


Re: [R] label plotting on nmds diagram

2007-11-14 Thread Antonio Olinto
Hi Gavin, thanks for your message.

You are right. The function nmds is from library labdsv and not from vegan. My
mistake.

Well, that's what I did:

library(vegan)
library(labdsv)
nms.euc <- nmds(dis.euc,4)
initial  value 13.457832 
iter   5 value 8.589716
iter  10 value 7.989993
iter  15 value 7.734503
iter  20 value 7.539719
iter  25 value 7.438721
final  value 7.397254 
converged

plot(nms.euc)
text(nms.euc$points[,1],nms.euc$points[,2],labels=row.names(nms.euc$points),pos=4,cex=0.7)
# got overwritten labels

identify(plot(nms.euc))
numeric(0)

It seems that ordiplot is not suitable to nmds but only to cca and rda.

But you understood exactly what I want. Your example showed for cca what I want
for nmds plot. 

Best regards,

Antonio


Citando Gavin Simpson <[EMAIL PROTECTED]>:

> On Wed, 2007-11-14 at 11:22 -0300, Antonio Olinto wrote:
> > Hi,
> > 
> > I'm using nmds command (library vegan) to analyze some fishing data.
> 
> I doubt it - there is no nmds command. Do you mean metaMDS() ? This uses
> isoMDS() from package MASS with some extra features.
> 
> > 
> > I'd like to plot not only points, but also the names of species and
> stations in
> > a specified position.
> > 
> > I used the command
> > text(nmds$points[,1], nmds $points[,2],labels=row.names(nmds
> > $points),pos=3,cex=0.5)
> > 
> > But the labels are sometimes overlapped.
> > 
> > Is there any way to use identify, or a similar command, to plot the row
> names in
> > a given position? Identify would be perfect but it indicates the row number
> and
> > I'd like to have the row name given in nmds$points.
> 
> Have a look at ?orditorp in vegan for an alternative approach, but yes
> you can use identify, but you have to capture the output of the plot
> call and use that in your identify call:
> 
> require(vegan)
> data(varespec)
> mod <- cca(varespec)
> mod.plt <- plot(mod, display = "species", type = "p")
> identify(mod.plt)
> 
> See ?identify.ordiplot for more info.
> 
> HTH
> 
> > 
> > Many thanks,
> > 
> > Antonio Olinto
> > Sao Paulo Fisheries Institute
> > Brasil
> > 
> > 
> > 
> > 
> > 
> > 
> > -
> > BCMG Internet Webmail
> > www.bcmg.com.br
> > 
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
> -- 
> %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
>  Dr. Gavin Simpson [t] +44 (0)20 7679 0522
>  ECRC, UCL Geography,  [f] +44 (0)20 7679 0565
>  Pearson Building, [e] gavin.simpsonATNOSPAMucl.ac.uk
>  Gower Street, London  [w] http://www.ucl.ac.uk/~ucfagls/
>  UK. WC1E 6BT. [w] http://www.freshwaters.org.uk
> %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
> 
> 




-
BCMG Internet Webmail
www.bcmg.com.br

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


Re: [R] Name collisions

2007-11-14 Thread Prof Brian Ripley
On Wed, 14 Nov 2007, Tim Fennell wrote:

> I am receiving a list of name collisions when I launch R as seen below.
> I'm new to R and any suggestion or help with how I can go about getting
> rid of these collisions would be greatly appreciated.

I don't think those messages are from R, but from your OS. What platform 
is this?  Googling came up with

http://osdir.com/ml/gnu.octave.general/2004-10/msg00228.html
http://forums.mysql.com/read.php?20,84522,84522

which look relevant, so is this Solaris?  (I have never seen it there.)

I'd suggest talking to whoever installed R.

>
> % R
>
> R version 2.5.1 (2007-06-27)
> Copyright (C) 2007 The R Foundation for Statistical Computing
> ISBN 3-900051-07-0
>
> R is free software and comes with ABSOLUTELY NO WARRANTY.
> You are welcome to redistribute it under certain conditions.
> Type 'license()' or 'licence()' for distribution details.
>
> R is a collaborative project with many contributors.
> Type 'contributors()' for more information and
> 'citation()' on how to cite R or R packages in publications.
>
> Type 'demo()' for some demos, 'help()' for on-line help, or
> 'help.start()' for an HTML browser interface to help.
> Type 'q()' to quit R.
>
> Name collision between M0 M0
> Name collision between M1 M1
> Name collision between M5 M5
> Name collision between MI MI
> Name collision between Mx Mx
> Name collision between My My
> Name collision between Mz Mz
> Name collision between Mu Mu Mu
> Name collision between Mu Mu
> Name collision between M2 M2
> Name collision between M3 M3
>
>
>

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

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


[R] About print a label in plot

2007-11-14 Thread affy snp
Dear list,

Hello! I have a question about how to print a label in the plot.
I am using the following code:

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


Re: [R] enumeration variable by groups

2007-11-14 Thread jim holtman
Here is a way to do it:

> x <- scan(textConnection("1 48  1 45  2 50  2 42  1 41  2 51  1 52  1 43  2 
> 52"), what=0L)
Read 18 items
> x <- matrix(x, ncol=2, byrow=TRUE)
> colnames(x) <- c('gender', 'score')
> x
  gender score
 [1,]  148
 [2,]  145
 [3,]  250
 [4,]  242
 [5,]  141
 [6,]  251
 [7,]  152
 [8,]  143
 [9,]  252
> # split out categories
> y <- split(seq_len(nrow(x)), x[,1])
> # combine into new matrix
> x.new <- do.call('rbind', lapply(y, function(.row) cbind(x[.row,], 
> index=seq_along(.row
> x.new
  gender score index
 [1,]  148 1
 [2,]  145 2
 [3,]  141 3
 [4,]  152 4
 [5,]  143 5
 [6,]  250 1
 [7,]  242 2
 [8,]  251 3
 [9,]  252 4
>
>



On Nov 14, 2007 12:58 PM, lamack lamack <[EMAIL PROTECTED]> wrote:
>
>
>
>
>
> Dear all, How can I create an enumeration variable by groups?
>
> I have:
> gender score  1 48  1 45  2 50  2 42  1 41  2 51  1 52  1 43  2 52
>
> and Y would like to get:
>
>genderscoreindex
> 148   1
> 145   2
> 141   3
> 152   4
> 143   5
> 250   1
> 242   2
> 251   3
> 252   4
> best regards
> _
>
> [[replacing trailing spam]]
>
>[[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



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

What is the problem you are trying to solve?

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


[R] reading tables from url

2007-11-14 Thread stubben
I'm trying to read some web tables directly into R.  These are both  
genome sequencing projects (eukaryotes and metagenomes) from NCBI and  
look very similar;  however, only the first one works.

http://www.ncbi.nlm.nih.gov/genomes/leuks.cgi
http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi

I added  ?dump=selected to the end of the url string to get a tab- 
delimited file (which is what happens if you click the Save button on  
either page).

 > options(internet.info=0)

## this one works

 > x1<-url("http://www.ncbi.nlm.nih.gov/genomes/leuks.cgi? 
dump=selected")
 > read.delim(x1, skip=1, nrows=5)[,1:3]

   X...Columns. ProjectID Organism.Name
120303 Acanthamoeba castellanii Neff  Protists
213657  Acyrthosiphon pisum LSR1   Animals
312434   Aedes aegypti Liverpool   Animals
412635 Ajellomyces capsulatus G186AR Fungi
512653  Ajellomyces capsulatus G217B Fungi

Warning messages:
1: connected to 'www.ncbi.nlm.nih.gov' on port 80. in: open.connection 
(file, "r")
2: -> GET /genomes/leuks.cgi?dump=selected HTTP/1.0

Host: www.ncbi.nlm.nih.gov

Pragma: no-cache

in: open.connection(file, "r")
3: <- HTTP/1.1 200 OK in: open.connection(file, "r")
4: <- Date: Wed, 14 Nov 2007 18:03:29 GMT in: open.connection(file, "r")
5: <- Server: Apache in: open.connection(file, "r")
6: <- Content-Disposition: attachment; filename="untitle.txt" in:  
open.connection(file, "r")
7: <- Content-Type: application/force-download in: open.connection 
(file, "r")
8: <- Vary: Accept-Encoding in: open.connection(file, "r")
9: <- Connection: close in: open.connection(file, "r")
10: Code 200, content-type 'application/force-download' in:  
open.connection(file, "r")


## this one fails to open a connection

 > x2<-url("http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi? 
dump=selected")
 > read.delim(x2, skip=1, nrows=5)[,1:3]

Error in open.connection(file, "r") : unable to open connection
In addition: Warning messages:
1: connected to 'www.ncbi.nlm.nih.gov' on port 80. in: open.connection 
(file, "r")
2: -> GET /genomes/lenvs.cgi?dump=selected HTTP/1.0

Host: www.ncbi.nlm.nih.gov

Pragma: no-cache

in: open.connection(file, "r")
3: <- HTTP/1.1 500 Internal Server Error in: open.connection(file, "r")
4: <- Date: Wed, 14 Nov 2007 18:04:26 GMT in: open.connection(file, "r")
5: <- Server: Apache in: open.connection(file, "r")
6: <- Content-Type: text/html; charset=ISO-8859-1 in: open.connection 
(file, "r")
7: <- Vary: Accept-Encoding in: open.connection(file, "r")
8: <- Connection: close in: open.connection(file, "r")
9: Code 500, content-type 'text/html; charset=ISO-8859-1' in:  
open.connection(file, "r")
10: cannot open: HTTP status was '500 Internal Server Error' in:  
open.connection(file, "r")

Also, I can't even read lines from the main page.

 > readLines("http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi";, n=10)
Error in file(con, "r") : unable to open connection
...
## now I'm just guessing...
 > readLines("http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi";, n=10,  
encoding="ISO-8859-1")
Error in file(con, "r") : unable to open connection
...


Download.file works fine, but I would like to avoid this if possible.

 > capabilities()[5]
http/ftp
 TRUE

 > download.file("http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi? 
dump=selected", "lenvs.tab")
 > read.delim("lenvs.tab", skip=1, nrows=5)[,1:3]
   X...Columns.  
Parent.ProjectID ProjectID
11973313694   Global Ocean Sampling  
Expedition Metagenome
22082313696  5-Way (CG) Acid Mine Drainage  
Biofilm Metagenome
3-13699Waseca County Farm  
Soil Metagenome
4-13702 Methane-Oxidizing Archaea from Deep- 
Sea Sediments
5-13729 Pacific Beach  
Sand Metagenome



Thanks for your help.  Hopefully this is something simple that I  
missed in the documentation/help.

Chris



--
---
Chris Stubben

Los Alamos National Lab
BioScience Division
MS M888
Los Alamos, NM 87545

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


Re: [R] label plotting on nmds diagram

2007-11-14 Thread Gavin Simpson
On Wed, 2007-11-14 at 15:19 -0300, Antonio Olinto wrote:
> Hi Gavin, thanks for your message.
> 
> You are right. The function nmds is from library labdsv and not from vegan. My
> mistake.
> 
> Well, that's what I did:
> 
> library(vegan)
> library(labdsv)
> nms.euc <- nmds(dis.euc,4)
> initial  value 13.457832 
> iter   5 value 8.589716
> iter  10 value 7.989993
> iter  15 value 7.734503
> iter  20 value 7.539719
> iter  25 value 7.438721
> final  value 7.397254 
> converged
> 
> plot(nms.euc)
> text(nms.euc$points[,1],nms.euc$points[,2],labels=row.names(nms.euc$points),pos=4,cex=0.7)
> # got overwritten labels
> 
> identify(plot(nms.euc))
> numeric(0)
> 
> It seems that ordiplot is not suitable to nmds but only to cca and rda.

Hi Antonio,

You haven't used ordiplot there or created an object of class "ordiplot"
so the problem is that you expected identify to work with an arbitrary
object. Someone would need to write an identify method that could deal
with whatever plot.nmds returns for what you have above to work.

ordiplot is suitable for lots of ordination types (effectively anything
that scores() can deal with) I was just not familiar with nmds() from
labdsv. Having looked at this, nmds() is an even simpler wrapper round
isoMDS() from MASS.

> 
> But you understood exactly what I want. Your example showed for cca what I 
> want
> for nmds plot.

Only because it is easier to write the cca code without having to wait
for metaMDS to finish the random starts ;-)

Here's what I think you need but using vegan and metaMDS for the random
starts and other NMDS sugar...:

require(vegan)
data(varespec)
nms.euc <- metaMDS(varespec, distance = "euclidean", k = 4)
plt <- plot(nms.euc, display = "sites", type = "p")
identify(plt, what = "sites")

I also forgot to mention that NMDS doesn't have species scores as such,
so your example has nothing to do with species as far as I can see - did
you omit something else?

vegan does allow you to get "species" scores via weighted averages
(see ?scores.metaMDS ), e.g.

plt2 <- plot(nms.euc, display = c("species","sites"), type = "p", shrink
= TRUE)
identify(plt2, what = "species", cex = 0.8, col = "red")
identify(plt2, what = "sites", cex = 0.8, col = "black")

Is this what you are looking for?

If you stick with labdsv:::nmds then take a look at plotid.nmds() which
is a wrapper around identify for nmds objects, which appears to do the
hard work for you. You could look at this to see how you could have gone
about crafting the identify call yourself - the function is very simple.

Alternatively, do look at orditorp() in vegan, as it is a nice function
that allows you to display only labels for species/sites that can fit:

plot(nms.euc, display = c("species","sites"), type = "n", shrink = TRUE)
points(nms.euc, display = "sites", cex = 0.7)
orditorp(nms.euc, display = "species", shrink = TRUE, col = "red")

The orditorp call will generate some warnings - they are harmless and
have been fixed in the latest development version of vegan on R-Forge.
Jari is planning a release of vegan real soon now that will contain this
fix.

Note that orditorp has a priority argument, which you can use to give
labelling priority to species/sites, say based on their maximal
abundance or number of occurrences or similar; see ?orditorp

HTH

G

> 
> Best regards,
> 
> Antonio
> 
> 
> Citando Gavin Simpson <[EMAIL PROTECTED]>:
> 
> > On Wed, 2007-11-14 at 11:22 -0300, Antonio Olinto wrote:
> > > Hi,
> > > 
> > > I'm using nmds command (library vegan) to analyze some fishing data.
> > 
> > I doubt it - there is no nmds command. Do you mean metaMDS() ? This uses
> > isoMDS() from package MASS with some extra features.
> > 
> > > 
> > > I'd like to plot not only points, but also the names of species and
> > stations in
> > > a specified position.
> > > 
> > > I used the command
> > > text(nmds$points[,1], nmds $points[,2],labels=row.names(nmds
> > > $points),pos=3,cex=0.5)
> > > 
> > > But the labels are sometimes overlapped.
> > > 
> > > Is there any way to use identify, or a similar command, to plot the row
> > names in
> > > a given position? Identify would be perfect but it indicates the row 
> > > number
> > and
> > > I'd like to have the row name given in nmds$points.
> > 
> > Have a look at ?orditorp in vegan for an alternative approach, but yes
> > you can use identify, but you have to capture the output of the plot
> > call and use that in your identify call:
> > 
> > require(vegan)
> > data(varespec)
> > mod <- cca(varespec)
> > mod.plt <- plot(mod, display = "species", type = "p")
> > identify(mod.plt)
> > 
> > See ?identify.ordiplot for more info.
> > 
> > HTH
> > 
> > > 
> > > Many thanks,
> > > 
> > > Antonio Olinto
> > > Sao Paulo Fisheries Institute
> > > Brasil
> > > 
> > > 
> > > 
> > > 
> > > 
> > > 
> > > -
> > > BCMG Internet Webmail
> > > www.bcmg.com.br
> > > 
> > > __
> > > R-help@r-projec

Re: [R] How to get row numbers of a subset of rows

2007-11-14 Thread Julian Burgos
One way to do this is

range(which(B[,2]==1))

Julian

affy snp wrote:
> Hello list,
> 
> I read in a txt file using
> 
>  
> by specifying the row.names=NULL so that the rows are numbered.
> Below is an example after how the table looks like using
>  
> 
>   SNPChromosome  PhysicalPosition
> 1 SNP_A-1909444  1   7924293
> 2 SNP_A-2237149  1   8173763
> 3 SNP_A-4303947  1   8191853
> 4 SNP_A-2236359  1   8323433
> 5 SNP_A-2205441  1   8393263
> 6 SNP_A-1909445  1   7924293
> 7 SNP_A-2237146  2   8173763
> 8 SNP_A-4303946  2   8191853
> 9 SNP_A-2236357  2   8323433
> 10 SNP_A-2205442 2   8393263
> 
> I am wondering if there is a way to return the start and end row numbers
> for a subset of rows.
> 
> For example, If I specify B[,2]=1, I would like to get
> start=1 and end=6
> 
> if B[,2]=2, then start=7 and end=10
> 
> Is there any way in R to quickly do this?
> 
> Thanks a bunch!
> 
> Allen
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] executable script

2007-11-14 Thread Peter Alspach
Kia ora unknown requestor

One option is to create a .bat file along the following lines:

path_to_R\bin\R CMD BATCH yourScript.R


Peter Alspach


> -Original Message-
> From: [EMAIL PROTECTED] 
> [mailto:[EMAIL PROTECTED] On Behalf Of mysimbaa
> Sent: Thursday, 15 November 2007 6:10 a.m.
> To: r-help@r-project.org
> Subject: [R] executable script
> 
> 
> Dear All,
> 
> Apologies for this simple question and thanks in advance for 
> any help given.
> 
> I want to make from my .R script an .exe file.
> 
> Is there any way to transfort my script to an autolaunch file?
> It means it runs the script by double clicking on it.
> 
> p.s.: I'm using windows
> --
> View this message in context: 
> http://www.nabble.com/executable-script-tf4806651.html#a13751752
> Sent from the R help mailing list archive at Nabble.com.
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 

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


Re: [R] How to get row numbers of a subset of rows

2007-11-14 Thread affy snp
Thank you very much, Julian. I got it.

Best,
Allen

On Nov 14, 2007 2:38 PM, Julian Burgos <[EMAIL PROTECTED]> wrote:
> One way to do this is
>
> range(which(B[,2]==1))
>
> Julian
>
>
> affy snp wrote:
> > Hello list,
> >
> > I read in a txt file using
> >
> >  >
> > by specifying the row.names=NULL so that the rows are numbered.
> > Below is an example after how the table looks like using
> >  >
> >
> >   SNPChromosome  PhysicalPosition
> > 1 SNP_A-1909444  1   7924293
> > 2 SNP_A-2237149  1   8173763
> > 3 SNP_A-4303947  1   8191853
> > 4 SNP_A-2236359  1   8323433
> > 5 SNP_A-2205441  1   8393263
> > 6 SNP_A-1909445  1   7924293
> > 7 SNP_A-2237146  2   8173763
> > 8 SNP_A-4303946  2   8191853
> > 9 SNP_A-2236357  2   8323433
> > 10 SNP_A-2205442 2   8393263
> >
> > I am wondering if there is a way to return the start and end row numbers
> > for a subset of rows.
> >
> > For example, If I specify B[,2]=1, I would like to get
> > start=1 and end=6
> >
> > if B[,2]=2, then start=7 and end=10
> >
> > Is there any way in R to quickly do this?
> >
> > Thanks a bunch!
> >
> > Allen
> >
>
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
>
>

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


[R] where is library file

2007-11-14 Thread Casey,Richard
Hello,

While trying to install R on 64-bit HP Proliant and RedHat Enterprise Linux v.4 
using R-2.6.0-3.rh4.x86_64.rpm, it says:

Failed dependencies: libg2c.so.0()(64bit) is needed by R

Does anyone know where to get libg2c.so or the rpm to install it?



Richard

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


Re: [R] Generating these matrices going backwards

2007-11-14 Thread Chris Stubben

Maybe something like this?

N<-c(1,2,1,3)
## create empty matrix
x<-diag(0,3)
## fill off diagonal
 x[row(x)==col(x)+1]<-N[1:2]
# fill 3rd column
 x[1:2,3]<-N[3:4]

Or create a function and return both x and y

mat3<-function(N)
{
 x<-diag(0,3)
# fill each element separately
x[2,1]<-N[1]
x[3,2]<-N[2]
x[1,3]<-N[3]
x[2,3]<-N[4]
dimnames(x)<-list(c("D1", "D2", "D3"), c("E1", "E2", "E3"))
y =sum(x) * x / (rowSums(x)%o%colSums(x)) 
list(x=x,y=y)
}
mat3(N)

$x
   E1 E2 E3
D1  0  0  1
D2  1  0  3
D3  0  2  0

$y
 E1  E2 E3
D1 0.00 0.0 1.7500
D2 1.75 0.0 1.3125
D3 0.00 3.5 0.




francogrex wrote:
> 
> I have generated the following:
> 
> x= 
> E1  E2  E3 
> D1  0   0   1 
> D2  1   0   3 
> D3  0   2   0 
> 
> y= 
> E1  E2  E3 
> D1  0   0   1.75 
> D2  1.750   1.3125 
> D3  0   3.5 0 
> 
> Where x and y are linked by: 
> y =sum(x) * x / (rowSums(x)%o%colSums(x)) 
> 
> N=x[x[1:3,]>0] 
> R=y[y[1:3,]>0] 
> 
> Now suppose I ONLY have N and R linked in this way below where each N 
> corresponds to an R 
> 
> N   R 
> 1   1.7500 
> 2   3.5000 
> 1   1.7500 
> 3   1.3125 
> 
> Is there a way to generate matrix "x" and matrix "y" having only the N 
> and the R as above? 
> 
> 
> 

-- 
View this message in context: 
http://www.nabble.com/Generating-these-matrices-going-backwards-tf4807447.html#a13756103
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] map - mapproj : problem of states localisation

2007-11-14 Thread Ray Brownrigg
On Wed, 14 Nov 2007, Amandine Chevalier wrote:
> Thank you very much for your help,
>
> I work about my code so as to give you the example (attached file) :
>
> In the first case (azequalarea projection), france from the "france" map
> and from the "world" map are well superimposed, whereas in the second case
> (lambert projection), it is not. However, due to projections I used to
> produce my data, it would be better to use the lambert projection...That is
> not only a problem of boundaries that do not line up...it is like the map
> of France (it would be the same with the Italy's map) was centered on the
> figure without taking into account the lat-long...Can I also trust the
> location of the world map ?
>
> Do you have any idea of the problem?
>
> Thank you veru much for your help,
>
> Regards,
>
> Amandine
>
It seems to me the problem doesn't need 200 lines of code and data to 
describe, it is merely:
library(mapproj)
map("france", proj="lambert", pa=c(30, 60))
map("world", proj="lambert", pa=c(30, 60),  lwd=1, col="red", add=TRUE)

So I've punted it to the maintainer of mapproj, who isn't me.

Ray

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


[R] Changing the text in the strips of lattice plots and y axis

2007-11-14 Thread Judith Flores
Dear R-helpers,

   I am sorry for asking something I know has been
asked before, I have tried different combinations in
the strip function without success...


I am using version 2.5.1 and work on a PC.

   I have  barcharts generated from the following
formula:

barchart(y1+y2+y3~x | g)

I need to change the names of the variables y1,y2 or
y3 that currently appear in the strips, I have two
strips per panel, one corresponds to the name of the
variable (y1 or y2) and the second strip to the value
of the conditional variable (g has two levels, this
means I have 6 panels total).

  Instead of y1,y2 and y3 I would like the strips to
read 'name of variable y1', 'name of variable y2',
'name of variable y3'; of course these are not the
real names I want to assign to those variables.

  And my second question is regarding the the limits
of the y axis. If I setup scales to relation='free', I
obtain 6 different scales, which is expected. But I
need the panels that correspond to each y variable
have the same scale. I tried something like this:

prepanel=function(x,y,...) {ylim=list(min(y),max(y)) }

   but didn't work...

What do I need to do?

Thank you very much in advance for your help,

Judith


  

Never miss a thing.  Make Yahoo your home page.

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


Re: [R] where is library file

2007-11-14 Thread Prof Brian Ripley
On Wed, 14 Nov 2007, Casey,Richard wrote:

> While trying to install R on 64-bit HP Proliant and RedHat Enterprise 
> Linux v.4 using R-2.6.0-3.rh4.x86_64.rpm, it says:
>
> Failed dependencies: libg2c.so.0()(64bit) is needed by R
>
> Does anyone know where to get libg2c.so or the rpm to install it?

libg2c is the runtime library for g77, the fortran compiler for gcc 3.x.

Loooking at the R spec file, this is likely to be in gcc-g77 (or a 
dependency of that), and that is a dependency of the R-devel RPM.  Have 
you tried installing the latter: you will need it to install packages 
later?

Googling got me to

http://www.redhat.com/archives/nahant-list/2007-April/msg00049.html

which suggests that you need the libf2c RPM on RHEL4.

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

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


[R] differing scales for lattice plot panels

2007-11-14 Thread Joe Bloggs
Hi,
   I want to compare densities across several groups with: 
densityplot(~values | groups, pch=".")
However one of the groups has an almost point-like distribution. 
This forces the y-axis scale to be so large that the other distribution curves 
are all flat in comparison and get lost amongst the data dots.
How can I produce a lattice plot with different y-axes (preferably 
automatically scaled) for each panel?

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


[R] Error while calculating AOV

2007-11-14 Thread apsawant

Hi ,

I am new to R.
I am trying to run a simple example script to calculate AOV.
Below is the script file (aov.R) I am trying to execute:
aov.R
--
data1<-c(49,47,46,47,48,47,41,46,43,47,46,45,48,46,47,45,49,44,44,45,42,45,45,40
,49,46,47,45,49,45,41,43,44,46,45,40,45,43,44,45,48,46,40,45,40,45,47,40)

matrix(data1, ncol= 4, dimnames = list(paste("subj", 1:12),
c("Shape1.Color1",
"Shape2.Color1", "Shape1.Color2", "Shape2.Color2")))

Hays.df<-data.frame(rt = data1,
subj=factor(rep(paste("subj", 1:12, sep=""), 4)),
shape=factor(rep(rep(c("shape1","shape2"), c(12, 12)), 2)),
color=factor(rep(c("color1","color2"), c(24, 24

aov(rt - shape * color + Error(subj/(shape * color)), data=Hays.df)

summary(aov(rt - shape * color + Error(subj/(shape * color)), data=Hays.df))


Here is the error I get:
> source("aov.R")
Error in terms(formula, "Error", data = data) : 
Object "shape" not found

I would appreciate your help as to why I am getting this error message.

Thanks,
Amit.
-- 
View this message in context: 
http://www.nabble.com/Error-while-calculating-AOV-tf4799288.html#a13730481
Sent from the R help mailing list archive at Nabble.com.

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


[R] Log random number

2007-11-14 Thread Tobias Schlottmann
 
   
   Dear R users,
   
  Simply my question is that how it is possible to generate some random numbers 
using rnorm( ) but in log transformed values.
   
  Thank you,
   
  Tobias

   
-

[[alternative HTML version deleted]]

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


[R] Generating these matrices going backwards

2007-11-14 Thread francogrex

I have generated the following:

x= 
E1  E2  E3 
D1  0   0   1 
D2  1   0   3 
D3  0   2   0 

y= 
E1  E2  E3 
D1  0   0   1.75 
D2  1.750   1.3125 
D3  0   3.5 0 

Where x and y are linked by: 
y =sum(x) * x / (rowSums(x)%o%colSums(x)) 

N=x[x[1:3,]>0] 
R=y[y[1:3,]>0] 

Now suppose I ONLY have N and R linked in this way below where each N 
corresponds to an R 

N   R 
1   1.7500 
2   3.5000 
1   1.7500 
3   1.3125 

Is there a way to generate matrix "x" and matrix "y" having only the N 
and the R as above? 


-- 
View this message in context: 
http://www.nabble.com/Generating-these-matrices-going-backwards-tf4807447.html#a13754521
Sent from the R help mailing list archive at Nabble.com.

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


[R] calc dist AB for B at levels of A

2007-11-14 Thread Melanie Ann Harsch
I have a question that I have been trying to figure out and I imagine 
there is a very simple answer to it.

I am trying to use the distAB function in the clim.pact package
I have a dataframe with 4 columns in which A is ref pt and B is site
1. longitude of A
2. Latitude of A
3. longitude of B
4. Latitude of B

Problem is the columns are unequal in length 
column 1 and 2 are length=29 and 3 and 4 length=53

What I want to do is create a matrix which are filled with values of 
distAB for B (site) at all levels of A (reference points)

Any ideas?
Melanie

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


Re: [R] Log random number

2007-11-14 Thread Wensui Liu
dont think so, unless i miss something here.
please do check the range of normal random number and the domain of
log function.

On 11/14/07, Tobias Schlottmann <[EMAIL PROTECTED]> wrote:
>
>
>Dear R users,
>
>   Simply my question is that how it is possible to generate some random 
> numbers using rnorm( ) but in log transformed values.
>
>   Thank you,
>
>   Tobias
>
>
> -
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>


-- 
===
WenSui Liu
Statistical Project Manager
ChoicePoint Precision Marketing
(http://spaces.msn.com/statcompute/blog)

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


Re: [R] Piecewise Linear Regression

2007-11-14 Thread Vadim Ogranovich
Hi, 

Let me pick up this old thread. How does one extract the locations of the knots 
(ends of the segments) from the fit object below? 

Thanks, 
Vadim 


>From : roger koenker < roger_at_ysidro.econ.uiuc.edu > 
Date : Tue 31 May 2005 - 10:23:19 EST 




It is conventional to fit piecewise linear models by assuming Gaussian error 
and 
using least squares methods, but one can argue that median regression provides 
a more robust approach to this problem. You might consider the following fit: 

x = c 
(6.25,6.25,12.50,12.50,18.75,25.00,25.00,25.00,31.25,31.25,37.50,37.50,5 
0.00,50.00,62.50,62.50,75.00,75.00,75.00,100.00,100.00) y = c 
(0.328,0.395,0.321,0.239,0.282,0.230,0.273,0.347,0.211,0.210,0.259,0.186 
,0.301,0.270,0.252,0.247,0.277,0.229,0.225,0.168,0.202) library(quantreg) 
plot(x,y) 
fit <- rqss(y ~ qss(x)) 
plot(fit) 

it gives 5 segments not 3, but this can be controlled by the choice of lambda 
in the qss 
function, for example, try: 

fit <- rqss(y ~ qss(x,lambda=3) 
plot(fit,col="red") 

which gives a fit like you suggest might be reasonable with only three 
segments. ... some text removed ... 




On May 30, 2005, at 6:38 PM, Abhyuday Mandal wrote: 

> Hi, 
> 
> I need to fit a piecewise linear regression. 
> 
> x = c 
> (6.25,6.25,12.50,12.50,18.75,25.00,25.00,25.00,31.25,31.25,37.50,37.50 
> ,50.00,50.00,62.50,62.50,75.00,75.00,75.00,100.00,100.00) 
> y = c 
> (0.328,0.395,0.321,0.239,0.282,0.230,0.273,0.347,0.211,0.210,0.259,0.1 
> 86,0.301,0.270,0.252,0.247,0.277,0.229,0.225,0.168,0.202) 
> 
> there are two change points. so the fitted curve should look like 


[[alternative HTML version deleted]]

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


Re: [R] Log random number

2007-11-14 Thread Ravi Varadhan
If you want to generate log-normally distributed random numbers, see:

?rlnorm

Note: X is said to be log-normally distributed when log(X) is normally
distributed.

Ravi.


---

Ravi Varadhan, Ph.D.

Assistant Professor, The Center on Aging and Health

Division of Geriatric Medicine and Gerontology 

Johns Hopkins University

Ph: (410) 502-2619

Fax: (410) 614-9625

Email: [EMAIL PROTECTED]

Webpage:  http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html

 





-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On
Behalf Of Wensui Liu
Sent: Wednesday, November 14, 2007 4:51 PM
To: Tobias Schlottmann
Cc: [EMAIL PROTECTED]
Subject: Re: [R] Log random number

dont think so, unless i miss something here.
please do check the range of normal random number and the domain of
log function.

On 11/14/07, Tobias Schlottmann <[EMAIL PROTECTED]> wrote:
>
>
>Dear R users,
>
>   Simply my question is that how it is possible to generate some random
numbers using rnorm( ) but in log transformed values.
>
>   Thank you,
>
>   Tobias
>
>
> -
>
> [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>


-- 
===
WenSui Liu
Statistical Project Manager
ChoicePoint Precision Marketing
(http://spaces.msn.com/statcompute/blog)

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

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


[R] Integrating a two variable function

2007-11-14 Thread Robert Schneider

Hi,

I would like to integrate over x the following function:

f1 <- function(x,int) x^(p+p.int*int) * (1-x)^(q+q.int)

where

- p, p.int, q, q.int: fixed parameters

- int is a binary variable (0 or 1), and specific to a given x, i.e.

 x=c(0,0.14,0.29,0.32,...,1)
 int=c(0,0,1,0,1,1,...,0)

I would like to calculate the area under the curve between each value of x. I 
can get this using the integrate function without taking into account the 
binary variable, i.e. fixing int to 0 or 1:

lim = length(x)
result.0 <- c()
result.1 <- c()
for (i in 1:lim-1) {
result.0[i] <- integrate(f1,x[i],x[i+1],int=0)
result.1[i] <- integrate(f1,x[i],x[i+1],int=1)
}

Is there a way to take into account the second variable, i.e. for the 
integration to take into account int for each x ?

Thanks in advance.

Robert Schneider

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


Re: [R] graphics in cocorresp package

2007-11-14 Thread Gavin Simpson
On Wed, 2007-11-14 at 16:46 +0100, Elisa Santi wrote:
> I'm using cocorresp package with predictive co-correspondence analysis
> model.
> i would need to visulize inside the graphic (plot function) the sample
> site numbers (is it possible to consider this a "label"?).
> is there any function to do this?
> 
> thanks for your help.
> 
> Elisa

Dear Elisa,

The posting guide asks you to direct questions about contributed
packages to the maintainer not the list. So that'd be me then...

The default plotting functions in cocorresp are very rudimentary and
reflect what I knocked up at the time to visualise the output simply.
I've somewhat neglected the package for a while due to work commitments
and a need to finish another package for a research project. Updating
the package is on my TODO list but that is a rather long list at the
moment.

In the interim, here's an example of how to alter the plotting which
will hopefully get you started. If you have other questions, please
direct them to me rather than the list.

require(cocorresp)
data(beetles)
data(plants)
bp.pred <- coca(beetles ~ ., data = plants, n.axes = 2)

## get the scores for plotting
bp.sc <- scores(bp.pred, display = c("species", "site"))

## beetles is the response, lets plot those scores; sites first
plot(bp.sc$site$X1[, 1], bp.sc$site$X1[, 2], type = "n", asp = 1)
## if you want labels then
text(bp.sc$site$X1[, 1], bp.sc$site$X1[, 2], 
 labels = rownames(bp.sc$site$X1), cex = 0.7)
## or just the points
points(bp.sc$site$X1[, 1], bp.sc$site$X1[, 2], cex = 0.7)

## or we can plot the species with labels
plot(bp.sc$species$U1[, 1], bp.sc$species$U1[, 2], type = "n", asp = 1)
## then draw species as labels
text(bp.sc$species$U1[, 1], bp.sc$species$U1[, 2], 
 labels = rownames(bp.sc$species$U1), cex = 0.5)

Problem is, these are terribly crowded. A better solution might be to
allow you the user to label the points. Here is a solution using
identify()

## solution using identify()
plot(bp.sc$species$U1[, 1], bp.sc$species$U1[, 2], asp = 1, cex = 0.5)
identify(bp.sc$species$U1[, 1], bp.sc$species$U1[, 2], 
 labels = rownames(bp.sc$species$U1), cex = 0.5, col = "red")
## click some points to label them, right click the plot window to end

Note that if your want the predictor data set (in my example here these
are the scores for plants), then switch the X1 for X2 or U1 for U2 in
bp.sc$site$X1 or bp.sc$species$U1 - yes, I know this is not user
friendly, but it matches the Matlab code from which cocorresp was
ported, but will be changed in a future version of the package.

Does this help?

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

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


Re: [R] reading tables from url

2007-11-14 Thread Prof Brian Ripley
On Thu, 15 Nov 2007, Duncan Temple Lang wrote:

> -BEGIN PGP SIGNED MESSAGE-
> Hash: SHA1
>
> Hi Chris.
>
>  Indeed, I cannot connect to that URL either.  So I did a bit of
> digging and experimentation to find out whether one needed to
> pass additional hidden options from the form or whether the problem was
> more to do with how we connect.
>
> It turns out that the script associated with NCBI leuks.cgi is being
> fussy and wants you tell it the user agent that is performing the
> request.  (Why the two behave differently is not clear after a very
> brief look, but it is probably not worth pursuing.)
>
> AFAIR, there is no way to tell R to include a UserAgent field in the
> header of the request using url(), etc. although it did come up at one
> point.

Which is presumably why download.file works: it does set a user agent (see 
option HTTPUserAgent). I don't understand why it was added to 
download.file but not url: but that is easy to rectify.

(I don't understand why the server is set up to give an internal server 
error for a missing field, though).


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

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


Re: [R] intersect character

2007-11-14 Thread Peter Alspach
Kevin

How about, for example,

plot(1, xlab=quote(B~intersect(A)))

Peter Alspach
 

> -Original Message-
> From: [EMAIL PROTECTED] 
> [mailto:[EMAIL PROTECTED] On Behalf Of Kevin J. Thompson
> Sent: Thursday, 15 November 2007 5:51 a.m.
> To: R-help
> Subject: [R] intersect character
> 
> hello,
> 
>   can anyone help me concatenate a string containing the 
> intersect character?  i need to use it in a figure legend.
> 
> thanks,
> kevin
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 

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


[R] graphics in cocorresp package

2007-11-14 Thread Elisa Santi

I'm using cocorresp package with predictive co-correspondence analysis model.
i would need to visulize inside the graphic (plot function) the sample site 
numbers (is it possible to consider this a "label"?).
is there any function to do this?
 
thanks for your help.
 
Elisa
 
 
 
Dr. Elisa Santi
Dipartimento di Scienze Ambientali "G. Sarfatti" Unità di Ricerca Biodiversità 
ed Ecologia delle Comunità Università degli Studi di SienaVia P.A. Mattioli 4, 
I-53100 Siena, Italy
Cell +39 338 4069230 - Email: [EMAIL PROTECTED]
 
  
_


[[alternative HTML version deleted]]

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


Re: [R] label plotting on nmds diagram

2007-11-14 Thread Dave Roberts
Using nmds from package labdsv will work, but you are confusing some 
elements with ordiplot from package vegan.  The plot function in labdsv 
does not produce a new object, it simply plots the requested ordination. 
  There is an accompanying command "plotid" that puts labels on the 
ordonation.  By default, it simply numbers them, but you can specify a 
vector of names of labels with the optional argument ids= (see 
?plotid.nmds), so

library(labdsv)
data(bryceveg)
data(brycesite)
dis.bc <- dsvdis(bryceveg,'bray') # Bray-Curtis dissimilarity
nmds.bc <- nmds(dis.bc)
plot(nmds.bc)
plotid(nmds.bc,ids=row.names(brycesite))

will put more complex row.names on items selected with a mouse.  It's 
tedious if you want to label ALL the points. but it will allow you to 
minimize overlapping labels.

You can also use it to print the values of categorical variables, rather 
than glyphs,

plotid(nmds.bc,ids=brycesite$depth)

Gavin is correct that NMDS doesn't have species coordinates, only sample 
coordinates.  package vegan addresses this by calculating the weighted 
centroid of all samples containing a species as the species location, 
using function 'wascores'.  This is a post-hoc calculation, and species 
don't influence the configuration of the ordination except as they 
influence the underlying dissimilarity/distance matrix.  wascores also 
works with nmds, as

spec.loc <- wascores(nmds.bc$points,bryceveg)
points(spec.loc,pch="+",col=2)

Alternatively, you could switch the entire operation to package vegan, 
using metaMDS and associated plotting routines.  metaMDS uses a number 
of arguments to specify whether to use step-acxross distances, and uses 
a number of random initial starts.  The default with nmds is to use a 
pco as the initial conditions, but random starts are also available in 
labdsv with bestnmds.


Good luck, Dave Roberts

Gavin Simpson wrote:
> On Wed, 2007-11-14 at 15:19 -0300, Antonio Olinto wrote:
>> Hi Gavin, thanks for your message.
>>
>> You are right. The function nmds is from library labdsv and not from vegan. 
>> My
>> mistake.
>>
>> Well, that's what I did:
>>
>> library(vegan)
>> library(labdsv)
>> nms.euc <- nmds(dis.euc,4)
>> initial  value 13.457832 
>> iter   5 value 8.589716
>> iter  10 value 7.989993
>> iter  15 value 7.734503
>> iter  20 value 7.539719
>> iter  25 value 7.438721
>> final  value 7.397254 
>> converged
>>
>> plot(nms.euc)
>> text(nms.euc$points[,1],nms.euc$points[,2],labels=row.names(nms.euc$points),pos=4,cex=0.7)
>> # got overwritten labels
>>
>> identify(plot(nms.euc))
>> numeric(0)
>>
>> It seems that ordiplot is not suitable to nmds but only to cca and rda.
> 
> Hi Antonio,
> 
> You haven't used ordiplot there or created an object of class "ordiplot"
> so the problem is that you expected identify to work with an arbitrary
> object. Someone would need to write an identify method that could deal
> with whatever plot.nmds returns for what you have above to work.
> 
> ordiplot is suitable for lots of ordination types (effectively anything
> that scores() can deal with) I was just not familiar with nmds() from
> labdsv. Having looked at this, nmds() is an even simpler wrapper round
> isoMDS() from MASS.
> 
>> But you understood exactly what I want. Your example showed for cca what I 
>> want
>> for nmds plot.
> 
> Only because it is easier to write the cca code without having to wait
> for metaMDS to finish the random starts ;-)
> 
> Here's what I think you need but using vegan and metaMDS for the random
> starts and other NMDS sugar...:
> 
> require(vegan)
> data(varespec)
> nms.euc <- metaMDS(varespec, distance = "euclidean", k = 4)
> plt <- plot(nms.euc, display = "sites", type = "p")
> identify(plt, what = "sites")
> 
> I also forgot to mention that NMDS doesn't have species scores as such,
> so your example has nothing to do with species as far as I can see - did
> you omit something else?
> 
> vegan does allow you to get "species" scores via weighted averages
> (see ?scores.metaMDS ), e.g.
> 
> plt2 <- plot(nms.euc, display = c("species","sites"), type = "p", shrink
> = TRUE)
> identify(plt2, what = "species", cex = 0.8, col = "red")
> identify(plt2, what = "sites", cex = 0.8, col = "black")
> 
> Is this what you are looking for?
> 
> If you stick with labdsv:::nmds then take a look at plotid.nmds() which
> is a wrapper around identify for nmds objects, which appears to do the
> hard work for you. You could look at this to see how you could have gone
> about crafting the identify call yourself - the function is very simple.
> 
> Alternatively, do look at orditorp() in vegan, as it is a nice function
> that allows you to display only labels for species/sites that can fit:
> 
> plot(nms.euc, display = c("species","sites"), type = "n", shrink = TRUE)
> points(nms.euc, display = "sites", cex = 0.7)
> orditorp(nms.euc, display = "species", shrink = TRUE, col = "red")
> 
> The orditorp call will generate some warnings - they are harmless and
> have been fixed in t

[R] Question about AGNES by Rousseeuw et al. in the package "cluster": How many clusters?

2007-11-14 Thread AUDE MARIE PLONTZ
Dear all,

I am no stat wiz and I am just trying to use the AGNES algorithm at my
very modest level of statistical of understanding. I have difficulties
understanding the ouput from AGNES.  My question is: how to interpret
the output, especially how do you I know which cluster solution is the
best? In SPSS, an Agglomeration Schedule table is produced and I used to
look at the biggest jump between the error coefficients for each
agglomerative steps in order to find where to stop clustering. But with
the Agnes output I don't know what I should be looking at.

Thanks so much for your help,

Aude

Aude Plontz
Refugee Health Research Centre, Melbourne, Australia
Phone:
at home:
  +61 (0)3 9917 2134
at La Trobe University, Bundoora campus (every other Thursday):
  +61 (0)3 9479 3161 (desk, Martin Building, room 475,
consultation hours: 9.00am-11.00am, 1.00pm-5.00pm)
  or +61 (0)3 9479 5874 (Estelle Purchas, RHRC secretary)
at the Victorian Foundation for Survivors of Torture, Brunswick (on some
Wednesdays):
 +61 (0)3 9940 1567 (desk, Brunswick Business Incubator)
 or +61 (0)3 9388 0022 (Brunswick VFST reception)
Email:
[EMAIL PROTECTED]  (email account checked once daily,
except when on leave)
[EMAIL PROTECTED]

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


Re: [R] adding in missing values in a sequence

2007-11-14 Thread Andrew Hoskins
That seems to have solved it.

Thanks Mark.

On 15/11/2007, at 9:57 AM, Marc Schwartz wrote:

>
> On Thu, 2007-11-15 at 09:41 +1100, Andrew Hoskins wrote:
>> Hi,
>>
>> I have a data frame with two columns of data, one an indexing column
>> and the other a data column. My issue is, this data frame is
>> incomplete and there are missing lines.  I want to know how I can
>> find and add data into these missing lines.  See example below
>>
>> ## Example data
>>
>> data <- data.frame(index=c(1:4,6:10), data=c
>> (1.5,4.3,5.6,6.7,7.1,12.5,14.5,16.8,3.4))
>>
>> index data
>> 1 1  1.5
>> 2 2  4.3
>> 3 3  5.6
>> 4 4  6.7
>> 5 6  7.1
>> 6 7 12.5
>> 7 8 14.5
>> 8 9 16.8
>> 910  3.4
>>
>> ## note: index number 5 is missing
>>
>> ## What I want
>>
>> index data
>> 1 1  1.5
>> 2 2  4.3
>> 3 3  5.6
>> 4 4  6.7
>> 55 NA
>> 6 6  7.1
>> 7 7 12.5
>> 8 8 14.5
>> 9 9 16.8
>> 1010  3.4
>>
>> I'm running R2.6.0 on Mac OSX.
>
> How about this:
>
>> DF
>   index data
> 1 1  1.5
> 2 2  4.3
> 3 3  5.6
> 4 4  6.7
> 5 6  7.1
> 6 7 12.5
> 7 8 14.5
> 8 9 16.8
> 910  3.4
>
>
> DF.NEW <- data.frame(index = seq(max(DF$index)))
>
>> DF.NEW
>index
> 1  1
> 2  2
> 3  3
> 4  4
> 5  5
> 6  6
> 7  7
> 8  8
> 9  9
> 1010
>
>
> DF.NEW <- merge(DF.NEW, DF, all.x = TRUE)
>
>> DF.NEW
>index data
> 1  1  1.5
> 2  2  4.3
> 3  3  5.6
> 4  4  6.7
> 5  5   NA
> 6  6  7.1
> 7  7 12.5
> 8  8 14.5
> 9  9 16.8
> 1010  3.4
>
>
> See ?merge for more information.
>
> HTH,
>
> Marc Schwartz
>

Andrew Hoskins
PhD Candidate
Deakin University, Australia
Email: [EMAIL PROTECTED]




[[alternative HTML version deleted]]

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


Re: [R] adding in missing values in a sequence

2007-11-14 Thread Achim Zeileis
On Thu, 15 Nov 2007, Andrew Hoskins wrote:

> Hi,
>
> I have a data frame with two columns of data, one an indexing column
> and the other a data column. My issue is, this data frame is
> incomplete and there are missing lines.  I want to know how I can
> find and add data into these missing lines.  See example below

You could use a "zoo" series (from the "zoo" package) which provides
infrastructure for indexed observations. With your example data:

data <- data.frame(index = c(1:4, 6:10),
  data = c(1.5,4.3,5.6,6.7,7.1,12.5,14.5,16.8,3.4))

you can create a series

z <- zoo(data$data, data$index)

end extend it to the grid 1:10

z <- merge(zoo(,1:10), z)

which then has an NA at index 5. Then you could use linear interpolation
to replace that NA

na.approx(z)

or replace it with some other number

z[is.na(z)] <- 42

See vignette("zoo", package = "zoo") for more details.
Z

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


[R] Help with Bartlett's test on linear model

2007-11-14 Thread Stephanie Bernard

Hi all,
I would like to test the homoegeneity of variances between several 
linear model for some analysis of covariance. It seems that the 
Bartlett's test is a good test to use but I am having problem using with 
linear model and I cannot find any examples on the internet. There are 
some examples for comparisons of variances but not linear models.
If I take the hellung data set, which is the example in Dalgaard's book. 
I know var.test works fine but I want to learn how to use the Bartlett's 
test.

> hellung$glucose <- factor(hellung$glucose, labels=c("Yes","No"))
> attach(hellung)
> tethym.gluc <- hellung[glucose=="Yes",]
> tethym.nogluc <- hellung[glucose=="No",]
> lm.nogluc <- lm(log10(diameter)~log10(conc), data=tethym.nogluc)
> lm.gluc <- lm(log10(diameter)~log10(conc), data=tethym.gluc)
I guess I have two questions. 1) How to use bartlett.test with linear 
model (using the model above) and 2) how to test for homogeneity of 
variances of linear models when there are more than two groups.

Thanks,
Stephanie




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


Re: [R] Help with Bartlett's test on linear model

2007-11-14 Thread Emmanuel Charpentier
Stephanie Bernard a écrit :
> Hi all,
> I would like to test the homoegeneity of variances between several
> linear model for some analysis of covariance. It seems that the
> Bartlett's test is a good test to use but I am having problem using with
> linear model and I cannot find any examples on the internet. There are
> some examples for comparisons of variances but not linear models.
> If I take the hellung data set, which is the example in Dalgaard's book.
> I know var.test works fine but I want to learn how to use the Bartlett's
> test.
>> hellung$glucose <- factor(hellung$glucose, labels=c("Yes","No"))
>> attach(hellung)
>> tethym.gluc <- hellung[glucose=="Yes",]
>> tethym.nogluc <- hellung[glucose=="No",]
>> lm.nogluc <- lm(log10(diameter)~log10(conc), data=tethym.nogluc)
>> lm.gluc <- lm(log10(diameter)~log10(conc), data=tethym.gluc)
> I guess I have two questions. 1) How to use bartlett.test with linear
> model (using the model above)

Please look at ?bartlett.test. This generic has methods both for
formulae (formulas ?) and for objects inheriting from the class "lm".
Therefore you might call :
bartlett.test(nogluc)
and
bartlett.test(gluc)
or (equivalently, I think)
bartlett.test(log10(diameter)~log10(conc), data=tethym.nogluc)
etc ... ad nauseam.

and 2) how to test for homogeneity of
> variances of linear models when there are more than two groups.


As the "?bartlett.test" will tell you, this function "[p]erforms
Bartlett's test of the null that the variances in each of the groups
(samples) are the same", which is exactly what you aim at...

HTH

Emmanuel Charpentier

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


Re: [R] calc dist AB for B at levels of A

2007-11-14 Thread Emmanuel Charpentier
Melanie Ann Harsch a écrit :
> I have a question that I have been trying to figure out and I imagine 
> there is a very simple answer to it.
> 
> I am trying to use the distAB function in the clim.pact package
> I have a dataframe with 4 columns in which A is ref pt and B is site
> 1. longitude of A
> 2. Latitude of A
> 3. longitude of B
> 4. Latitude of B
> 
> Problem is the columns are unequal in length 
> column 1 and 2 are length=29 and 3 and 4 length=53

Curious representation... I'd rather use two dataframes (those are
different sets, of course ?

> What I want to do is create a matrix which are filled with values of 
> distAB for B (site) at all levels of A (reference points)
> 
> Any ideas?

OTTOMM :

D1<-df[1:29,1:3]
D2<-df[,3:4]
t<-outer(1:29,1:53,FUN=function(i,j)distAB(D1[i,1],D1[i,2],D2[j,1],D2[j,2]))

This is probably the most straightforward, and might be the more
efficient (no loops). However, you'll note that, if one set is a subset
of the other, each distance between points of the inner set is computed
*twice*. If the 'cost" of computing a distance is much higher than the
cost of an ifelse(), treating this special case might be worthy (but
that's not easy).

HTH,

Emmanuel Charpentier

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


Re: [R] Help with Bartlett's test on linear model

2007-11-14 Thread Ben Bolker



Stephanie Bernard wrote:
> 
> Hi all,
> I would like to test the homoegeneity of variances between several 
> linear model for some analysis of covariance. It seems that the 
> Bartlett's test is a good test to use but I am having problem using with 
> linear model and I cannot find any examples on the internet. There are 
> some examples for comparisons of variances but not linear models.
> If I take the hellung data set, which is the example in Dalgaard's book. 
> I know var.test works fine but I want to learn how to use the Bartlett's 
> test.
>  > hellung$glucose <- factor(hellung$glucose, labels=c("Yes","No"))
>  > attach(hellung)
>  > tethym.gluc <- hellung[glucose=="Yes",]
>  > tethym.nogluc <- hellung[glucose=="No",]
>  > lm.nogluc <- lm(log10(diameter)~log10(conc), data=tethym.nogluc)
>  > lm.gluc <- lm(log10(diameter)~log10(conc), data=tethym.gluc)
> I guess I have two questions. 1) How to use bartlett.test with linear 
> model (using the model above) and 2) how to test for homogeneity of 
> variances of linear models when there are more than two groups.
> 

1.

bartlett.test(list(lm.gluc,lm.nogluc))

2.

x = rnorm(40,sd=1)
y = rnorm(40,sd=2)
z = rnorm(40,sd=3)
bartlett.test(list(x,y,z))

-- 
View this message in context: 
http://www.nabble.com/Help-with-Bartlett%27s-test-on-linear-model-tf4807461.html#a13756626
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Log random number

2007-11-14 Thread Julian Burgos
Hi Tobias,

You'll have to explain what you mean with "log transformed values".  If 
your question is if it is possible to generate random numbers from a 
normal distribution (using rnorm()) and then getting their logarithm, 
then you can do that with the obvious caveat that the logarithm is not 
defined for negative numbers.

Julian


Wensui Liu wrote:
> dont think so, unless i miss something here.
> please do check the range of normal random number and the domain of
> log function.
> 
> On 11/14/07, Tobias Schlottmann <[EMAIL PROTECTED]> wrote:
>>
>>Dear R users,
>>
>>   Simply my question is that how it is possible to generate some random 
>> numbers using rnorm( ) but in log transformed values.
>>
>>   Thank you,
>>
>>   Tobias
>>
>>
>> -
>>
>> [[alternative HTML version deleted]]
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
> 
>

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


[R] adding in missing values in a sequence

2007-11-14 Thread Andrew Hoskins
Hi,

I have a data frame with two columns of data, one an indexing column  
and the other a data column. My issue is, this data frame is  
incomplete and there are missing lines.  I want to know how I can  
find and add data into these missing lines.  See example below

## Example data

data <- data.frame(index=c(1:4,6:10), data=c 
(1.5,4.3,5.6,6.7,7.1,12.5,14.5,16.8,3.4))

index data
1 1  1.5
2 2  4.3
3 3  5.6
4 4  6.7
5 6  7.1
6 7 12.5
7 8 14.5
8 9 16.8
910  3.4

## note: index number 5 is missing

## What I want

index data
1 1  1.5
2 2  4.3
3 3  5.6
4 4  6.7
5   5 NA
6 6  7.1
7 7 12.5
8 8 14.5
9 9 16.8
1010  3.4

I'm running R2.6.0 on Mac OSX.

Thanks in advance,

Andrew Hoskins
PhD Candidate
Deakin University, Australia
Email: [EMAIL PROTECTED]




[[alternative HTML version deleted]]

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


[R] Quantile Regression Question

2007-11-14 Thread Arti Mann
Hi,

Could you please explain what is non-positive fis error? I have been trying to 
use quantile regression (rq) procedure and I keep ending up with this error. I 
haven't been able to find an explanation for the same.

Best Regards,
Arti

Arti Mann
Ph.D. Student
Department of Information Systems
W.P. Carey School of Business
Arizona State University
Email : [EMAIL PROTECTED]

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


Re: [R] Quantile Regression Question

2007-11-14 Thread roger koenker
Package questions to package maintainers please.

"non-positive fis"  is not an error, it is a warning -- if the number  
of negative fis
is large relative to the sample  size then there is some reason to  
doubt the
plausibility of the specification of the model specified by the rq  
formula.  So you
can consider this _warning_ to be a diagnostic.


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 Nov 14, 2007, at 6:43 PM, Arti Mann wrote:

> Hi,
>
> Could you please explain what is non-positive fis error? I have  
> been trying to use quantile regression (rq) procedure and I keep  
> ending up with this error. I haven't been able to find an  
> explanation for the same.
>
> Best Regards,
> Arti
>
> Arti Mann
> Ph.D. Student
> Department of Information Systems
> W.P. Carey School of Business
> Arizona State University
> Email : [EMAIL PROTECTED]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting- 
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] adding in missing values in a sequence

2007-11-14 Thread Marc Schwartz

On Thu, 2007-11-15 at 09:41 +1100, Andrew Hoskins wrote:
> Hi,
> 
> I have a data frame with two columns of data, one an indexing column  
> and the other a data column. My issue is, this data frame is  
> incomplete and there are missing lines.  I want to know how I can  
> find and add data into these missing lines.  See example below
> 
> ## Example data
> 
> data <- data.frame(index=c(1:4,6:10), data=c 
> (1.5,4.3,5.6,6.7,7.1,12.5,14.5,16.8,3.4))
> 
> index data
> 1 1  1.5
> 2 2  4.3
> 3 3  5.6
> 4 4  6.7
> 5 6  7.1
> 6 7 12.5
> 7 8 14.5
> 8 9 16.8
> 910  3.4
> 
> ## note: index number 5 is missing
> 
> ## What I want
> 
> index data
> 1 1  1.5
> 2 2  4.3
> 3 3  5.6
> 4 4  6.7
> 5 5 NA
> 6 6  7.1
> 7 7 12.5
> 8 8 14.5
> 9 9 16.8
> 1010  3.4
> 
> I'm running R2.6.0 on Mac OSX.

How about this:

> DF
  index data
1 1  1.5
2 2  4.3
3 3  5.6
4 4  6.7
5 6  7.1
6 7 12.5
7 8 14.5
8 9 16.8
910  3.4


DF.NEW <- data.frame(index = seq(max(DF$index)))

> DF.NEW
   index
1  1
2  2
3  3
4  4
5  5
6  6
7  7
8  8
9  9
1010


DF.NEW <- merge(DF.NEW, DF, all.x = TRUE)

> DF.NEW
   index data
1  1  1.5
2  2  4.3
3  3  5.6
4  4  6.7
5  5   NA
6  6  7.1
7  7 12.5
8  8 14.5
9  9 16.8
1010  3.4


See ?merge for more information.

HTH,

Marc Schwartz

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


Re: [R] geotiff calculations

2007-11-14 Thread Dylan Beaudette
On Wednesday 14 November 2007, Monica Pisica wrote:
> Laura,
>
> As far as i know RGDAL and maybe proj4 will help you with geotiff, but it
> may be tricky if you want to preserve the geotiff projection. Your safest
> bet is to use a GIS software and ESRI ArcGIS can do everything you want
> without using R. To calculate differences, mean, and so on is a simple task
> in ArcGIS, no need for R. Evem more advanced statistics can be done in
> ArcGIS in their "expanded" calculator. You also can do some geostatistics,
> analysis of variances, correlations, "hot spots" and local outliers.
>
> I hope this helps,
>
> Monica

Or, if you like open source software, GRASS will do all of the above, has a 
great connector to R, and is free!

cheers,

Dylan

>
> -- Message: 110Date: Wed, 14 Nov 2007 09:22:52
> +0100From: "Laura Poggio" <[EMAIL PROTECTED]>Subject: [R] geotiff
> calculationsTo:
> [EMAIL PROTECTED]:[EMAIL PROTECTED]>Content-Type: text/plain Dear list,I have to compare two
> digital elevation models in raster format (geotiff).I then have to
> calculate the differences in altitude for each cell and makesome statistics
> (basic as mean, median, std, range but also more "advanced"as RMSE) on
> that.I do not know very much how to proceed:1) is it possible to import the
> geotiff in R? If so with which package? ifnot which is the best way to
> import such files?2) is it better to perform the calculations of the
> differences in a GISsoftware and then to use R only for statistical
> analysis? or it is better todo everything in R?3) is there any specific
> package for doing this kind of analysis? Thank you very much in advance
> Laura [[alternative HTML version del! eted]]  
> --
> _
> [[replacing trailing spam]]
>
>   [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html and provide commented, minimal,
> self-contained, reproducible code.



-- 
Dylan Beaudette
Soil Resource Laboratory
http://casoilresource.lawr.ucdavis.edu/
University of California at Davis
530.754.7341

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


[R] Moving to a Mac environment - quick question

2007-11-14 Thread David Kaplan

   Hi all,
   I'm moving to a Mac platform and have installed R and all seems fine.
   My  question  concerns moving my R objects over.  On my pc, which file
   contains my R objects and once I copy that, where do I copy it within the
   Mac?
   I hope that was clear.
   Thanks in advance,
   David
-- 
===
David Kaplan, Ph.D.
Professor
Department of Educational Psychology
University of Wisconsin - Madison
Educational Sciences, Room, 1061
1025 W. Johnson Street
Madison, WI 53706

email: [EMAIL PROTECTED]
homepage:
[2]http://www.education.wisc.edu/edpsych/facstaff/kaplan/kaplan.htm
Phone: 608-262-0836
===

References

   1. mailto:[EMAIL PROTECTED]
   2. http://www.education.wisc.edu/edpsych/facstaff/kaplan/kaplan.htm
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] reading tables from url

2007-11-14 Thread Duncan Temple Lang
-BEGIN PGP SIGNED MESSAGE-
Hash: SHA1

Hi Chris.

  Indeed, I cannot connect to that URL either.  So I did a bit of
digging and experimentation to find out whether one needed to
pass additional hidden options from the form or whether the problem was
more to do with how we connect.

It turns out that the script associated with NCBI leuks.cgi is being
fussy and wants you tell it the user agent that is performing the
request.  (Why the two behave differently is not clear after a very
brief look, but it is probably not worth pursuing.)

AFAIR, there is no way to tell R to include a UserAgent field in the
header of the request using url(), etc. although it did come up at one
point.

So here is an alternative. Use the RCurl package and this allows you
a great deal of control over the composition of the request and how
to read it back.

 getURL("http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi?dump=selected";,
 useragent = "curl", verbose = TRUE)

(The verbose is there to show the header of the request being made to
see the UserAgent field.)

One could do the same with sockets directly or use the
httpRequest for simple HTTP requests.

 D.

stubben wrote:
> I'm trying to read some web tables directly into R.  These are both  
> genome sequencing projects (eukaryotes and metagenomes) from NCBI and  
> look very similar;  however, only the first one works.
> 
> http://www.ncbi.nlm.nih.gov/genomes/leuks.cgi
> http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi
> 
> I added  ?dump=selected to the end of the url string to get a tab- 
> delimited file (which is what happens if you click the Save button on  
> either page).
> 
>  > options(internet.info=0)
> 
> ## this one works
> 
>  > x1<-url("http://www.ncbi.nlm.nih.gov/genomes/leuks.cgi? 
> dump=selected")
>  > read.delim(x1, skip=1, nrows=5)[,1:3]
> 
>X...Columns. ProjectID Organism.Name
> 120303 Acanthamoeba castellanii Neff  Protists
> 213657  Acyrthosiphon pisum LSR1   Animals
> 312434   Aedes aegypti Liverpool   Animals
> 412635 Ajellomyces capsulatus G186AR Fungi
> 512653  Ajellomyces capsulatus G217B Fungi
> 
> Warning messages:
> 1: connected to 'www.ncbi.nlm.nih.gov' on port 80. in: open.connection 
> (file, "r")
> 2: -> GET /genomes/leuks.cgi?dump=selected HTTP/1.0
> 
> Host: www.ncbi.nlm.nih.gov
> 
> Pragma: no-cache
> 
> in: open.connection(file, "r")
> 3: <- HTTP/1.1 200 OK in: open.connection(file, "r")
> 4: <- Date: Wed, 14 Nov 2007 18:03:29 GMT in: open.connection(file, "r")
> 5: <- Server: Apache in: open.connection(file, "r")
> 6: <- Content-Disposition: attachment; filename="untitle.txt" in:  
> open.connection(file, "r")
> 7: <- Content-Type: application/force-download in: open.connection 
> (file, "r")
> 8: <- Vary: Accept-Encoding in: open.connection(file, "r")
> 9: <- Connection: close in: open.connection(file, "r")
> 10: Code 200, content-type 'application/force-download' in:  
> open.connection(file, "r")
> 
> 
> ## this one fails to open a connection
> 
>  > x2<-url("http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi? 
> dump=selected")
>  > read.delim(x2, skip=1, nrows=5)[,1:3]
> 
> Error in open.connection(file, "r") : unable to open connection
> In addition: Warning messages:
> 1: connected to 'www.ncbi.nlm.nih.gov' on port 80. in: open.connection 
> (file, "r")
> 2: -> GET /genomes/lenvs.cgi?dump=selected HTTP/1.0
> 
> Host: www.ncbi.nlm.nih.gov
> 
> Pragma: no-cache
> 
> in: open.connection(file, "r")
> 3: <- HTTP/1.1 500 Internal Server Error in: open.connection(file, "r")
> 4: <- Date: Wed, 14 Nov 2007 18:04:26 GMT in: open.connection(file, "r")
> 5: <- Server: Apache in: open.connection(file, "r")
> 6: <- Content-Type: text/html; charset=ISO-8859-1 in: open.connection 
> (file, "r")
> 7: <- Vary: Accept-Encoding in: open.connection(file, "r")
> 8: <- Connection: close in: open.connection(file, "r")
> 9: Code 500, content-type 'text/html; charset=ISO-8859-1' in:  
> open.connection(file, "r")
> 10: cannot open: HTTP status was '500 Internal Server Error' in:  
> open.connection(file, "r")
> 
> Also, I can't even read lines from the main page.
> 
>  > readLines("http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi";, n=10)
> Error in file(con, "r") : unable to open connection
> ...
> ## now I'm just guessing...
>  > readLines("http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi";, n=10,  
> encoding="ISO-8859-1")
> Error in file(con, "r") : unable to open connection
> ...
> 
> 
> Download.file works fine, but I would like to avoid this if possible.
> 
>  > capabilities()[5]
> http/ftp
>  TRUE
> 
>  > download.file("http://www.ncbi.nlm.nih.gov/genomes/lenvs.cgi? 
> dump=selected", "lenvs.tab")
>  > read.delim("lenvs.tab", skip=1, nrows=5)[,1:3]
>X...Columns.  
> Parent.ProjectID ProjectID
> 11973313694   Global Ocean Sampling  
> Expedition Metagenome
> 220823   

[R] kalman filter estimation

2007-11-14 Thread adschai
Hi,

Following convention below:
y(t) = Ax(t)+Bu(t)+eps(t) # observation eq
x(t) = Cx(t-1)+Du(t)+eta(t) # state eq

I modified the following routine (which I copied from: 
http://www.stat.pitt.edu/stoffer/tsa2/Rcode/Kall.R) to accommodate u(t), an 
exogenous input to the system. 

for (i in 2:N){
 xp[[i]]=C%*%xf[[i-1]]
 Pp[[i]]=C%*%Pf[[i-1]]%*%t(C)+Q
   siginv=A[[i]]%*%Pp[[i]]%*%t(A[[i]])+R
 sig[[i]]=(t(siginv)+siginv)/2 # make sure sig is symmetric
   siginv=solve(sig[[i]])  # now siginv is sig[[i]]^{-1}
 K=Pp[[i]]%*%t(A[[i]])%*%siginv
 innov[[i]]=as.matrix(yobs[i,])-A[[i]]%*%xp[[i]]
 xf[[i]]=xp[[i]]+K%*%innov[[i]]
 Pf[[i]]=Pp[[i]]-K%*%A[[i]]%*%Pp[[i]]
 like= like + log(det(sig[[i]])) + t(innov[[i]])%*%siginv%*%innov[[i]]
 }
   like=0.5*like
   list(xp=xp,Pp=Pp,xf=xf,Pf=Pf,like=like,innov=innov,sig=sig,Kn=K)
}

I tried to fit my problem and observe that I got positive log likelihood mainly 
because the log of determinant of my variance matrix is largely negative. 
That's not good because they should be positive. Have anyone experience this 
kind of instability?

Also, I realize that I have about 800 sample points. The above routine when 
being plugged to optim becomes very slow. Could anyone share a faster way to 
compute kalman filter? 

And my last problem is, optim with my defined feasible space does not converge. 
I have about 20 variables that I need to identify using MLE method. Is there 
any other way that I can try out? I tried most of the methods available in 
optim already. They do not converge at all.. Thank you.

- adschai

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


Re: [R] Moving to a Mac environment - quick question

2007-11-14 Thread Prof Brian Ripley
On Wed, 14 Nov 2007, David Kaplan wrote:

>
>   Hi all,
>   I'm moving to a Mac platform and have installed R and all seems fine.
>   My  question  concerns moving my R objects over.  On my pc, which file
>   contains my R objects and once I copy that, where do I copy it within the
>   Mac?

It is the .RData file in the R working directory.  Use getwd() in an R 
session to see where that is.  (That applies whatever OS 'my pc' runs.)

You put it in the working directory on the Mac.  Same advice to find what 
it is (I use command-line R on the Mac, but presume you mean for the GUI 
console).

Please use R-sig-mac for Mac questions (in the future, not that this is 
one).

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

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


[R] generate combination set

2007-11-14 Thread aiminy
I have a set data={A,B,C,D,E,F,G}
I want to choose 2 letter from 8 letters, i.e. generate the combination set 
for choose 2 letters from 8 letters.
I want to get the liking:
combination set={AB,AC,AD,}
Does anyone konw how to do in R.

thanks,

Aimin

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


[R] what is the "right" way to obtain frequencies of vector values?

2007-11-14 Thread Vlad Skvortsov
Hi!

Let's say I have vector x with positive integer values ranging from 1 to 
N. I need to obtain another vector y of size N where y[i] contains the 
number of times value i occurs in x. It is in a sense similar to hist() 
(with appropriate number of breaks) or table() with numeric "factors".

Currenlty I use a custom function for that, but thought maybe there is a 
more "direct" way in R.

Thanks!

-- 
Vlad Skvortsov, [EMAIL PROTECTED], http://vss.73rus.com

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


Re: [R] differing scales for lattice plot panels

2007-11-14 Thread Frede Aakmann Tøgersen
See the scales argument of xyplot: e.g.

densityplot(~values | groups, pch=".", scales = list(y = "free")).


Best regards

Frede Aakmann Tøgersen
Scientist


UNIVERSITY OF AARHUS
Faculty of Agricultural Sciences
Dept. of Genetics and Biotechnology
Blichers Allé 20, P.O. BOX 50
DK-8830 Tjele

Phone:   +45 8999 1900
Direct:  +45 8999 1878

E-mail:  [EMAIL PROTECTED]
Web:   http://www.agrsci.org

This email may contain information that is confidential.
Any use or publication of this email without written permission from Faculty of 
Agricultural Sciences is not allowed.
If you are not the intended recipient, please notify Faculty of Agricultural 
Sciences immediately and delete this email.



 

> -Oprindelig meddelelse-
> Fra: [EMAIL PROTECTED] 
> [mailto:[EMAIL PROTECTED] På vegne af Joe Bloggs
> Sendt: 14. november 2007 03:24
> Til: [EMAIL PROTECTED]
> Emne: [R] differing scales for lattice plot panels
> 
> Hi,
>I want to compare densities across several groups with: 
> densityplot(~values | groups, pch=".")
> However one of the groups has an almost point-like distribution. 
> This forces the y-axis scale to be so large that the other 
> distribution curves are all flat in comparison and get lost 
> amongst the data dots.
> How can I produce a lattice plot with different y-axes 
> (preferably automatically scaled) for each panel?
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 

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


Re: [R] Generating these matrices going backwards

2007-11-14 Thread francogrex

Hi, thanks but the way you are doing it is to assign the values of N in the x
matrix, knowing from the example I have given where they are supposed to be.
While the assumption is, you ONLY have values of N and R and do NOT know
where they would be placed in the x and y matrix a-priori, but their
position has to be derived from only the (N and R) dataframe you have.


Chris Stubben wrote:
> 
> Maybe something like this?
> 
> N<-c(1,2,1,3)
> ## create empty matrix
> x<-diag(0,3)
> ## fill off diagonal
>  x[row(x)==col(x)+1]<-N[1:2]
> # fill 3rd column
>  x[1:2,3]<-N[3:4]
> 
> Or create a function and return both x and y
> 
> mat3<-function(N)
> {
>  x<-diag(0,3)
> # fill each element separately
> x[2,1]<-N[1]
> x[3,2]<-N[2]
> x[1,3]<-N[3]
> x[2,3]<-N[4]
> dimnames(x)<-list(c("D1", "D2", "D3"), c("E1", "E2", "E3"))
> y =sum(x) * x / (rowSums(x)%o%colSums(x)) 
> list(x=x,y=y)
> }
> mat3(N)
> 

-- 
View this message in context: 
http://www.nabble.com/Generating-these-matrices-going-backwards-tf4807447.html#a13763015
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Changing the text in the strips of lattice plots and y axis

2007-11-14 Thread Deepayan Sarkar
On Nov 14, 2007 9:06 PM, Judith Flores <[EMAIL PROTECTED]> wrote:
> Dear R-helpers,
>
>I am sorry for asking something I know has been
> asked before, I have tried different combinations in
> the strip function without success...
>
>
> I am using version 2.5.1 and work on a PC.
>
>I have  barcharts generated from the following
> formula:
>
> barchart(y1+y2+y3~x | g)
>
> I need to change the names of the variables y1,y2 or
> y3 that currently appear in the strips, I have two
> strips per panel, one corresponds to the name of the
> variable (y1 or y2) and the second strip to the value
> of the conditional variable (g has two levels, this
> means I have 6 panels total).
>
>   Instead of y1,y2 and y3 I would like the strips to
> read 'name of variable y1', 'name of variable y2',
> 'name of variable y3'; of course these are not the
> real names I want to assign to those variables.

The easiest way is to set the dimnames; e.g., something like

p <- barchart(...)
dimanames(p)[[1]] <- c("A", "B", "C")
p

>   And my second question is regarding the the limits
> of the y axis. If I setup scales to relation='free', I
> obtain 6 different scales, which is expected. But I
> need the panels that correspond to each y variable
> have the same scale. I tried something like this:
>
> prepanel=function(x,y,...) {ylim=list(min(y),max(y)) }
>
>but didn't work...
>
> What do I need to do?

I'm not sure what you are looking for. Please provide a reproducible example.

-Deepayan

>
> Thank you very much in advance for your help,
>
> Judith
>
>
>   
> 
> Never miss a thing.  Make Yahoo your home page.
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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