Re: [R] package update

2011-05-10 Thread eric
ok, how do I do root permissions from the RStudio GUI ?  What is the specific
command that allows me to do that when I type ...update.packages()  ?

Also, why are the packages installing in the first place if I can't write to
that location ?

Currently running linux ubuntu 10.04 by the way.

--
View this message in context: 
http://r.789695.n4.nabble.com/package-update-tp3507479p3510977.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] spdep::poly2nb problem

2011-05-10 Thread Roger Bivand
This message should have had a more informative subject, and, given its
specialist character, been sent to R-sig-geo, as is recommended on the
Spatial Task View. From initial trials, the file you are using includes far
too much boundary detail, and appears to have been assembled from disparate
state/county maps (Florida is added late, not in sequence). With a malformed
input map intended for display only, you shouldn't be surprised that the
topology is wrong. I'll investigate whether there is a code issue, but the
real problem is the choice of a very inappropriate input map of an object
size of 64Mb for 56 objects, where you do not need all of the coastal
detail.

Roger


Philip A. Viton wrote:
> 
> I downloaded a US states+territories shapefile from
> http://www.weather.gov/geodata/catalog/national/html/us-state.htm
> and unzipped it into my working directory.
> Then in R (version 12.2.0 for ms-win32):
> 
> library(rgdal) # Version: 0.6-28
> library(spdep) # Version: 0.5-24 (loads sp: Version: 0.9-72)
> 
> # read in the data
> states<-readOGR(".","s_01de10")
> summary(states)
> 
> # now create a neighbors list
> states_nb<-poly2nb(states)
> # create a 0/1 matrix of the neighbors
> nb_mat<-nb2mat(states_nb,style="B",zero.policy=TRUE)
> # add names
> nms<-as.character(states$NAMES)
> # set the matrix names
> dimnames(nb_mat)<-list(nms,nms)
> 
> # now see which states have no neighbors
> rowSums(nb_mat)   
> # it turns out that Florida   seems to have zero neighbors.
> 
> I then repeated it with
> states_nb<-poly2nb(states,useC=FALSE)
> 
> and this time -- the computation takes a lot longer, of course
> --- Florida has 2 neighbors (correctly: AL and GA).
> 
> 
> Can anyone tell me what's going wrong here?
> 
> 
> Thanks!
> 
> 
> 
> Philip A. Viton
> City Planning, Ohio State University
> 275 West Woodruff Avenue, Columbus OH 43210
> vito...@osu.edu
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 


-
Roger Bivand
Economic Geography Section
Department of Economics
Norwegian School of Economics and Business Administration
Helleveien 30
N-5045 Bergen, Norway

--
View this message in context: 
http://r.789695.n4.nabble.com/spdep-poly2nb-problem-tp3506042p3511103.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] Scale time series in a way that 90% of the data is in the -0.-9/ +0.9 range

2011-05-10 Thread Mr Copper
Hello,

please excuse my ignorance as i am new to R and this might seem trivial to
you

How can i scale my time series in a way that 90% of the data is in the
-0.-9/ +0.9 range?

I think i first have to select the 90% closest to the mean and then
calculate my scaling parameter on it, then scale the whole series.

I think i could use outlier() from outliers package to identify the outer
10%, extract them and then do my scaling, but i assume there is a more
elegant /efficient way to do so.

Does anyone have a hint for me?

Thanks in advance,

Warren

[[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] Creating binary variable depending on strings of two dataframes

2011-05-10 Thread noxyp...@gmail.com
On Fri, May 6, 2011 at 7:41 PM, David Winsemius  wrote:
>
> On May 6, 2011, at 11:35 AM, Pete Pete wrote:
>
>>
>> Gabor Grothendieck wrote:
>>>
>>> On Tue, Dec 7, 2010 at 11:30 AM, Pete Pete 
>>> wrote:

 Hi,
 consider the following two dataframes:
 x1=c("232","3454","3455","342","13")
 x2=c("1","1","1","0","0")
 data1=data.frame(x1,x2)

 y1=c("232","232","3454","3454","3455","342","13","13","13","13")
 y2=c("E1","F3","F5","E1","E2","H4","F8","G3","E1","H2")
 data2=data.frame(y1,y2)

 I need a new column in dataframe data1 (x3), which is either 0 or 1
 depending if the value "E1" in y2 of data2 is true while x1=y1. The
 result
 of data1 should look like this:
  x1     x2 x3
 1 232   1   1
 2 3454 1   1
 3 3455 1   0
 4 342   0   0
 5 13     0   1

 I think a SQL command could help me but I am too inexperienced with it
 to
 get there.

>>>
>>> Try this:
>>>
 library(sqldf)
 sqldf("select x1, x2, max(y2 = 'E1') x3 from data1 d1 left join data2 d2
 on (x1 = y1) group by x1, x2 order by d1.rowid")
>>>
>>>   x1 x2 x3
>>> 1  232  1  1
>>> 2 3454  1  1
>>> 3 3455  1  0
>>> 4  342  0  0
>>> 5   13  0  1
>>>
>>>
> snipped Gabor's sig
>>
>> That works pretty cool but I need to automate this a bit more. Consider
>> the
>> following example:
>>
>> list1=c("A01","B04","A64","G84","F19")
>>
>> x1=c("232","3454","3455","342","13")
>> x2=c("1","1","1","0","0")
>> data1=data.frame(x1,x2)
>>
>> y1=c("232","232","3454","3454","3455","342","13","13","13","13")
>> y2=c("E13","B04","F19","A64","E22","H44","F68","G84","F19","A01")
>> data2=data.frame(y1,y2)
>>
>> I want now to creat a loop, which creates for every value in list1 a new
>> binary variable in data1. Result should look like:
>> x1      x2      A01     B04     A64     G84     F19
>> 232     1       0       1       0       0       0
>> 3454    1       0       0       1       0       1
>> 3455    1       0       0       0       0       0
>> 342     0       0       0       0       0       0
>> 13      0       1       0       0       1       1
>
> Loops!?! We don't nee no steenking loops!
>
>> xtb <-  with(data2, table(y1,y2))
>> cbind(data1, xtb[match(data1$x1, rownames(xtb)), ] )
>       x1 x2 A01 A64 B04 E13 E22 F19 F68 G84 H44
> 232   232  1   0   0   1   1   0   0   0   0   0
> 3454 3454  1   0   1   0   0   0   1   0   0   0
> 3455 3455  1   0   0   0   0   1   0   0   0   0
> 342   342  0   0   0   0   0   0   0   0   0   1
> 13     13  0   1   0   0   0   0   1   1   1   0
>
> I am guessing that you were to ... er, busy? ... to complete the table?
>
> --
>
> David Winsemius, MD
> West Hartford, CT
>
>

Thanks a lot! Pretty simple. I am so much used to SQLDF right now.

So how would you handle more complicated strings like that:
y1=c("232","232", "232", "3454","3454","3455","342","13","13","13","13")
y2=c("E13","B04 A01 F19","B04","F19","A64 G84 A05","E22","H44
C35","F68","G84","F19","A01")
data2=data.frame(y1,y2)

Where you want to extract for instance all "A01" from the strings?

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

2011-05-10 Thread Komine
Hi, 
I try always to do my nomogram. I change little my code in order to definite
the range of my linear predictor. 
But I do not succeed to have 0 until 1. I tried this code.   
library(rms)
d <- datadist(Fire)
options(datadist='d')
Fire$Ignition=c(0,1,by=0.1)
f<-lrm(Ignition~FMC+Charge,data=Fire)
print(f)
summary(f,FMC=c(0,506.40),Charge=c(45,646))
plot(nomogram(f, interact=list(FMC,Charge))) 

The result is: 

> print(f)

Logistic Regression Model

lrm(formula = Ignition ~ FMC + Charge, data = Fire)

 Model Likelihood DiscriminationRank Discrim.
Ratio TestIndexes  Indexes   

Obs   231LR chi2  0.04R2   0.000C   0.501
 0 77d.f.2g0.026Dxy 0.002
 0.1   77Pr(> chi2) 0.9796gr   1.026gamma   0.003
 1 77 gp   0.006tau-a   0.002
max |deriv| 7e-13 Brier0.222 


   CoefS.E.   Wald Z Pr(>|Z|)
y>=0.1  0.6832 0.3018  2.26  0.0236  
y>=1   -0.7033 0.3020 -2.33  0.0199  
FMC 0.0002 0.0011  0.18  0.8601  
Charge -0.0001 0.0009 -0.07  0.9424  


> summary(f,FMC=c(0,506.40),Charge=c(45,646)) 
 Effects  Response : Ignition 

 Factor  Low High  Diff. Effect S.E. Lower 0.95 Upper 0.95
 FMC  0  506.4 506.4  0.10  0.58 -1.04  1.24  
  Odds Ratio  0  506.4 506.4  1.11NA  0.36  3.45  
 Charge  45  646.0 601.0 -0.04  0.55 -1.11  1.03  
  Odds Ratio 45  646.0 601.0  0.96NA  0.33  2.81  
> plot(nomogram(f, interact=list(FMC,Charge)))   
http://r.789695.n4.nabble.com/file/n3511218/nomogram_2.jpeg 
Please can I do to have linear predictor between 0 to 1. 
Thanks very much for your Help
Komine

--
View this message in context: 
http://r.789695.n4.nabble.com/Draw-a-nomogram-after-glm-tp3498144p3511218.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Saving multiple 3x3 TIFF graphics inside a loop

2011-05-10 Thread Uwe Ligges



On 10.05.2011 03:36, John S wrote:

Dear Friends,

I have been trying to save multiple 3x3 (mfrow=c(3,3) graphics inside a loop
using tiff figure format (not using PDF or savePlot functions) with no
success. Could you please help?

Here is a simplified example code:



dat=data.frame (ID=rep(1:10,each=10),IDV=rep(seq(1:10),times=10))

dat$DV<- with(dat, 50+15*IDV)

dat=dat[order(dat$ID,dat$IDV),]

for(i in 1:10){

dt1 = dat[dat$ID==i,]

pagei=as.integer(ceiling(i/9))

#This does not work

if(i==1) {tiff(file=paste("Outfile",pagei,".tiff",sep="")

  ,width=17.15,height=17.15,units="cm", res=1200, pointsize=10,
compression = "lzw")}


You may not be able to allocate space for this huge tiffs (at least I am 
not able to do so under Windows).





   if(i==10) {tiff(file=paste("Outfile",pagei,".tiff",sep="")

  ,width=17.15,height=17.15,units="cm", res=1200, pointsize=10,
compression = "lzw")}

par(mar=c(5,5,4,1),mfrow=c(3,3))

plot(dt1[, "IDV"],dt1[, "DV"], type="n")

lines(dt1[, "IDV"], dt1[, "DV"], type="b", col=1,lwd=3,pch=1)

title(paste(i),col=8,cex=.8)

  if(i==9) dev.off ## close file


Your forgot parentheses after dev.off.

Uwe Ligges



}

dev.off()

  The idea is to save 2 files with 3x3 layouts but I only get one subject in
Outfile1.



OS: Windows 7 64 bits

R version 2.13.0

Many Thanks,

John

[[alternative HTML version deleted]]

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


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


Re: [R] Uniroot - error

2011-05-10 Thread peter dalgaard

On May 6, 2011, at 16:37 , Joshua Wiley wrote:

> On Fri, May 6, 2011 at 6:33 AM, CarJabo  wrote:
>> sorry I am not asking someone to do my homework, as I have finished all the
>> procedure. I am just wondering why this technical error occurs, so I can fix
>> it myself.
> 

> If whatever institution you are taking this at simply gives
> assignments, grades them and penalizes for plagiarism without having
> any instructors or teachers, I recommend moving to an institution
> where classes are taught by someone who can answer questions etc.  As
> Dr. Murdoc (and the posting guide) said and say, R-help is not for
> homework.

The occasional hint might slip though the iron-cladding once in a while though.
In this case, learning how to use options(error=recover) might be relevant. 
Turn it on, let the error occur and do the post mortem. The culprit is 
obviously that f(2) isn't positive. Unless F has turned out all-zero or has 
negative elements, you need a larger value of 2, otherwise, find out why F is 
weird.

-- 
Peter Dalgaard
Center for Statistics, Copenhagen Business School
Solbjerg Plads 3, 2000 Frederiksberg, Denmark
Phone: (+45)38153501
Email: pd@cbs.dk  Priv: pda...@gmail.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] Question Using Old R Packages

2011-05-10 Thread Duncan Murdoch

On 11-05-09 8:59 PM, c.hulmelowe wrote:

In case anyone else is having similar problems, it turns out it may have been
an issue with the interface I was using (R Studio). When I close R Studio
and tried it again through R it installed without the error message.


You should report this as a serious bug to the R Studio developers.  It 
sounds as though it is looking in the wrong location for packages to 
install.  You were "lucky" in that rgl wasn't in that location, but in 
other cases it may be finding the wrong version.


Duncan Murdoch

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


Re: [R] Installing package

2011-05-10 Thread Uwe Ligges



On 10.05.2011 06:34, Hurr wrote:

Using install packages pull-down, I chose a mirror and then package rJava.
Then it said that it cannot remove the rJava package already installed.
Then when I typed library(rJava) it said there was no rjava package.
Isn't that rather contradictory for R?



Start out without loading rJava and then try to reinstall. You probably 
tried to update when rJava was loaded and not al files could be removed 
due to file locking.


Uwe Ligges





--
View this message in context: 
http://r.789695.n4.nabble.com/Installing-package-tp3510987p3510987.html
Sent from the R help mailing list archive at Nabble.com.

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


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


[R] Join and Meet on Lattice

2011-05-10 Thread Vahid
Hi,

I am using a package to calculate Join and Meet on lattice. 
I have written a function that generates partitions and then try to see if 
their 
meet is the power set.
When I give values to this function greater than lets say 6 (n=k=6) and delta 
=0 
, it slows down.
I believe it is something due to commands that I use. But I can't see it. 
Is there any one who can help me?

vahid


* The code:

#library(partitions)
library(clue)

# First clear all
rm(list=ls(all=TRUE)) 


# n: number of balls
# k: number of urns
# this is for the case of no difficulty
# delta shows the difference between number of urns. delta=0 shows that urns 
are 
the same in both cases.
# 

urn <- function(n,k,delta){

# trivial partitions or power set: all elements are in one partition. i.e. (1 1 
1 1 1 ) for n = 5
trivial <- as.numeric(cbind(rep(1,n)))
t <- as.cl_class_ids(trivial)

# write Trivial in the matrice form (hard partition)
# set counter and number of rounds.
counter <- 0
rounds  <- 10

for(i in 1:rounds){
x<- as.numeric(cbind(sample.int( k , n , replace = TRUE, prob = NULL)))
y<- as.numeric(cbind(sample.int( k + delta , n , replace = TRUE, prob = NULL)))


# A matrice which shows each object belongs to which class/partition. This is a 
Hard partition; membership id shown by {0,1}
A <- as.cl_partition(x)
B <- as.cl_partition(y)


#join of them: The out put is a vector of the class ids 
m <- cl_join(A,B)


# check if the join of two partition is the trivial partition 
if(m==t){counter<-counter+1}
}
probability <- counter/rounds

return(probability)}
[[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] Cumulative dose in survival models

2011-05-10 Thread klaus

Dear List,
 
I want to fit a CPH survival model with a time-dependent cumulative exposure 
variable. More specifically, I have exposure measurements for all cohort 
members at discrete times (yearly measurements) during follow-up.
 
I suppose that when estimating the model parameters the most correct would be 
to have the cumulative exposure for everybody at risk at each death time, as I 
cannot assume the cumulative exposure to be constant between the yearly 
measurements . Thus, I need to calculate (interpolate) cumulative exposure for 
everybody at risk at each death time and pass this to the estimation routine.
 
Does anyone know how to do this in R? Any help would be much appreciated.
 
Klaus Andersen

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

2011-05-10 Thread Frank Harrell
by is not an argument to c( )
datadist needs to know the final coding of variables so run it after doing
any re-codes.
You are using the wrong syntax for interact=
Please take a while to read the full documentation before responding
further, and read my earlier posts.
Frank


Komine wrote:
> 
> Hi, 
> I try always to do my nomogram. I change little my code in order to
> definite the range of my linear predictor. 
> But I do not succeed to have 0 until 1. I tried this code.   
> library(rms)
> d <- datadist(Fire)
> options(datadist='d')
> Fire$Ignition=c(0,1,by=0.1)
> f<-lrm(Ignition~FMC+Charge,data=Fire)
> print(f)
> summary(f,FMC=c(0,506.40),Charge=c(45,646))
> plot(nomogram(f, interact=list(FMC,Charge))) 
> 
> The result is: 
> 
>> print(f)
> 
> Logistic Regression Model
> 
> lrm(formula = Ignition ~ FMC + Charge, data = Fire)
> 
>  Model Likelihood DiscriminationRank Discrim.
> Ratio TestIndexes  Indexes   
> 
> Obs   231LR chi2  0.04R2   0.000C   0.501
>  0 77d.f.2g0.026Dxy 0.002
>  0.1   77Pr(> chi2) 0.9796gr   1.026gamma   0.003
>  1 77 gp   0.006tau-a   0.002
> max |deriv| 7e-13 Brier0.222 
> 
> 
>CoefS.E.   Wald Z Pr(>|Z|)
> y>=0.1  0.6832 0.3018  2.26  0.0236  
> y>=1   -0.7033 0.3020 -2.33  0.0199  
> FMC 0.0002 0.0011  0.18  0.8601  
> Charge -0.0001 0.0009 -0.07  0.9424  
> 
> 
>> summary(f,FMC=c(0,506.40),Charge=c(45,646)) 
>  Effects  Response : Ignition 
> 
>  Factor  Low High  Diff. Effect S.E. Lower 0.95 Upper 0.95
>  FMC  0  506.4 506.4  0.10  0.58 -1.04  1.24  
>   Odds Ratio  0  506.4 506.4  1.11NA  0.36  3.45  
>  Charge  45  646.0 601.0 -0.04  0.55 -1.11  1.03  
>   Odds Ratio 45  646.0 601.0  0.96NA  0.33  2.81  
>> plot(nomogram(f, interact=list(FMC,Charge)))   
> http://r.789695.n4.nabble.com/file/n3511218/nomogram_2.jpeg 
> Please can I do to have linear predictor between 0 to 1. 
> Thanks very much for your Help
> Komine
> 


-
Frank Harrell
Department of Biostatistics, Vanderbilt University
--
View this message in context: 
http://r.789695.n4.nabble.com/Draw-a-nomogram-after-glm-tp3498144p3511758.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Creating Observation ID

2011-05-10 Thread Robert Baer
Thanks Bill and William.  Both solutions work, and both have introduced me 
to whole new R functions which I'm still wrapping my head around.


What a great resource this list is.

Thanks again,
Rob

--
Robert W. Baer, Ph.D.
Professor of Physiology
Kirksville College of Osteopathic Medicine
A. T. Still University of Health Sciences
800 W. Jefferson St.
Kirksville, MO 63501
660-626-2322
FAX 660-626-2965


--
From: "William Dunlap" 
Sent: Monday, May 09, 2011 5:17 PM
To: "Robert Baer" ; 
Subject: RE: [R] Creating Observation ID


Does the following work for you?

> df2 <- transform(df, ObsID=ave(rep(0,length(Group)), Group,
FUN=seq_along))
> head(df2)
   Value Group ObsID
1 -0.0025132 B 1
2 -1.2456156 A 1
3 -2.0531704 B 2
4  1.5861770 B 3
5  0.1900908 A 2
6  0.7197067 B 4

Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com


-Original Message-
From: r-help-boun...@r-project.org
[mailto:r-help-boun...@r-project.org] On Behalf Of Robert Baer
Sent: Monday, May 09, 2011 2:22 PM
To: R-help@r-project.org
Subject: [R] Creating Observation ID

If I have a data frame something like:
Value=rnorm(30)
Group = sample(c('A','B','C'), 30, replace=TRUE)
df = data.frame(Value, Group)

It seems like it should be simple to create an 'ObsID' column
which indicates the observation order of each Value within
each of the 3 groups.  Somehow, I can't quite see how to do
it without manually sub-setting the parent data frame and
then putting it back together again.

Anyone able to get me started on a cleaner (more R-like) approach?

Thanks,

Rob

--
Robert W. Baer, Ph.D.
Professor of Physiology
Kirksville College of Osteopathic Medicine
A. T. Still University of Health Sciences
800 W. Jefferson St.
Kirksville, MO 63501
660-626-2322
FAX 660-626-2965

[[alternative HTML version deleted]]

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





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


Re: [R] Creating binary variable depending on strings of two dataframes

2011-05-10 Thread David Winsemius


On May 10, 2011, at 3:18 AM, noxyp...@gmail.com wrote:

On Fri, May 6, 2011 at 7:41 PM, David Winsemius > wrote:


On May 6, 2011, at 11:35 AM, Pete Pete wrote:



Gabor Grothendieck wrote:


On Tue, Dec 7, 2010 at 11:30 AM, Pete Pete  


wrote:


Hi,
consider the following two dataframes:
x1=c("232","3454","3455","342","13")
x2=c("1","1","1","0","0")
data1=data.frame(x1,x2)

y1=c("232","232","3454","3454","3455","342","13","13","13","13")
y2=c("E1","F3","F5","E1","E2","H4","F8","G3","E1","H2")
data2=data.frame(y1,y2)

I need a new column in dataframe data1 (x3), which is either 0  
or 1
depending if the value "E1" in y2 of data2 is true while x1=y1.  
The

result
of data1 should look like this:
 x1 x2 x3
1 232   1   1
2 3454 1   1
3 3455 1   0
4 342   0   0
5 13 0   1

I think a SQL command could help me but I am too inexperienced  
with it

to
get there.



Try this:


library(sqldf)
sqldf("select x1, x2, max(y2 = 'E1') x3 from data1 d1 left join  
data2 d2

on (x1 = y1) group by x1, x2 order by d1.rowid")


  x1 x2 x3
1  232  1  1
2 3454  1  1
3 3455  1  0
4  342  0  0
5   13  0  1



snipped Gabor's sig


That works pretty cool but I need to automate this a bit more.  
Consider

the
following example:

list1=c("A01","B04","A64","G84","F19")

x1=c("232","3454","3455","342","13")
x2=c("1","1","1","0","0")
data1=data.frame(x1,x2)

y1=c("232","232","3454","3454","3455","342","13","13","13","13")
y2=c("E13","B04","F19","A64","E22","H44","F68","G84","F19","A01")
data2=data.frame(y1,y2)

I want now to creat a loop, which creates for every value in list1  
a new

binary variable in data1. Result should look like:
x1  x2  A01 B04 A64 G84 F19
232 1   0   1   0   0   0
34541   0   0   1   0   1
34551   0   0   0   0   0
342 0   0   0   0   0   0
13  0   1   0   0   1   1


Loops!?! We don't nee no steenking loops!


xtb <-  with(data2, table(y1,y2))
cbind(data1, xtb[match(data1$x1, rownames(xtb)), ] )

  x1 x2 A01 A64 B04 E13 E22 F19 F68 G84 H44
232   232  1   0   0   1   1   0   0   0   0   0
3454 3454  1   0   1   0   0   0   1   0   0   0
3455 3455  1   0   0   0   0   1   0   0   0   0
342   342  0   0   0   0   0   0   0   0   0   1
13 13  0   1   0   0   0   0   1   1   1   0

I am guessing that you were to ... er, busy? ... to complete the  
table?


--

David Winsemius, MD
West Hartford, CT




Thanks a lot! Pretty simple. I am so much used to SQLDF right now.

So how would you handle more complicated strings like that:
y1=c("232","232", "232",  
"3454","3454","3455","342","13","13","13","13")

y2=c("E13","B04 A01 F19","B04","F19","A64 G84 A05","E22","H44
C35","F68","G84","F19","A01")
data2=data.frame(y1,y2)

Where you want to extract for instance all "A01" from the strings?


I think you need either to explain what you want in more words of the  
English language or to offer an example of the desired output. I  
suspect you did not want something as simple as this:


> A01.instances <- grep("A01" , data2$y2)
> A01.instances
[1]  2 11
> data2[A01.instances, ]
y1  y2
2  232 B04 A01 F19
11  13 A01

Or maybe you did?

--
David Winsemius, MD
West Hartford, CT

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Saving multiple 3x3 TIFF graphics inside a loop

2011-05-10 Thread John S
Thanks Uwe . The parentheses were a typo. Please note that if I have only 9
subjects in my data then the following works (creates 493 KB data file)



#windows()

{tiff(file=paste("Outfile",1,".tiff",sep=""),width=17.15

  ,height=17.15,units="cm", res=1200, pointsize=10, compression =
"lzw")}

par(mar=c(5,5,4,1),mfrow=c(3,3))

for(i in 1:9){

   dt1 = dat[dat$ID==i,]

   #pagei=as.integer(ceiling(i/9))

   #This does not work

   #if(i==1) {tiff(file=paste("Outfile",pagei,".tiff",sep=""),width=17.15

   ##,height=17.15,units="cm", res=1200, pointsize=10, compression = "lzw")}


   #par(mar=c(5,5,4,1),mfrow=c(3,3))

   plot(dt1[, "IDV"],dt1[, "DV"], type="n",ylim=range(dt1[, "DV"]),
xlim=range(dt1[, "IDV"]))

   lines(dt1[, "IDV"], dt1[, "DV"], type="b", col=1,lwd=3,pch=1)

   title(paste(i),col=8,cex=.8)

   #if(i==5) dev.off ## close file

   ## open new file

   #if(i==5) {tiff(file=paste("Outfile",pagei,".tiff",sep=""),width=17.15

   ##,height=17.15,units="cm", res=1200, pointsize=10, compression = "lzw")}

   #if(i==9) dev.off () ## close file

   #if(i==9) {tiff(file=paste("Outfile",pagei,".tiff",sep=""),width=17.15

   ##,height=17.15,units="cm", res=1200, pointsize=10, compression = "lzw")}

}

dev.off()





My workaround around was to break the dataset into small data each with 9
subjects but if I change the mfrow then I have to change no. of subjects.
Any smart suggestion to get around this? I tried to open and close devices
based on iterations (i) but it does not work for me




2011/5/10 Uwe Ligges 

>
>
> On 10.05.2011 03:36, John S wrote:
>
>> Dear Friends,
>>
>> I have been trying to save multiple 3x3 (mfrow=c(3,3) graphics inside a
>> loop
>> using tiff figure format (not using PDF or savePlot functions) with no
>> success. Could you please help?
>>
>> Here is a simplified example code:
>>
>>
>>
>> dat=data.frame (ID=rep(1:10,each=10),IDV=rep(seq(1:10),times=10))
>>
>> dat$DV<- with(dat, 50+15*IDV)
>>
>> dat=dat[order(dat$ID,dat$IDV),]
>>
>> for(i in 1:10){
>>
>>dt1 = dat[dat$ID==i,]
>>
>>pagei=as.integer(ceiling(i/9))
>>
>>#This does not work
>>
>>if(i==1) {tiff(file=paste("Outfile",pagei,".tiff",sep="")
>>
>>  ,width=17.15,height=17.15,units="cm", res=1200, pointsize=10,
>> compression = "lzw")}
>>
>
> You may not be able to allocate space for this huge tiffs (at least I am
> not able to do so under Windows).
>
>
>
>
>   if(i==10) {tiff(file=paste("Outfile",pagei,".tiff",sep="")
>>
>>  ,width=17.15,height=17.15,units="cm", res=1200, pointsize=10,
>> compression = "lzw")}
>>
>>par(mar=c(5,5,4,1),mfrow=c(3,3))
>>
>>plot(dt1[, "IDV"],dt1[, "DV"], type="n")
>>
>>lines(dt1[, "IDV"], dt1[, "DV"], type="b", col=1,lwd=3,pch=1)
>>
>>title(paste(i),col=8,cex=.8)
>>
>>  if(i==9) dev.off ## close file
>>
>
> Your forgot parentheses after dev.off.
>
> Uwe Ligges
>
>
>  }
>>
>> dev.off()
>>
>>  The idea is to save 2 files with 3x3 layouts but I only get one subject
>> in
>> Outfile1.
>>
>>
>>
>> OS: Windows 7 64 bits
>>
>> R version 2.13.0
>>
>> Many Thanks,
>>
>> John
>>
>>[[alternative HTML version deleted]]
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>

[[alternative HTML version deleted]]

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


[R] Plotting confidence Intervals

2011-05-10 Thread Shubham Saxena

Hey guys,
 
Suppose i have
 
X ~ LN(0,1)
 
A = rlnorm(5)
 
How can I plot the 95% Confidence Interval of the ecdf of A in R?
 
Thank You 
[[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] Creating binary variable depending on strings of two dataframes

2011-05-10 Thread noxyp...@gmail.com
On Tue, May 10, 2011 at 3:09 PM, David Winsemius 
wrote:
>
> On May 10, 2011, at 3:18 AM, noxyp...@gmail.com wrote:
>
>> On Fri, May 6, 2011 at 7:41 PM, David Winsemius 
>> wrote:
>>>
>>> On May 6, 2011, at 11:35 AM, Pete Pete wrote:
>>>

 Gabor Grothendieck wrote:
>
> On Tue, Dec 7, 2010 at 11:30 AM, Pete Pete 
> wrote:
>>
>> Hi,
>> consider the following two dataframes:
>> x1=c("232","3454","3455","342","13")
>> x2=c("1","1","1","0","0")
>> data1=data.frame(x1,x2)
>>
>> y1=c("232","232","3454","3454","3455","342","13","13","13","13")
>> y2=c("E1","F3","F5","E1","E2","H4","F8","G3","E1","H2")
>> data2=data.frame(y1,y2)
>>
>> I need a new column in dataframe data1 (x3), which is either 0 or 1
>> depending if the value "E1" in y2 of data2 is true while x1=y1. The
>> result
>> of data1 should look like this:
>>  x1 x2 x3
>> 1 232   1   1
>> 2 3454 1   1
>> 3 3455 1   0
>> 4 342   0   0
>> 5 13 0   1
>>
>> I think a SQL command could help me but I am too inexperienced with
it
>> to
>> get there.
>>
>
> Try this:
>
>> library(sqldf)
>> sqldf("select x1, x2, max(y2 = 'E1') x3 from data1 d1 left join data2
>> d2
>> on (x1 = y1) group by x1, x2 order by d1.rowid")
>
>  x1 x2 x3
> 1  232  1  1
> 2 3454  1  1
> 3 3455  1  0
> 4  342  0  0
> 5   13  0  1
>
>
>>> snipped Gabor's sig

 That works pretty cool but I need to automate this a bit more. Consider
 the
 following example:

 list1=c("A01","B04","A64","G84","F19")

 x1=c("232","3454","3455","342","13")
 x2=c("1","1","1","0","0")
 data1=data.frame(x1,x2)

 y1=c("232","232","3454","3454","3455","342","13","13","13","13")
 y2=c("E13","B04","F19","A64","E22","H44","F68","G84","F19","A01")
 data2=data.frame(y1,y2)

 I want now to creat a loop, which creates for every value in list1 a
new
 binary variable in data1. Result should look like:
 x1  x2  A01 B04 A64 G84 F19
 232 1   0   1   0   0   0
 34541   0   0   1   0   1
 34551   0   0   0   0   0
 342 0   0   0   0   0   0
 13  0   1   0   0   1   1
>>>
>>> Loops!?! We don't nee no steenking loops!
>>>
 xtb <-  with(data2, table(y1,y2))
 cbind(data1, xtb[match(data1$x1, rownames(xtb)), ] )
>>>
>>>  x1 x2 A01 A64 B04 E13 E22 F19 F68 G84 H44
>>> 232   232  1   0   0   1   1   0   0   0   0   0
>>> 3454 3454  1   0   1   0   0   0   1   0   0   0
>>> 3455 3455  1   0   0   0   0   1   0   0   0   0
>>> 342   342  0   0   0   0   0   0   0   0   0   1
>>> 13 13  0   1   0   0   0   0   1   1   1   0
>>>
>>> I am guessing that you were to ... er, busy? ... to complete the table?
>>>
>>> --
>>>
>>> David Winsemius, MD
>>> West Hartford, CT
>>>
>>>
>>
>> Thanks a lot! Pretty simple. I am so much used to SQLDF right now.
>>
>> So how would you handle more complicated strings like that:
>> y1=c("232","232", "232", "3454","3454","3455","342","13","13","13","13")
>> y2=c("E13","B04 A01 F19","B04","F19","A64 G84 A05","E22","H44
>> C35","F68","G84","F19","A01")
>> data2=data.frame(y1,y2)
>>
>> Where you want to extract for instance all "A01" from the strings?
>
> I think you need either to explain what you want in more words of the
> English language or to offer an example of the desired output. I suspect
you
> did not want something as simple as this:
>
>> A01.instances <- grep("A01" , data2$y2)
>> A01.instances
> [1]  2 11
>> data2[A01.instances, ]
>y1  y2
> 2  232 B04 A01 F19
> 11  13 A01
>
> Or maybe you did?
>
> --
> David Winsemius, MD
> West Hartford, CT
>
>

No, that was not my intention. Consider the following example:

list1=c("A01","B04","A64","G84","F19") # My "substrings" to screen for in
> data2
>
>
> x1=c("232","3454","3455","342","13")
> x2=c("1","1","1","0","0")
> data1=data.frame(x1,x2) # Target dataframe where the 5 new binary variables
> (namely from list1) are added
>
>
> y1=c("232","232", "232", "3454","3454","3455","342","13","13","13","13")
> y2=c("E133","B04 A01A F194","B04","F19","A642 G84 A05","E223","H44
> C35","F68","G84","F19","A01")
> data2=data.frame(y1,y2) # Dataframe to be screen by list1
>

Result should look like this:

x1  x2  A01 B04 A64 G84 F19
> 232 1   1   1   0   0   0
> 34541   0   0   1   0   1
> 34551   0   0   0   0   0
> 342 0   0   0   0   0   0
> 13  0   1   0   0   1   1
>

[[alternative HTML version deleted]]

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

Re: [R] Cumulative dose in survival models

2011-05-10 Thread peter dalgaard

On May 10, 2011, at 10:52 , kl...@kaae-consult.dk wrote:

> 
> Dear List,
>  
> I want to fit a CPH survival model with a time-dependent cumulative exposure 
> variable. More specifically, I have exposure measurements for all cohort 
> members at discrete times (yearly measurements) during follow-up.
>  
> I suppose that when estimating the model parameters the most correct would be 
> to have the cumulative exposure for everybody at risk at each death time, as 
> I cannot assume the cumulative exposure to be constant between the yearly 
> measurements . Thus, I need to calculate (interpolate) cumulative exposure 
> for everybody at risk at each death time and pass this to the estimation 
> routine.
>  
> Does anyone know how to do this in R? Any help would be much appreciated.

You didn't say what your data layout is. One record per yearly measurement?

Generic time splitting is discussed in the 2nd ed. of ISwR (or peek into the 
scripts in the ISwR package). That is for a simplified situation (5 year age 
groups) but should generalize to arbitrary cutpoints with some elaboration. 
Then you "just" need to figure out the exposure. If you can things narrowed 
down to 

- cumulative exposure at last measurement
- current exposure
- time since last measurement at end of observation interval

then it would be obvious how to proceed. One worry would be whether you have 
the current exposure for those who die. I.e., you do not want to use 
interpolations that depend on _next_ year's measurements. 


-- 
Peter Dalgaard
Center for Statistics, Copenhagen Business School
Solbjerg Plads 3, 2000 Frederiksberg, Denmark
Phone: (+45)38153501
Email: pd@cbs.dk  Priv: pda...@gmail.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] Reading a large netCDF file using R

2011-05-10 Thread David William Pierce
On Mon, May 9, 2011 at 6:29 PM, sulochandhungel
wrote:

> I got a netCDF file from a climate group and wanted to aggregate the data.
> The file summary looks like this
> [...]
>
> I didnt think it worked... can u suggest me some other way ... I could not
> use the start and count properly


 Could you be a little more explicit? Why don't you think it worked?

In general, you seem to be using count and start correctly.

--Dave


-- 
David W. Pierce
Division of Climate, Atmospheric Science, and Physical Oceanography
Scripps Institution of Oceanography
(858) 534-8276 (voice)  /  (858) 534-8561 (fax)dpie...@ucsd.edu

[[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] Creating binary variable depending on strings of two dataframes

2011-05-10 Thread noxyp...@gmail.com
On Tue, May 10, 2011 at 3:09 PM, David Winsemius wrote:

>
> On May 10, 2011, at 3:18 AM, noxyp...@gmail.com wrote:
>
>  On Fri, May 6, 2011 at 7:41 PM, David Winsemius 
>> wrote:
>>
>>>
>>> On May 6, 2011, at 11:35 AM, Pete Pete wrote:
>>>
>>>
 Gabor Grothendieck wrote:

>
> On Tue, Dec 7, 2010 at 11:30 AM, Pete Pete 
> wrote:
>
>>
>> Hi,
>> consider the following two dataframes:
>> x1=c("232","3454","3455","342","13")
>> x2=c("1","1","1","0","0")
>> data1=data.frame(x1,x2)
>>
>> y1=c("232","232","3454","3454","3455","342","13","13","13","13")
>> y2=c("E1","F3","F5","E1","E2","H4","F8","G3","E1","H2")
>> data2=data.frame(y1,y2)
>>
>> I need a new column in dataframe data1 (x3), which is either 0 or 1
>> depending if the value "E1" in y2 of data2 is true while x1=y1. The
>> result
>> of data1 should look like this:
>>  x1 x2 x3
>> 1 232   1   1
>> 2 3454 1   1
>> 3 3455 1   0
>> 4 342   0   0
>> 5 13 0   1
>>
>> I think a SQL command could help me but I am too inexperienced with it
>> to
>> get there.
>>
>>
> Try this:
>
>  library(sqldf)
>> sqldf("select x1, x2, max(y2 = 'E1') x3 from data1 d1 left join data2
>> d2
>> on (x1 = y1) group by x1, x2 order by d1.rowid")
>>
>
>  x1 x2 x3
> 1  232  1  1
> 2 3454  1  1
> 3 3455  1  0
> 4  342  0  0
> 5   13  0  1
>
>
>  snipped Gabor's sig
>>>

 That works pretty cool but I need to automate this a bit more. Consider
 the
 following example:

 list1=c("A01","B04","A64","G84","F19")

 x1=c("232","3454","3455","342","13")
 x2=c("1","1","1","0","0")
 data1=data.frame(x1,x2)

 y1=c("232","232","3454","3454","3455","342","13","13","13","13")
 y2=c("E13","B04","F19","A64","E22","H44","F68","G84","F19","A01")
 data2=data.frame(y1,y2)

 I want now to creat a loop, which creates for every value in list1 a new
 binary variable in data1. Result should look like:
 x1  x2  A01 B04 A64 G84 F19
 232 1   0   1   0   0   0
 34541   0   0   1   0   1
 34551   0   0   0   0   0
 342 0   0   0   0   0   0
 13  0   1   0   0   1   1

>>>
>>> Loops!?! We don't nee no steenking loops!
>>>
>>>  xtb <-  with(data2, table(y1,y2))
 cbind(data1, xtb[match(data1$x1, rownames(xtb)), ] )

>>>  x1 x2 A01 A64 B04 E13 E22 F19 F68 G84 H44
>>> 232   232  1   0   0   1   1   0   0   0   0   0
>>> 3454 3454  1   0   1   0   0   0   1   0   0   0
>>> 3455 3455  1   0   0   0   0   1   0   0   0   0
>>> 342   342  0   0   0   0   0   0   0   0   0   1
>>> 13 13  0   1   0   0   0   0   1   1   1   0
>>>
>>> I am guessing that you were to ... er, busy? ... to complete the table?
>>>
>>> --
>>>
>>> David Winsemius, MD
>>> West Hartford, CT
>>>
>>>
>>>
>> Thanks a lot! Pretty simple. I am so much used to SQLDF right now.
>>
>> So how would you handle more complicated strings like that:
>> y1=c("232","232", "232", "3454","3454","3455","342","13","13","13","13")
>> y2=c("E13","B04 A01 F19","B04","F19","A64 G84 A05","E22","H44
>> C35","F68","G84","F19","A01")
>> data2=data.frame(y1,y2)
>>
>> Where you want to extract for instance all "A01" from the strings?
>>
>
> I think you need either to explain what you want in more words of the
> English language or to offer an example of the desired output. I suspect you
> did not want something as simple as this:
>
> > A01.instances <- grep("A01" , data2$y2)
> > A01.instances
> [1]  2 11
> > data2[A01.instances, ]
>y1  y2
> 2  232 B04 A01 F19
> 11  13 A01
>
> Or maybe you did?
>
> --
> David Winsemius, MD
> West Hartford, CT
>
>
With sqldf I could do it manually:

> > data1=sqldf("SELECT data1.*, max(data2.y2 LIKE '% A01%') OR max(data2.y2
> LIKE 'A01%') A01 FROM data1 left join data2 on (data1.x1 = data2.y1) group
> by data1.x1, data2.y1")
> > data1=sqldf("SELECT data1.*, max(data2.y2 LIKE '% B04%') OR max(data2.y2
> LIKE 'B04%') B04 FROM data1 left join data2 on (data1.x1 = data2.y1) group
> by data1.x1, data2.y1")
> > data1
> x1 x2 A01 B04
> 1   13  0   1   0
> 2  232  1   1   1
> 3  342  0   0   0
> 4 3454  1   0   0
> 5 3455  1   0   0
> >
>

But I need to automate this for some thousand "substrings".
Any suggestion?

[[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] Creating binary variable depending on strings of two dataframes

2011-05-10 Thread David Winsemius


On May 10, 2011, at 9:49 AM, noxyp...@gmail.com wrote:

On Tue, May 10, 2011 at 3:09 PM, David Winsemius >

wrote:


On May 10, 2011, at 3:18 AM, noxyp...@gmail.com wrote:

On Fri, May 6, 2011 at 7:41 PM, David Winsemius >

wrote:


On May 6, 2011, at 11:35 AM, Pete Pete wrote:



Gabor Grothendieck wrote:


On Tue, Dec 7, 2010 at 11:30 AM, Pete Pete  


wrote:


Hi,
consider the following two dataframes:
x1=c("232","3454","3455","342","13")
x2=c("1","1","1","0","0")
data1=data.frame(x1,x2)

y1=c("232","232","3454","3454","3455","342","13","13","13","13")
y2=c("E1","F3","F5","E1","E2","H4","F8","G3","E1","H2")
data2=data.frame(y1,y2)

I need a new column in dataframe data1 (x3), which is either 0  
or 1
depending if the value "E1" in y2 of data2 is true while  
x1=y1. The

result
of data1 should look like this:
x1 x2 x3
1 232   1   1
2 3454 1   1
3 3455 1   0
4 342   0   0
5 13 0   1

I think a SQL command could help me but I am too inexperienced  
with

it

to
get there.



Try this:


library(sqldf)
sqldf("select x1, x2, max(y2 = 'E1') x3 from data1 d1 left  
join data2

d2
on (x1 = y1) group by x1, x2 order by d1.rowid")


x1 x2 x3
1  232  1  1
2 3454  1  1
3 3455  1  0
4  342  0  0
5   13  0  1



snipped Gabor's sig


That works pretty cool but I need to automate this a bit more.  
Consider

the
following example:

list1=c("A01","B04","A64","G84","F19")

x1=c("232","3454","3455","342","13")
x2=c("1","1","1","0","0")
data1=data.frame(x1,x2)

y1=c("232","232","3454","3454","3455","342","13","13","13","13")
y2=c("E13","B04","F19","A64","E22","H44","F68","G84","F19","A01")
data2=data.frame(y1,y2)

I want now to creat a loop, which creates for every value in  
list1 a

new

binary variable in data1. Result should look like:
x1  x2  A01 B04 A64 G84 F19
232 1   0   1   0   0   0
34541   0   0   1   0   1
34551   0   0   0   0   0
342 0   0   0   0   0   0
13  0   1   0   0   1   1


Loops!?! We don't nee no steenking loops!


xtb <-  with(data2, table(y1,y2))
cbind(data1, xtb[match(data1$x1, rownames(xtb)), ] )


   x1 x2 A01 A64 B04 E13 E22 F19 F68 G84 H44
232   232  1   0   0   1   1   0   0   0   0   0
3454 3454  1   0   1   0   0   0   1   0   0   0
3455 3455  1   0   0   0   0   1   0   0   0   0
342   342  0   0   0   0   0   0   0   0   0   1
13 13  0   1   0   0   0   0   1   1   1   0

I am guessing that you were to ... er, busy? ... to complete the  
table?


--

David Winsemius, MD
West Hartford, CT




Thanks a lot! Pretty simple. I am so much used to SQLDF right now.

So how would you handle more complicated strings like that:
y1=c("232","232", "232",  
"3454","3454","3455","342","13","13","13","13")

y2=c("E13","B04 A01 F19","B04","F19","A64 G84 A05","E22","H44
C35","F68","G84","F19","A01")
data2=data.frame(y1,y2)

Where you want to extract for instance all "A01" from the strings?


I think you need either to explain what you want in more words of the
English language or to offer an example of the desired output. I  
suspect

you

did not want something as simple as this:


A01.instances <- grep("A01" , data2$y2)
A01.instances

[1]  2 11

data2[A01.instances, ]

 y1  y2
2  232 B04 A01 F19
11  13 A01

Or maybe you did?

--
David Winsemius, MD
West Hartford, CT




No, that was not my intention. Consider the following example:

list1=c("A01","B04","A64","G84","F19") # My "substrings" to screen  
for in

data2


x1=c("232","3454","3455","342","13")
x2=c("1","1","1","0","0")
data1=data.frame(x1,x2) # Target dataframe where the 5 new binary  
variables

(namely from list1) are added


y1=c("232","232", "232",  
"3454","3454","3455","342","13","13","13","13")

y2=c("E133","B04 A01A F194","B04","F19","A642 G84 A05","E223","H44
C35","F68","G84","F19","A01")
data2=data.frame(y1,y2) # Dataframe to be screen by list1



Result should look like this:

x1  x2  A01 B04 A64 G84 F19

232 1   1   1   0   0   0
34541   0   0   1   0   1
34551   0   0   0   0   0
342 0   0   0   0   0   0
13  0   1   0   0   1   1


And how were we supposed to figure out that 3454/G84 was not supposed  
to be counted?


OK. let's assume you were just sloppy ... then build a new data.frame:

> data4 <- data.frame(y1=rep(data3[,1],
sapply (strsplit(gsub("\\\n"," ",data3$y2), " "),  
length) ),
 y2 = unlist (strsplit(gsub("\\\n"," ",data3$y2),  
" ") ) + )

> data4
y1  y2
1   232 E13
2   232 B04
3   232 A01
4   232 F19
5   232 B04
6  3454 F19
7  3454 A64
8  3454 G84
9  3454 A05
10 3455 E22
11  342 H44
12  342 C35
13   13 F68
14   13 G84
15   13 F19
16   13 A01



--
David Winsemius, MD
West Hartford, CT

__
R-help@

Re: [R] Saving multiple 3x3 TIFF graphics inside a loop

2011-05-10 Thread nidal . alhuniti
Try this

dat=data.frame (ID=rep(1:27,each=10),IDV=rep(seq(1:10),times=27))
dat$DV <- with(dat, 50+15*IDV)
dat=dat[order(dat$ID,dat$IDV),]
for(i in 1:27){
   dt1 = dat[dat$ID==i,]
   pagei=as.integer(ceiling(i/9))

   #open device for subjects 1-9
   if(i==1) {tiff(file=paste("Outfile",pagei,".tiff",sep=""),width=17.15
 ,height=17.15,units="cm", res=1200, pointsize=10, compression =
"lzw")}
   if(i==1) {par(mar=c(5,5,4,1),mfrow=c(3,3))}

   ## open new file for subjects 10-18
   if(i==10) {tiff(file=paste("Outfile",pagei,".tiff",sep=""),width=17.15
   ,height=17.15,units="cm", res=1200, pointsize=10, compression = "lzw")}
   if(i==10) {par(mar=c(5,5,4,1),mfrow=c(3,3))}

   ## open new file for subjects 19-27
   if(i==19) {tiff(file=paste("Outfile",pagei,".tiff",sep=""),width=17.15
   ,height=17.15,units="cm", res=1200, pointsize=10, compression = "lzw")}
   if(i==19) {par(mar=c(5,5,4,1),mfrow=c(3,3))}

   plot(dt1[, "IDV"],dt1[, "DV"], type="n",ylim=range(dt1[, "DV"]),
xlim=range(dt1[, "IDV"]))
   lines(dt1[, "IDV"], dt1[, "DV"], type="b", col=1,lwd=3,pch=1)
   title(paste(i),col=8,cex=.8)

   if(i==9) dev.off ()  ## close file no 1
   if(i==18) dev.off () ## close file no 2
   if(i==27) dev.off () ## close file no 3

}




On Tue, May 10, 2011 at 9:25 AM, John S  wrote:

> Thanks Uwe . The parentheses were a typo. Please note that if I have only 9
> subjects in my data then the following works (creates 493 KB data file)
>
>
>
> #windows()
>
> {tiff(file=paste("Outfile",1,".tiff",sep=""),width=17.15
>
>  ,height=17.15,units="cm", res=1200, pointsize=10, compression =
> "lzw")}
>
> par(mar=c(5,5,4,1),mfrow=c(3,3))
>
> for(i in 1:9){
>
>   dt1 = dat[dat$ID==i,]
>
>   #pagei=as.integer(ceiling(i/9))
>
>   #This does not work
>
>   #if(i==1) {tiff(file=paste("Outfile",pagei,".tiff",sep=""),width=17.15
>
>   ##,height=17.15,units="cm", res=1200, pointsize=10, compression = "lzw")}
>
>
>   #par(mar=c(5,5,4,1),mfrow=c(3,3))
>
>   plot(dt1[, "IDV"],dt1[, "DV"], type="n",ylim=range(dt1[, "DV"]),
> xlim=range(dt1[, "IDV"]))
>
>   lines(dt1[, "IDV"], dt1[, "DV"], type="b", col=1,lwd=3,pch=1)
>
>   title(paste(i),col=8,cex=.8)
>
>   #if(i==5) dev.off ## close file
>
>   ## open new file
>
>   #if(i==5) {tiff(file=paste("Outfile",pagei,".tiff",sep=""),width=17.15
>
>   ##,height=17.15,units="cm", res=1200, pointsize=10, compression = "lzw")}
>
>   #if(i==9) dev.off () ## close file
>
>   #if(i==9) {tiff(file=paste("Outfile",pagei,".tiff",sep=""),width=17.15
>
>   ##,height=17.15,units="cm", res=1200, pointsize=10, compression = "lzw")}
>
> }
>
> dev.off()
>
>
>
>
>
> My workaround around was to break the dataset into small data each with 9
> subjects but if I change the mfrow then I have to change no. of subjects.
> Any smart suggestion to get around this? I tried to open and close devices
> based on iterations (i) but it does not work for me
>
>
>
>
> 2011/5/10 Uwe Ligges 
>
> >
> >
> > On 10.05.2011 03:36, John S wrote:
> >
> >> Dear Friends,
> >>
> >> I have been trying to save multiple 3x3 (mfrow=c(3,3) graphics inside a
> >> loop
> >> using tiff figure format (not using PDF or savePlot functions) with no
> >> success. Could you please help?
> >>
> >> Here is a simplified example code:
> >>
> >>
> >>
> >> dat=data.frame (ID=rep(1:10,each=10),IDV=rep(seq(1:10),times=10))
> >>
> >> dat$DV<- with(dat, 50+15*IDV)
> >>
> >> dat=dat[order(dat$ID,dat$IDV),]
> >>
> >> for(i in 1:10){
> >>
> >>dt1 = dat[dat$ID==i,]
> >>
> >>pagei=as.integer(ceiling(i/9))
> >>
> >>#This does not work
> >>
> >>if(i==1) {tiff(file=paste("Outfile",pagei,".tiff",sep="")
> >>
> >>  ,width=17.15,height=17.15,units="cm", res=1200, pointsize=10,
> >> compression = "lzw")}
> >>
> >
> > You may not be able to allocate space for this huge tiffs (at least I am
> > not able to do so under Windows).
> >
> >
> >
> >
> >   if(i==10) {tiff(file=paste("Outfile",pagei,".tiff",sep="")
> >>
> >>  ,width=17.15,height=17.15,units="cm", res=1200, pointsize=10,
> >> compression = "lzw")}
> >>
> >>par(mar=c(5,5,4,1),mfrow=c(3,3))
> >>
> >>plot(dt1[, "IDV"],dt1[, "DV"], type="n")
> >>
> >>lines(dt1[, "IDV"], dt1[, "DV"], type="b", col=1,lwd=3,pch=1)
> >>
> >>title(paste(i),col=8,cex=.8)
> >>
> >>  if(i==9) dev.off ## close file
> >>
> >
> > Your forgot parentheses after dev.off.
> >
> > Uwe Ligges
> >
> >
> >  }
> >>
> >> dev.off()
> >>
> >>  The idea is to save 2 files with 3x3 layouts but I only get one subject
> >> in
> >> Outfile1.
> >>
> >>
> >>
> >> OS: Windows 7 64 bits
> >>
> >> R version 2.13.0
> >>
> >> Many Thanks,
> >>
> >> John
> >>
> >>[[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
> 

Re: [R] Weird read.xls behavior

2011-05-10 Thread Jun Shen
Hi, Kenneth,

I would like a bit more advice on how to use gsub to remove the spaces at
the beginning or the end for a data frame and still keep the data frame in
the same format. The way I see is something like

gsub(' *$', '' , unlist(dataframe)). But the result is not a data frame
anymore. Thanks.

Jun

On Mon, May 9, 2011 at 11:01 PM, Kenneth Roy Cabrera Torres <
krcab...@une.net.co> wrote:

> Maybe in the original there are some "ghost" spaces after ng/ml,
> verify them.
>
> Any way you can erase the first and the last white spaces with gsub
> function.
>
> HTH
>
> Kenneth
> El lun, 09-05-2011 a las 22:44 -0500, Jun Shen escribió:
> > Dear list,
> >
> > I used read.xls to read in an excel file, which has both character and
> > numeric columns. Everything seems fine except in the last column, I have
> > this character value "ng/ml". When reading in, read.xls seems to add a
> space
> > at the end of it, became "ng/ml ". How can I prevent read.xls doing so?
> > Thanks.
> >
> > >read.xls(data, header=T, as.is=T)
> >
> > Jun
> >
> >   [[alternative HTML version deleted]]
> >
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
>
>
>

[[alternative HTML version deleted]]

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


Re: [R] Weird read.xls behavior

2011-05-10 Thread David Winsemius


On May 10, 2011, at 10:54 AM, Jun Shen wrote:


Hi, Kenneth,

I would like a bit more advice on how to use gsub to remove the  
spaces at
the beginning or the end for a data frame and still keep the data  
frame in

the same format. The way I see is something like

gsub(' *$', '' , unlist(dataframe)). But the result is not a data  
frame

anymore. Thanks.



1) Please do not attach unrelated questions to existing threads.
2) Please do include reproducible code that creates an example to work  
with.
3) Please do read the help page for gsub with careful attention to the  
datatypes it expects as input and the values it will return.


--
David.

Jun

On Mon, May 9, 2011 at 11:01 PM, Kenneth Roy Cabrera Torres <
krcab...@une.net.co> wrote:


Maybe in the original there are some "ghost" spaces after ng/ml,
verify them.

Any way you can erase the first and the last white spaces with gsub
function.

HTH

Kenneth
El lun, 09-05-2011 a las 22:44 -0500, Jun Shen escribió:

Dear list,

I used read.xls to read in an excel file, which has both character  
and
numeric columns. Everything seems fine except in the last column,  
I have
this character value "ng/ml". When reading in, read.xls seems to  
add a

space
at the end of it, became "ng/ml ". How can I prevent read.xls  
doing so?

Thanks.


read.xls(data, header=T, as.is=T)


Jun

 [[alternative HTML version deleted]]

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

http://www.R-project.org/posting-guide.html

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






[[alternative HTML version deleted]]

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


David Winsemius, MD
West Hartford, CT

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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 symbol in evaluate expression function in axis labels

2011-05-10 Thread fernando espindola
Hi R-helps

I have the follow problem when try to put degree symbol in axis labels via
expression function in the next example code, but when put infinity in
expression put the degree symbol

x<-sample(10);y<-sample(10)
plot(x,y,pch="+",axes=FALSE,ylab="",xlab="")
ax<-seq(1,10,by=2)
axis(1,at=ax,labels=parse(text=paste(ax,"*degree~W",sep="")))
axis(2,at=ax,labels=parse(text=paste(ax,"*infinity~S,sep="")),las=2)

Thank for all

This is my R version...

R> R.Version()
$platform
[1] "i686-pc-linux-gnu"

$arch
[1] "i686"

$os
[1] "linux-gnu"

$system
[1] "i686, linux-gnu"

$status
[1] ""

$major
[1] "2"

$minor
[1] "13.0"

$year
[1] "2011"

$month
[1] "04"

$day
[1] "13"

$`svn rev`
[1] "55427"

$language
[1] "R"

$version.string
[1] "R version 2.13.0 (2011-04-13)"

[[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] need to delete by time, not date

2011-05-10 Thread Schatzi
I have a matrix where one column has a date and another column has a time. I
need to delete all times before 6am. I had combined the Date and Time column
into DateTime.

Mat1:
Weight Date  Time
7.6 04/28/11 09:03 
8.4 04/29/11 03:11 
8.6 04/29/11 05:32 
8.6 04/29/11 09:53 
1.4 05/01/11 19:52 

I tried just picking up the time:

as.POSIXct(Mat1$Time, format = "%H:%M")

but this assigns the time value today's date.

The final matrix should look like:
Weight Date  Time   
7.6 04/28/11 09:03 
8.6 04/29/11 09:53 
1.4 05/01/11 19:52 


-
In theory, practice and theory are the same. In practice, they are not - Albert 
Einstein
--
View this message in context: 
http://r.789695.n4.nabble.com/need-to-delete-by-time-not-date-tp3512162p3512162.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] need to delete by time, not date

2011-05-10 Thread Clint Bowman

You may be happy using:

strptime(paste(Date,Time,sep=" "),"%m/%d/%y %H:%M")$hour

--
Clint BowmanINTERNET:   cl...@ecy.wa.gov
Air Quality Modeler INTERNET:   cl...@math.utah.edu
Department of Ecology   VOICE:  (360) 407-6815
PO Box 47600FAX:(360) 407-7534
Olympia, WA 98504-7600


USPS:   PO Box 47600, Olympia, WA 98504-7600
Parcels:300 Desmond Drive, Lacey, WA 98503-1274


On Tue, 10 May 2011, Schatzi wrote:


I have a matrix where one column has a date and another column has a time. I
need to delete all times before 6am. I had combined the Date and Time column
into DateTime.

Mat1:
Weight Date  Time
7.6 04/28/11 09:03
8.4 04/29/11 03:11
8.6 04/29/11 05:32
8.6 04/29/11 09:53
1.4 05/01/11 19:52

I tried just picking up the time:

as.POSIXct(Mat1$Time, format = "%H:%M")

but this assigns the time value today's date.

The final matrix should look like:
Weight Date  Time
7.6 04/28/11 09:03
8.6 04/29/11 09:53
1.4 05/01/11 19:52


-
In theory, practice and theory are the same. In practice, they are not - Albert 
Einstein
--
View this message in context: 
http://r.789695.n4.nabble.com/need-to-delete-by-time-not-date-tp3512162p3512162.html
Sent from the R help mailing list archive at Nabble.com.

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



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


Re: [R] need to delete by time, not date

2011-05-10 Thread Adele_Thompson
I am already using that to combine the date and time (I had it in the code, but 
took it out for the post to simplify it), but I still have the problem where I 
want to remove any time before 6am across all dates.

-Original Message-
From: cl...@ecy.wa.gov [mailto:cl...@ecy.wa.gov] 
Sent: Tuesday, May 10, 2011 10:15 AM
To: Thompson, Adele - adele_thomp...@cargill.com
Cc: r-help@r-project.org
Subject: Re: [R] need to delete by time, not date

You may be happy using:

strptime(paste(Date,Time,sep=" "),"%m/%d/%y %H:%M")$hour

-- 
Clint BowmanINTERNET:   cl...@ecy.wa.gov
Air Quality Modeler INTERNET:   cl...@math.utah.edu
Department of Ecology   VOICE:  (360) 407-6815
PO Box 47600FAX:(360) 407-7534
Olympia, WA 98504-7600


 USPS:   PO Box 47600, Olympia, WA 98504-7600
 Parcels:300 Desmond Drive, Lacey, WA 98503-1274


On Tue, 10 May 2011, Schatzi wrote:

> I have a matrix where one column has a date and another column has a time. I
> need to delete all times before 6am. I had combined the Date and Time column
> into DateTime.
>
> Mat1:
> Weight Date  Time
> 7.6 04/28/11 09:03
> 8.4 04/29/11 03:11
> 8.6 04/29/11 05:32
> 8.6 04/29/11 09:53
> 1.4 05/01/11 19:52
>
> I tried just picking up the time:
>
> as.POSIXct(Mat1$Time, format = "%H:%M")
>
> but this assigns the time value today's date.
>
> The final matrix should look like:
> Weight Date  Time
> 7.6 04/28/11 09:03
> 8.6 04/29/11 09:53
> 1.4 05/01/11 19:52
>
>
> -
> In theory, practice and theory are the same. In practice, they are not - 
> Albert Einstein
> --
> View this message in context: 
> http://r.789695.n4.nabble.com/need-to-delete-by-time-not-date-tp3512162p3512162.html
> Sent from the R help mailing list archive at Nabble.com.
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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


Re: [R] Weird read.xls behavior

2011-05-10 Thread Gabor Grothendieck
On Tue, May 10, 2011 at 10:59 AM, David Winsemius
 wrote:
>
> On May 10, 2011, at 10:54 AM, Jun Shen wrote:
>
>> Hi, Kenneth,
>>
>> I would like a bit more advice on how to use gsub to remove the spaces at
>> the beginning or the end for a data frame and still keep the data frame in
>> the same format. The way I see is something like
>>
>> gsub(' *$', '' , unlist(dataframe)). But the result is not a data frame
>> anymore. Thanks.
>>
>
> 1) Please do not attach unrelated questions to existing threads.
> 2) Please do include reproducible code that creates an example to work with.
> 3) Please do read the help page for gsub with careful attention to the

To be fair this was related and he did send me something reproducible
privately although you would not know that.

-- 
Statistics & Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.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] need to delete by time, not date

2011-05-10 Thread David Winsemius


On May 10, 2011, at 11:06 AM, Schatzi wrote:

I have a matrix where one column has a date and another column has a  
time. I
need to delete all times before 6am. I had combined the Date and  
Time column

into DateTime.

Mat1:
Weight Date  Time
7.6 04/28/11 09:03
8.4 04/29/11 03:11
8.6 04/29/11 05:32
8.6 04/29/11 09:53
1.4 05/01/11 19:52

I tried just picking up the time:

as.POSIXct(Mat1$Time, format = "%H:%M")

but this assigns the time value today's date.


Since these are either character strings (or factors) at the  
moment ... and if they are matrices they are character...


Mat1[Mat1$Time >= "06:00" , ]

If actually a dataframe then

Mat1[as.character(Mat1$Time) >= "06:00" , ]

In the future please include the output of dput(Mat1).



The final matrix should look like:
Weight Date  Time
7.6 04/28/11 09:03
8.6 04/29/11 09:53
1.4 05/01/11 19:52




David Winsemius, MD
West Hartford, CT

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] need to delete by time, not date

2011-05-10 Thread Adele_Thompson
Thanks, that did exactly what I wanted. I realize that I didn't put it in the 
best form, but was having trouble coming up with a replicable example. I see 
that many people here quickly create examples and I am learning how to do that, 
but still struggle, especially when using as.POSIXct and other more complex 
classes.
Thanks for the help.

-Original Message-
From: dwinsem...@comcast.net [mailto:dwinsem...@comcast.net] 
Sent: Tuesday, May 10, 2011 10:19 AM
To: Thompson, Adele - adele_thomp...@cargill.com
Cc: r-help@r-project.org
Subject: Re: [R] need to delete by time, not date


On May 10, 2011, at 11:06 AM, Schatzi wrote:

> I have a matrix where one column has a date and another column has a  
> time. I
> need to delete all times before 6am. I had combined the Date and  
> Time column
> into DateTime.
>
> Mat1:
> Weight Date  Time
> 7.6 04/28/11 09:03
> 8.4 04/29/11 03:11
> 8.6 04/29/11 05:32
> 8.6 04/29/11 09:53
> 1.4 05/01/11 19:52
>
> I tried just picking up the time:
>
> as.POSIXct(Mat1$Time, format = "%H:%M")
>
> but this assigns the time value today's date.

Since these are either character strings (or factors) at the  
moment ... and if they are matrices they are character...

Mat1[Mat1$Time >= "06:00" , ]

If actually a dataframe then

Mat1[as.character(Mat1$Time) >= "06:00" , ]

In the future please include the output of dput(Mat1).

>
> The final matrix should look like:
> Weight Date  Time
> 7.6 04/28/11 09:03
> 8.6 04/29/11 09:53
> 1.4 05/01/11 19:52
>
>

David Winsemius, MD
West Hartford, CT

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

2011-05-10 Thread Pavan G
Thanks David and Steve. This is exactly what I was looking for. tapply seems
to be one robust command.
(Humor of any potency is always a welcome thing in academia :)
WSS

[[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] Weird read.xls behavior

2011-05-10 Thread Gabor Grothendieck
On Tue, May 10, 2011 at 12:12 AM, Jun Shen  wrote:
> Kenneth,
>
> Thanks for the reply. I checked the original data. There is no space. I even
> manually added a space to one value. After reading in with read.xls, the
> value has two spaces. The reason I don't like it is I am going to do some
> comparison with another dataset, which is supposed to be the same as this
> one. Now I am getting a bunch of false negatives.

It seems that the perl program underlying gdata's read.xls puts out
lines like this:

|"KAI-4169-002","830","5 mg" |
where | mark the beginning and end and are not part of the line.
read.csv includes the space after the last double quote in the last
field even though its outside of the double quote.

As an interim fix, edit the file at this location:

   system.file("perl", "xls2csv.pl", package = "gdata")

removing the space before the \n in this line:
   print OutFile "$outputLine \n"
so it becomes this:
   print OutFile "$outputLine\n"

Now it should work.

-- 
Statistics & Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.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] maximize portfolio return under risk target

2011-05-10 Thread xavier abulker
Hello,

I'm trying to build a simple example to maximize the return on a portfolio with 
a risk target (12%) and long  positions only (0 <= weights <=1).
I've read all the previous posts on fPortfolio and PortfolioAnalytics but I 
cannot find anything simple.
Could you please help me? My code currently looks like this:


library(fPortfolio)
library(PortfolioAnalytics)
TargetRisk<-0.12 #my target risk = sd
## return data 
Data <- SMALLCAP.RET
Data <- Data[, c("BKE", "GG", "GYMB", "KRON")]



Regards
[[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] fitting non-intercept model with lrm

2011-05-10 Thread Clarissa Valim
I would appreciate if someone could tell me how to fit a non-intercept model
using lrm (and not glm). The -1 in the formula of the glm does not work with
lrm.

Thanks,

Clarissa

[[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] need to delete by time, not date

2011-05-10 Thread David Winsemius


On May 10, 2011, at 11:23 AM, adele_thomp...@cargill.com wrote:

Thanks, that did exactly what I wanted. I realize that I didn't put  
it in the best form, but was having trouble coming up with a  
replicable example. I see that many people here quickly create  
examples and I am learning how to do that, but still struggle,  
especially when using as.POSIXct and other more complex classes.


Excuse not accepted. All you needed to do with Mat1 was dput(Mat1). Or  
if your real data object was a lot longer, then just dput(head(Mat1,  
10)).


--
David.

Thanks for the help.

-Original Message-
From: dwinsem...@comcast.net [mailto:dwinsem...@comcast.net]
Sent: Tuesday, May 10, 2011 10:19 AM
To: Thompson, Adele - adele_thomp...@cargill.com
Cc: r-help@r-project.org
Subject: Re: [R] need to delete by time, not date


On May 10, 2011, at 11:06 AM, Schatzi wrote:


I have a matrix where one column has a date and another column has a
time. I
need to delete all times before 6am. I had combined the Date and
Time column
into DateTime.

Mat1:
Weight Date  Time
7.6 04/28/11 09:03
8.4 04/29/11 03:11
8.6 04/29/11 05:32
8.6 04/29/11 09:53
1.4 05/01/11 19:52

I tried just picking up the time:

as.POSIXct(Mat1$Time, format = "%H:%M")

but this assigns the time value today's date.


Since these are either character strings (or factors) at the
moment ... and if they are matrices they are character...

Mat1[Mat1$Time >= "06:00" , ]

If actually a dataframe then

Mat1[as.character(Mat1$Time) >= "06:00" , ]

In the future please include the output of dput(Mat1).



The final matrix should look like:
Weight Date  Time
7.6 04/28/11 09:03
8.6 04/29/11 09:53
1.4 05/01/11 19:52




David Winsemius, MD
West Hartford, CT



David Winsemius, MD
West Hartford, CT

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


Re: [R] Multiple comparisons for coxph survival analysis model

2011-05-10 Thread jul 2-pom
Hi ExpeRts,

I find this post very interesting... as I have the same problem. I found on
other webpages that a solution would be to run simple comparison and then to
adjust for the P-value accordingly; However, I am not very satisfied by such
approach, and i would very greatful to R if a multiple comparison would
exist.

Thank You
Julien Pompon, a biologist that difficulty try to use R

--
View this message in context: 
http://r.789695.n4.nabble.com/Multiple-comparisons-for-coxph-survival-analysis-model-tp906080p3512079.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Using $ accessor in GAM formula

2011-05-10 Thread Gene Leynes
Gavin,

Yes, I thought I was responding to Berwin, whoops!

Since I'm often mixing different data frames, I rarely use this formulation:
gam(y ~ s(x), dat)

Actually, I often forget about that syntax, and it's a good reminder.  It
really is much easier to read.  In my "real" model that people actually use,
I do use the formulation that you recommend.

Thank you though.  I believe that syntax makes a huge difference in
readability and long term error reduction.

Gene

On Fri, May 6, 2011 at 12:18 PM, Gavin Simpson wrote:

> On Fri, 2011-05-06 at 11:20 -0500, Gene Leynes wrote:
> > Hmmm
> >
> > After reading that email four times, I think I see what you mean.
> >
> > Checking for variables within particular scopes is probably one of the
> most
> > challenging things in R, and I would guess in other languages too.  In R
> > it's compounded by situations when you're writing a function to accept
> > variables as either a stand alone global variable, or as an element of a
> > data.frame or list.
>
> Dear Gene,
>
> > I think this is a new problem, and I'll switch to the lengthier
> > data.frame[,'x'] syntax in gam for now.
>
> No, No, No, No, No!
>
> If there is one thing you **should** take from this thread is that there
> is no need to perform subsetting like that in a model formula.
>
> Why would you want (or prefer) to do:
>
> gam(dat$y ~ s(dat$x))
>
> or
>
> gam(dat[, "y"] ~ s(dat[, "x"]))
>
> when
>
> gam(y ~ s(x), dat)
>
> will suffice?
>
> > By the way, about the $ accessors.  I can see why some people don't like
> > them, but they are a part of the language.  And, I think they're a good
> > part.  They make the code much more readable, and I have yet to make a
> > mistake using them (which makes me think that they can be used
> > responsibly).  Making code harder to read is its own source of error, and
> is
> > not something that I take lightly!
>
> If you use R's formula notation properly, you'll get cleaner code than
> either of your suggestions.
>
> > Thanks for the replies.  And, thank you Rolf for the detailed analysis.
>  Do
> > you think that your or I should submit a bug report to the package
> > maintainer?  I'm not sure how that works.  Very few of my challenges turn
> > out to be actual bugs, but I think this one is.
>
> I would consider this a bug - but in the sense that Simon didn't foresee
> the strange formulas that users of his software might concoct.
>
> By the way, I think you perhaps meant Berwin (re the detailed analysis)?
>
> HTH
>
> G
>
> > On Fri, May 6, 2011 at 4:53 AM, Berwin A Turlach
> > wrote:
> >
> > > G'day Rolf,
> > >
> > > On Fri, 06 May 2011 09:58:50 +1200
> > > Rolf Turner  wrote:
> > >
> > > > but it's strange that the dodgey code throws an error with gam(dat1$y
> > > > ~ s(dat1$x))  but not with gam(dat2$cf ~ s(dat2$s))
> > >
> > > > Something a bit subtle is going on; it would be nice to be able to
> > > > understand it.
> > >
> > > Well,
> > >
> > > R> traceback()
> > > 3: eval(expr, envir, enclos)
> > > 2: eval(inp, data, parent.frame())
> > > 1: gam(dat$y ~ s(dat$x))
> > >
> > > So the lines leading up to the problem seem to be the following from
> > > the gam() function:
> > >
> > >vars <- all.vars(gp$fake.formula[-2])
> > >inp <- parse(text = paste("list(", paste(vars, collapse = ","),
> > >")"))
> > >if (!is.list(data) && !is.data.frame(data))
> > >data <- as.data.frame(data)
> > >
> > >
> > >
> > > Setting
> > >
> > > R> options(error=recover)
> > >
> > > running the code until the error occurs, and then examining the frame
> > > number for the gam() call shows that "inp" is
> > > "expression(list( dat1,x ))" in your first example and
> > > "expression(list( dat2,s ))" in your second example.  In both
> > > examples, "data" is "list()" (not unsurprisingly).  When,
> > >
> > >dl <- eval(inp, data, parent.frame())
> > >
> > > is executed, it tries to eval "inp", in both cases "dat1" and "dat2"
> > > are found, obviously, in the parent frame.  In your first example "x"
> is
> > > (typically) not found and an error is thrown, in your second example an
> > > object with name "s" is found in "package:mgcv" and the call to eval
> > > succeeds.  "dl" becomes a list with two components, the first being,
> > > respectively, "dat1" or "dat2", and the second the body of the function
> > > "s".  (To verify that, you should probably issue the command
> > > "debug(gam)" and step through those first few lines of the function
> > > until you reach the above command.)
> > >
> > > The corollary is that you can use the name of any object that R will
> > > find in the parent frame, if it is another data set, then that data
> > > set will become the second component of "inp".  E.g.:
> > >
> > > R> dat=data.frame(min=1:100,cf=sin(1:100/50)+rnorm(100,0,.05))
> > > R> gam(dat$cf ~ s(dat$min))
> > >
> > > Family: gaussian
> > > Link function: identity
> > >
> > > Formula:
> > > dat$cf ~ s(dat$min)
> > >
> > > Estimated degr

[R] Default Colors in Script

2011-05-10 Thread Kishorenalluri
Dear Friends,
   I would appreciate one could help me, to prescribe the colors as
default if i have 10 more lines in my Y-Axis. 
I don't wanna set for all the column color like this


n_red <- subset(ts_data, my data>0, select=c(Fraction, column)
points(n_red$Fraction, n_red$Column, type = "h", col = "red")

I need default settings.

Thanks in advance


--
View this message in context: 
http://r.789695.n4.nabble.com/Default-Colors-in-Script-tp3512139p3512139.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] how to calculate the mean of a group in a table

2011-05-10 Thread tornanddesperate
Just wanted to thank everybody here for their great help. I tried the
"tapply" and it worked brilliantly. By changing its parameters I can easily
compute the values for the other treatments as well. Without help, I
probably would just have despaired.

Have a great time

mat

--
View this message in context: 
http://r.789695.n4.nabble.com/how-to-calculate-the-mean-of-a-group-in-a-table-tp3505986p3512145.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Weird read.xls behavior

2011-05-10 Thread Jun Shen
Hi, Gabor,

The interim fix solved the problem! Now it's working perfectly. And thanks
for everyone's input.

Jun

On Tue, May 10, 2011 at 10:39 AM, Gabor Grothendieck <
ggrothendi...@gmail.com> wrote:

> On Tue, May 10, 2011 at 12:12 AM, Jun Shen  wrote:
> > Kenneth,
> >
> > Thanks for the reply. I checked the original data. There is no space. I
> even
> > manually added a space to one value. After reading in with read.xls, the
> > value has two spaces. The reason I don't like it is I am going to do some
> > comparison with another dataset, which is supposed to be the same as this
> > one. Now I am getting a bunch of false negatives.
>
> It seems that the perl program underlying gdata's read.xls puts out
> lines like this:
>
> |"KAI-4169-002","830","5 mg" |
> where | mark the beginning and end and are not part of the line.
> read.csv includes the space after the last double quote in the last
> field even though its outside of the double quote.
>
> As an interim fix, edit the file at this location:
>
>   system.file("perl", "xls2csv.pl", package = "gdata")
>
> removing the space before the \n in this line:
>   print OutFile "$outputLine \n"
> so it becomes this:
>   print OutFile "$outputLine\n"
>
> Now it should work.
>
> --
> Statistics & Software Consulting
> GKX Group, GKX Associates Inc.
> tel: 1-877-GKX-GROUP
> email: ggrothendieck at gmail.com
>

[[alternative HTML version deleted]]

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


Re: [R] Weird read.xls behavior

2011-05-10 Thread David Winsemius


On May 10, 2011, at 11:39 AM, Gabor Grothendieck wrote:

On Tue, May 10, 2011 at 12:12 AM, Jun Shen   
wrote:

Kenneth,

Thanks for the reply. I checked the original data. There is no  
space. I even
manually added a space to one value. After reading in with  
read.xls, the
value has two spaces. The reason I don't like it is I am going to  
do some
comparison with another dataset, which is supposed to be the same  
as this

one. Now I am getting a bunch of false negatives.


It seems that the perl program underlying gdata's read.xls puts out
lines like this:


While we are on the topic of gdata functions I just looked at the trim  
function and find that it does return a data.frame when one is offered  
to it. (It was not clear from the documentation that a dataframe fell  
under the classification of "character strings and other related  
objects.")  Using a dataframe in my workspace from an earlier question:


> require(gdata)
> Mat1$Time[1] <- "09:30 "
> Mat1
  Weight Date   Time
17.6 04/28/11 09:30
28.4 04/29/11  03:11
38.6 04/29/11  05:32
48.6 04/29/11  09:53
51.4 05/01/11  19:52
> trim(Mat1)
  Weight Date  Time
17.6 04/28/11 09:30   # no space on my console
28.4 04/29/11 03:11
38.6 04/29/11 05:32
48.6 04/29/11 09:53
51.4 05/01/11 19:52
> nchar(trim(Mat1)$Time[1])
[1] 5
> nchar(Mat1$Time[1])
[1] 6

--
David.


|"KAI-4169-002","830","5 mg" |
where | mark the beginning and end and are not part of the line.
read.csv includes the space after the last double quote in the last
field even though its outside of the double quote.

As an interim fix, edit the file at this location:

  system.file("perl", "xls2csv.pl", package = "gdata")

removing the space before the \n in this line:
  print OutFile "$outputLine \n"
so it becomes this:
  print OutFile "$outputLine\n"

Now it should work.

--
Statistics & Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.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.


David Winsemius, MD
West Hartford, CT

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] need to delete by time, not date

2011-05-10 Thread Adele_Thompson
Good deal. I will do that in the future. Thanks.

-Original Message-
From: dwinsem...@comcast.net [mailto:dwinsem...@comcast.net] 
Sent: Tuesday, May 10, 2011 10:45 AM
To: Thompson, Adele - adele_thomp...@cargill.com
Cc: r-help@r-project.org
Subject: Re: [R] need to delete by time, not date


On May 10, 2011, at 11:23 AM, adele_thomp...@cargill.com wrote:

> Thanks, that did exactly what I wanted. I realize that I didn't put  
> it in the best form, but was having trouble coming up with a  
> replicable example. I see that many people here quickly create  
> examples and I am learning how to do that, but still struggle,  
> especially when using as.POSIXct and other more complex classes.

Excuse not accepted. All you needed to do with Mat1 was dput(Mat1). Or  
if your real data object was a lot longer, then just dput(head(Mat1,  
10)).

-- 
David.
> Thanks for the help.
>
> -Original Message-
> From: dwinsem...@comcast.net [mailto:dwinsem...@comcast.net]
> Sent: Tuesday, May 10, 2011 10:19 AM
> To: Thompson, Adele - adele_thomp...@cargill.com
> Cc: r-help@r-project.org
> Subject: Re: [R] need to delete by time, not date
>
>
> On May 10, 2011, at 11:06 AM, Schatzi wrote:
>
>> I have a matrix where one column has a date and another column has a
>> time. I
>> need to delete all times before 6am. I had combined the Date and
>> Time column
>> into DateTime.
>>
>> Mat1:
>> Weight Date  Time
>> 7.6 04/28/11 09:03
>> 8.4 04/29/11 03:11
>> 8.6 04/29/11 05:32
>> 8.6 04/29/11 09:53
>> 1.4 05/01/11 19:52
>>
>> I tried just picking up the time:
>>
>> as.POSIXct(Mat1$Time, format = "%H:%M")
>>
>> but this assigns the time value today's date.
>
> Since these are either character strings (or factors) at the
> moment ... and if they are matrices they are character...
>
> Mat1[Mat1$Time >= "06:00" , ]
>
> If actually a dataframe then
>
> Mat1[as.character(Mat1$Time) >= "06:00" , ]
>
> In the future please include the output of dput(Mat1).
>
>>
>> The final matrix should look like:
>> Weight Date  Time
>> 7.6 04/28/11 09:03
>> 8.6 04/29/11 09:53
>> 1.4 05/01/11 19:52
>>
>>
>
> David Winsemius, MD
> West Hartford, CT
>

David Winsemius, MD
West Hartford, CT

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] fitting non-intercept model with lrm

2011-05-10 Thread David Winsemius


On May 10, 2011, at 11:38 AM, Clarissa Valim wrote:

I would appreciate if someone could tell me how to fit a non- 
intercept model
using lrm (and not glm). The -1 in the formula of the glm does not  
work with

lrm.


Can you describe why one might want to have such a facility? I cannot  
think that it should change the value or interpretation of the  
remaining coefficients.


--

David Winsemius, MD
West Hartford, CT

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

2011-05-10 Thread Alaios
Dear all,
today I have writted the following code,
to plot the contents of some matrices I have

plot_shad_f
function(f){
  library(ggplot2)
  dev.new()
  plotdata<-melt(f)
  names(plotdata)<-c('x','y','z')
  v<-ggplot(plotdata, aes(x, y, z = z))
  print(v + geom_tile(aes(fill=z))) 
}

I would like to ask your help add a small circle in this plotting. What would 
be the easiest way to do that in ggplot2?

Best Regards
Alex

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

2011-05-10 Thread Scott Chamberlain
You should provide reproducible data in addition to your code. 

S
On Tuesday, May 10, 2011 at 11:54 AM, Alaios wrote: 
> Dear all,
> today I have writted the following code,
> to plot the contents of some matrices I have
> 
> plot_shad_f
> function(f){
>  library(ggplot2)
>  dev.new()
>  plotdata<-melt(f)
>  names(plotdata)<-c('x','y','z')
>  v<-ggplot(plotdata, aes(x, y, z = z))
>  print(v + geom_tile(aes(fill=z))) 
> }
> 
> I would like to ask your help add a small circle in this plotting. What would 
> be the easiest way to do that in ggplot2?
> 
> Best Regards
> Alex
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 

[[alternative HTML version deleted]]

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


Re: [R] ggplot2 and add circle

2011-05-10 Thread Alaios
Here you are :)


plot_shad_f<-function(f){
  library(ggplot2)
  dev.new()
  plotdata<-melt(f)
  names(plotdata)<-c('x','y','z')
  v<-ggplot(plotdata, aes(x, y, z = z))
  print(v + geom_tile(aes(fill=z))) 

}


f<-matrix(data=seq(1:1),nrow=100,ncol=100)
plot_shad_f(f)


I would like to add a circle at the middle of this region with a range of 10.

Best Regards

Alex

--- On Tue, 5/10/11, Scott Chamberlain  wrote:

From: Scott Chamberlain 
Subject: Re: [R] ggplot2 and add circle
To: "Alaios" 
Cc: R-help@r-project.org
Date: Tuesday, May 10, 2011, 5:59 PM



You should provide reproducible data in addition to your code. 
S





On Tuesday, May 10, 2011 at 11:54 AM, Alaios wrote:

Dear all,
today I have writted the following code,
to plot the contents of some matrices I have

plot_shad_f
function(f){
  library(ggplot2)
  dev.new()
  plotdata<-melt(f)
  names(plotdata)<-c('x','y','z')
  v<-ggplot(plotdata, aes(x, y, z = z))
  print(v + geom_tile(aes(fill=z))) 
}

I would like to ask your help add a small circle in this plotting. What would 
be the easiest way to do that in ggplot2?

Best Regards
Alex

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













[[alternative HTML version deleted]]

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


[R] Converting ordinal dates and time into sensible formats

2011-05-10 Thread Sam Albers
Hello all,

I am having a little trouble working with "strptime" and I was hoping
someone might be able to give me a hand. I have an instrument that outputs
an ordinal date and time in two columns something like this:

  day.hour min.sec
1125252050
2125182029
3125242023
4125242028
5125072035

Now the problem I am having is converting these numbers into dates and
times. I am able to convert these into their respective POSIXlt formats but
I am left with two columns where hour is left with the date (data$Only.Date)
and the date is left with the time (data$Only.Time). Can anyone recommend a
good way to convert these ordinal dates into something like the following?

day.hour min.secDateTime
12511 20332011-05-05   11:20:33

## A trivial example
##a data frame
day.hour <-as.integer(runif(5, 12500, 12523)) #First 3 digits are the day of
the year, last 2 are the hour of the day
data <- as.data.frame(day.hour)
data$min.sec <-as.integer(runif(5, 2000, 2060)) #First 2 digits are the
minute, last 2 are the seconds

##example of how things get a little jumbled. strptime was easy enough to
use.
data$Date <- strptime(data$day.hour, format="%j%H")
data$Time <- strptime(data$min.sec, format="%M%S")
data

Using Ubuntu 10.10 and R 2.11.1. Thanks in advance

Sam

-- 
*
Sam Albers
Geography Program
University of Northern British Columbia
 University Way
Prince George, British Columbia
Canada, V2N 4Z9
*

[[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] Dual colour ramps based on pos/neg values

2011-05-10 Thread Tyler Hayes
My apologies for the late reply but I was out of town for a while. The
solution I wound up using is below. Sorry about the text if it didn't
wrap well. You should be able to pick out the code from the comments
though.

Thanks for all the help!

Cheers,

t.


## > Start hack
##
##*** OLD COLOUR STYLES ***
## r2b <- c("#0571B0", "#92C5DE", "#F7F7F7", "#F4A582", "#CA0020") #red to blue
## r2g <- c("Red", "DarkRed", "Green", "Chartreuse")
## w2b <- c("#045A8D", "#2B8CBE", "#74A9CF", "#BDC9E1", "#F1EEF6")
#white to blue
## assign("col.sty", get(color))
## calendar.pal <- colorRampPalette((col.sty), space = "Lab")
##
## *** NEW METHOD FOR ZERO CROSSING RAMP SCALES ***
##
## First low colour is the MIN value; second is closest to 0.0
## First high is the closest to 0.0; second is the MAX value.
## If names are not known to the graphics driver, use HEX values
## I tend to put the extremes as DARK and the low values as MUTED/LIGHT colours
lowColFun  <- colorRampPalette(c("#80","#FF","#FF82AB","#FFE4E1"),
space = "Lab")
highColFun <- colorRampPalette(c("#BDFCC9","#7FFF00","#00EE00","#008000"),
space = "Lab")
## These are hard coded, but should be made available to tune
scaleFac   <- 1.001
ncolors<- 99  # Should be odd for now b/c of tmpXEnd calculations below
## Define middle cutoff values, this could also be arbitrary for
## highlighting specific regions
tmpMid <- c(-0.0001,0.0001)
## Cuts where to put the colors
tmpLowEnd  <- 
seq(from=min(values)*scaleFac,to=tmpMid[1]*scaleFac,length=((ncolors-1)/2))
tmpHighEnd <- 
seq(from=tmpMid[2]*scaleFac,to=max(values)*scaleFac,length=((ncolors-1)/2))
tmpAtVals  <- c(tmpLowEnd,tmpMid,tmpHighEnd)
## Create final values for levelplot
my.at  <- c(tmpLowEnd,tmpMid,tmpHighEnd)
my.col.reg <- c(lowColFun(length(tmpLowEnd)),
rep("black",length(tmpMid)), highColFun(length(tmpHighEnd)) )
my.cuts<- length(my.col.reg)-1
##
## > End hack

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

2011-05-10 Thread Scott Chamberlain
Here is a thread that might help. 

http://groups.google.com/group/ggplot2/browse_thread/thread/7f4238b5658e85bb/25d4cfb6215b7544?lnk=gst&q=circle+size#25d4cfb6215b7544
On Tuesday, May 10, 2011 at 12:15 PM, Alaios wrote: 
> Here you are :)
> 
> 
> plot_shad_f<-function(f){
>  library(ggplot2)
>  dev.new()
>  plotdata<-melt(f)
>  names(plotdata)<-c('x','y','z')
>  v<-ggplot(plotdata, aes(x, y, z = z))
>  print(v + geom_tile(aes(fill=z))) 
> 
> }
> 
> 
> f<-matrix(data=seq(1:1),nrow=100,ncol=100)
> plot_shad_f(f)
> 
> 
> I would like to add a circle at the middle of this region with a range of 10.
> 
> Best Regards
> 
> Alex
> 
> --- On Tue, 5/10/11, Scott Chamberlain  wrote:
> > 
> > From: Scott Chamberlain 
> > Subject: Re: [R] ggplot2 and add circle
> > To: "Alaios" 
> > Cc: R-help@r-project.org
> > Date: Tuesday, May 10, 2011, 5:59 PM
> > 
> > You should provide reproducible data in addition to your code. 
> > 
> > S
> > On Tuesday, May 10, 2011 at 11:54 AM, Alaios wrote:
> > > Dear all,
> > > today I have writted the following code,
> > > to plot the contents of some matrices I have
> > > 
> > > plot_shad_f
> > > function(f){
> > >  library(ggplot2)
> > >  dev.new()
> > >  plotdata<-melt(f)
> > >  names(plotdata)<-c('x','y','z')
> > >  v<-ggplot(plotdata, aes(x, y, z = z))
> > >  print(v + geom_tile(aes(fill=z))) 
> > > }
> > > 
> > > I would like to ask your help add a small circle in this plotting. What 
> > > would be the easiest way to do that in ggplot2?
> > > 
> > > Best Regards
> > > Alex
> > > 
> > > __
> > > R-help@r-project.org mailing list
> > > https://stat.ethz.ch/mailman/listinfo/r-help
> > > PLEASE do read the posting guide 
> > > http://www.R-project.org/posting-guide.html
> > > and provide commented, minimal, self-contained, reproducible code.
> > > 
> > 

[[alternative HTML version deleted]]

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


Re: [R] ggplot2 and add circle

2011-05-10 Thread Dennis Murphy
Hi:

Here's one way:

plot_shad <- function(d, r) {
   require('ggplot2')
   plotdata <- melt(d)
   names(plotdata)<-c('x','y','z')
   xc <- mean(range(plotdata$x))
   yc <- mean(range(plotdata$y))
   theta <- seq(-pi, pi, length = 200)
   circ <- data.frame(xv = xc + r * cos(theta),
  yv = yc + r * sin(theta))
   v <- ggplot(plotdata)
   print(v + geom_tile(aes(x = x, y = y, fill = z)) +
 geom_path(data = circ, aes(x = xv, y = yv), color = 'white',
size = 1) +
 coord_equal()
)
 }

plot_shad(f, 10)

HTH,
Dennis

On Tue, May 10, 2011 at 10:15 AM, Alaios  wrote:
> Here you are :)
>
>
> plot_shad_f<-function(f){
>   library(ggplot2)
>   dev.new()
>   plotdata<-melt(f)
>   names(plotdata)<-c('x','y','z')
>   v<-ggplot(plotdata, aes(x, y, z = z))
>   print(v + geom_tile(aes(fill=z)))
>
> }
>
>
> f<-matrix(data=seq(1:1),nrow=100,ncol=100)
> plot_shad_f(f)
>
>
> I would like to add a circle at the middle of this region with a range of 10.
>
> Best Regards
>
> Alex
>
> --- On Tue, 5/10/11, Scott Chamberlain  wrote:
>
> From: Scott Chamberlain 
> Subject: Re: [R] ggplot2 and add circle
> To: "Alaios" 
> Cc: R-help@r-project.org
> Date: Tuesday, May 10, 2011, 5:59 PM
>
>
>
>                You should provide reproducible data in addition to your code.
> S
>
>
>
>
>
>                On Tuesday, May 10, 2011 at 11:54 AM, Alaios wrote:
>
>                    Dear all,
> today I have writted the following code,
> to plot the contents of some matrices I have
>
> plot_shad_f
> function(f){
>  library(ggplot2)
>  dev.new()
>  plotdata<-melt(f)
>  names(plotdata)<-c('x','y','z')
>  v<-ggplot(plotdata, aes(x, y, z = z))
>  print(v + geom_tile(aes(fill=z)))
> }
>
> I would like to ask your help add a small circle in this plotting. What would 
> be the easiest way to do that in ggplot2?
>
> Best Regards
> Alex
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>
>
>
>
>
>
>
>
>
>
>
>
>        [[alternative HTML version deleted]]
>
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>

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


[R] for loop populate matrix

2011-05-10 Thread Chris English

Dear List,
I am trying to populate a matrix using a for loop.This works.
> four<- matrix(nrow=4,ncol=5)> for (j in 1:ncol(four)){+ for (i in 
> 1:nrow(four)){+ four[i,j]<-i-1}+ }> print(four)     [,1] [,2] [,3] [,4] 
> [,5][1,]    0    0    0    0    0[2,]    1    1    1    1    1[3,]    2    2  
>   2    2    2[4,]    3    3    3    3    3
What I would like is the inverse ofthe above, i.e.: [,1] [,2] [,3] [,4] 
[,5][1,]     3    3    3    3    3[2,]     2    2    2    2    2[3,]     1    1 
   1    1    1[4,]     0    0    0    0    0
and can't seem to get there.  Thanks,Chris  
  
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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 a large netCDF file using R

2011-05-10 Thread Sulochan Dhungel
Hi Dave,

I asked for the climate group for one particular day's data so that I could
match it up with the data I got from the codes for the same day. It did not
match. So I knew it was not working. Also time is the first thing on the
list of summary for "nc" , but when I saw other netCDF files, time was the
last thing.

Thanks again for looking into my problem.

Sulochan

On Tue, May 10, 2011 at 7:50 AM, David William Pierce wrote:

> On Mon, May 9, 2011 at 6:29 PM, sulochandhungel  > wrote:
>
>> I got a netCDF file from a climate group and wanted to aggregate the data.
>> The file summary looks like this
>> [...]
>>
>> I didnt think it worked... can u suggest me some other way ... I could not
>> use the start and count properly
>
>
>  Could you be a little more explicit? Why don't you think it worked?
>
> In general, you seem to be using count and start correctly.
>
> --Dave
>
>
> --
> David W. Pierce
> Division of Climate, Atmospheric Science, and Physical Oceanography
> Scripps Institution of Oceanography
> (858) 534-8276 (voice)  /  (858) 534-8561 (fax)dpie...@ucsd.edu
>
>

[[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] for loop populate matrix

2011-05-10 Thread Dennis Murphy
Why loop?

matrix(rep(3:0, each = 5), ncol = 5, byrow = TRUE)
 [,1] [,2] [,3] [,4] [,5]
[1,]33333
[2,]22222
[3,]11111
[4,]00000

HTH,
Dennis

On Tue, May 10, 2011 at 11:10 AM, Chris English  wrote:
>
> Dear List,
> I am trying to populate a matrix using a for loop.This works.
>> four<- matrix(nrow=4,ncol=5)> for (j in 1:ncol(four)){+ for (i in 
>> 1:nrow(four)){+ four[i,j]<-i-1}+ }> print(four)     [,1] [,2] [,3] [,4] 
>> [,5][1,]    0    0    0    0    0[2,]    1    1    1    1    1[3,]    2    2 
>>    2    2    2[4,]    3    3    3    3    3
> What I would like is the inverse ofthe above, i.e.:     [,1] [,2] [,3] [,4] 
> [,5][1,]     3    3    3    3    3[2,]     2    2    2    2    2[3,]     1    
> 1    1    1    1[4,]     0    0    0    0    0
> and can't seem to get there.  Thanks,Chris
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] Dual colour ramps based on pos/neg values

2011-05-10 Thread Achim Zeileis

On Tue, 10 May 2011, Tyler Hayes wrote:


My apologies for the late reply but I was out of town for a while. The
solution I wound up using is below. Sorry about the text if it didn't
wrap well. You should be able to pick out the code from the comments
though.


Note that despite computing the ramp in LAB space, you still get very 
flashy and unbalanced palettes that are hard to decode for color-blind 
viewers. The reason is that you start out from very flashy and unbalanced 
colors as input values for your color ramps. This is illustrated in the 
following:


  ## your two parts of the diverging palette
  lowColFun <- colorRampPalette(
c("#80","#FF","#FF82AB","#FFE4E1"), space = "Lab")
  highColFun <- colorRampPalette(
c("#BDFCC9","#7FFF00","#00EE00","#008000"), space = "Lab")

  ## choose small set of values
  cols <- c(lowColFun(4), highColFun(4))

  ## mapped to dichromatic vision
  library("dichromat")
  cols1 <- dichromat(cols, "deutan")
  cols2 <- dichromat(cols, "protan")

  ## desaturated version (computed in HCL space)
  library("colorspace")
  desaturate <- function(x) {
x <- as(hex2RGB(cols), "polarLUV")
x@coords[,2] <- 0
hex(x)
  }
  cols0 <- desaturate(cols)

  ## display all four sets of colors
  pal <- function(col, border = "light gray") {
n <- length(col)
plot(0, 0, type="n", xlim = c(0, 1), ylim = c(0, 1), axes = FALSE,
  xlab = "", ylab = "")
rect(0:(n-1)/n, 0, 1:n/n, 1, col = col, border = border)
  }
  par(mfrow = c(4, 1), oma = rep(1, 4), mar = rep(0, 4))
  pal(cols)
  pal(cols0)
  pal(cols1)
  pal(cols2)

The first panel of the original colors shows that the colors are very 
flashy and hence hard to look at for a longer time. The second panel shows 
that the colors are unbalanced, i.e., the red colors are darker than the 
green colors. And the third and fourth panel shows that the colors will be 
hard to decode properly for dichromatic viewers.


The "Escaping RGBland" paper that I posted in my previous reply discusses 
how all of these problems can be mitigated. The result is provided in 
diverge_hcl() in "colorspace".



Thanks for all the help!

Cheers,

t.


## > Start hack
##
##*** OLD COLOUR STYLES ***
## r2b <- c("#0571B0", "#92C5DE", "#F7F7F7", "#F4A582", "#CA0020") #red to blue
## r2g <- c("Red", "DarkRed", "Green", "Chartreuse")
## w2b <- c("#045A8D", "#2B8CBE", "#74A9CF", "#BDC9E1", "#F1EEF6")
#white to blue
## assign("col.sty", get(color))
## calendar.pal <- colorRampPalette((col.sty), space = "Lab")
##
## *** NEW METHOD FOR ZERO CROSSING RAMP SCALES ***
##
## First low colour is the MIN value; second is closest to 0.0
## First high is the closest to 0.0; second is the MAX value.
## If names are not known to the graphics driver, use HEX values
## I tend to put the extremes as DARK and the low values as MUTED/LIGHT colours
lowColFun  <- colorRampPalette(c("#80","#FF","#FF82AB","#FFE4E1"),
space = "Lab")
highColFun <- colorRampPalette(c("#BDFCC9","#7FFF00","#00EE00","#008000"),
space = "Lab")
## These are hard coded, but should be made available to tune
scaleFac   <- 1.001
ncolors<- 99  # Should be odd for now b/c of tmpXEnd calculations below
## Define middle cutoff values, this could also be arbitrary for
## highlighting specific regions
tmpMid <- c(-0.0001,0.0001)
## Cuts where to put the colors
tmpLowEnd  <- 
seq(from=min(values)*scaleFac,to=tmpMid[1]*scaleFac,length=((ncolors-1)/2))
tmpHighEnd <- 
seq(from=tmpMid[2]*scaleFac,to=max(values)*scaleFac,length=((ncolors-1)/2))
tmpAtVals  <- c(tmpLowEnd,tmpMid,tmpHighEnd)
## Create final values for levelplot
my.at  <- c(tmpLowEnd,tmpMid,tmpHighEnd)
my.col.reg <- c(lowColFun(length(tmpLowEnd)),
rep("black",length(tmpMid)), highColFun(length(tmpHighEnd)) )
my.cuts<- length(my.col.reg)-1
##
## > End hack

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] [caret package] [trainControl] supplying predefined partitions to train with cross validation

2011-05-10 Thread Fabon Dzogang
Here is an answer from Max Khun thank you !

Fabon,

If I understand the problem, there are two ways of doing it. First, if you
are using caret's trian(), rfe() or sbf(), if you set the seed right before
you call the models, they end up using the same resampled data sets. (btw,
if you use the resamples() function in caret, it checks for the same
resampling indices)

If you want to manually fix the data sets, there is an example in section
5.2 of

 http://cran.r-project.org/web/packages/caret/vignettes/caretTrain.pdf

Using LGOCV. For 10-fold CV, you can use createFolds() with an additional
argument:

> createFolds(1:10, returnTrain = TRUE)
$Fold01
[1] 2 3 4 5 6 7 8 9 10

$Fold02
[1]  1  3  4  5  6  7  8  9 10

$Fold03
[1]  1  2  4  5  6  7  8  9 10

$Fold04
[1]  1  2  3  5  6  7  8  9 10

$Fold05
[1]  1  2  3  4  6  7  8  9 10

$Fold06
[1]  1  2  3  4  5  7  8  9 10

$Fold07
[1]  1  2  3  4  5  6  8  9 10

$Fold08
[1]  1  2  3  4  5  6  7  9 10

$Fold09
[1]  1  2  3  4  5  6  7  8 10

$Fold10
[1] 1 2 3 4 5 6 7 8 9

For the trainControl() function, the index argument should be a list of
samples indices for each resample. So if I give it the the above results of
createFolds(), it will do 10-fold cv.

MAx

On Fri, May 6, 2011 at 12:32 PM, Fabon Dzogang  wrote:
> Hello,
>
> Thank you for your reply but I'm not sure your code answers my needs,
> from what I read it creates a 10-fold partition and then extracts the
> kth partition for future processing.
>
> My question was rather: once I have a 10-fold partition of my data,
> how to supply it to the "train" function of the caret package. Here's
> some sample code :
>
> folds <- createFolds(my_dataset_classes, 10)
>
> # I can't use index=folds on this one, it will train on the 1/k and test on 
> k-1
> t_control <- trainControl(method="cv", number=10)
>
> # here I would like train to take account of my predefined folds
> model <- train(my_dataset_predictors, my_dataset_classes,
> method="svmLinear", trControl = t_control)
>
> Cheers,
> Fabon.
>
> On Fri, May 6, 2011 at 10:59 AM, neetika nath  wrote:
>> Hi,
>> I did the similar experiment with my data. may be following code will give
>> you some idea. It might not be the best solution but for me it worked.
>> please do share if you get other idea.
>> Thank you
>>  CODE###
>>
>> library(dismo)
>>
>> set.seed(111)
>>
>> dd<-read.delim("yourfile.csv",sep=",",header=T)
>>
>> # To keep a check on error
>>
>> options(error=utils::recover)
>>
>> # dd- data to be split for 10 Fold CV, this will split complete data into 10
>> fold
>>
>> number<-kfold(dd, k=10)
>>
>> case 1: if k ==1
>>
>> x<-NULL;
>>
>> #retrieve all the index (from your data) for 1st fold in x, such that you
>> can use it as a test set and remaining can be used as train set for #1st
>> iteration.
>>
>> x<-which(number==k)
>>
>> On Thu, May 5, 2011 at 11:43 PM, Fabon Dzogang 
>> wrote:
>>>
>>> Hi all,
>>>
>>> I run R 2.11.1 under ubuntu 10.10 and caret version 2.88.
>>>
>>> I use the caret package to compare different models on a dataset. In
>>> order to compare their different performances I would like to use the
>>> same data partitions for every models. I understand that using a LGOCV
>>> or a boot type re-sampling method along with the "index" argument of
>>> the trainControl function, one is able to supply a training partition
>>> to the train function.
>>>
>>> However, I would like to apply a 10-fold cross validation to validate
>>> the models and I did not find any way to supply some predefined
>>> partition (created with createFolds) in this setting. Any help ?
>>>
>>> Thank you and great package by the way !
>>>
>>> Fabon Dzogang.
>>>
>>> __
>>> R-help@r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>> PLEASE do read the posting guide
>>> http://www.R-project.org/posting-guide.html
>>> and provide commented, minimal, self-contained, reproducible code.
>>
>>
>
>
>
> --
> Fabon Dzogang
>



-- 
Fabon Dzogang

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] maximum likelihood convergence reproducing Anderson Blundell 1982 Econometrica R vs Stata

2011-05-10 Thread Berend Hasselman

Alex Olssen wrote:
> 
> Wow that is really interesting,
> 
> Sorry I was asleep when you emailed these.
> 
> And yes, of course, I had been trying to implement model 18, not 18s,
> that was a typo, sorry.
> 
> I will have a look at the code you posted.
> 
> Thanks,
> 
> Alex
> 

I have run nlm with a zero starting vector on model 18.


> start.zero <- rep(0,8)
> 
> q.zero <- nlm(p=start.zero, f=lnl, hessian=TRUE, gradtol=1e-15,
+  y1=y1, y2=y2, x1=x1, x2=x2, x3=x3)
> 
> q.zero
$minimum
[1] -76.73878

$estimate
[1] -0.68124025  0.28686058 -0.15394903 -0.06173779  1.30304998 -0.25446946
[7]  0.14531856  0.05349810

$gradient
[1]  0.0001356761  0.0014407493  0.0013240913  0.0006367333 -0.0005172628
[6] -0.0051911383 -0.0044429086 -0.0019081432

$hessian
  [,1] [,2] [,3]  [,4]  [,5] [,6] [,7]
[1,]  197993.1  1982843  1746582  807261.9  135628.5  1355944  1194854
[2,] 1982842.8 19766730 17430209 8084282.7 1357883.8 13506353 11915748
[3,] 1746581.8 17430209 15376070 7138194.2 1196134.9 11911149 10512490
[4,]  807261.9  8084283  7138194 3323204.7  552932.9  5526861  4882154
[5,]  135628.5  1357884  1196135  552932.9  271571.7  2713953  2391658
[6,] 1355943.5 13506353 11911149 5526860.9 2713953.0 26956003 23785984
[7,] 1194854.2 11915748 10512490 4882154.2 2391658.1 23785984 20996480
[8,]  552860.7  5534816  4887292 2275676.3 1106862.5 11057031  9768028
   [,8]
[1,]   552860.7
[2,]  5534816.0
[3,]  4887292.4
[4,]  2275676.3
[5,]  1106862.5
[6,] 11057030.7
[7,]  9768027.6
[8,]  4554570.6

$code
[1] 2

$iterations
[1] 68

> 
> q.zero$estimate
[1] -0.68124025  0.28686058 -0.15394903 -0.06173779  1.30304998 -0.25446946
[7]  0.14531856  0.05349810

Quite different from what you get if one starts with AB's model 18 results
of Table II.
I have also used Gretl's FIML procedure to estimate the model (I don't know
what starting values Gretl uses).
The results are very similar to those obtained with nlm as above.

Berend

--
View this message in context: 
http://r.789695.n4.nabble.com/maximum-likelihood-convergence-reproducing-Anderson-Blundell-1982-Econometrica-R-vs-Stata-tp3502516p3512807.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] specifying scales in lattice xyplot makes the lines disappear?

2011-05-10 Thread Christopher W Ryan
I have a dataframe concerning manner of death from death certificates,
from 2005 to 2009 inclusive, with the following structure:

> str(MannerYoung.plot.data)
'data.frame':   245 obs. of  4 variables:
 $ year : Factor w/ 5 levels "2005","2006",..: 1 2 3 4 5 1 2 3 4 5 ...
 $ manner   : Factor w/ 7 levels "Accident","Homicide",..: 1 1 1 1
1 2 2 2 2 2 ...
 $ new.age.group: Factor w/ 7 levels "1 to 4","5 to 9",..: 1 1 1 1 1 1
1 1 1 1 ...
 $ Freq : int  [omitted from this post]...
>

new.age.group categorizes age of death in 5-year age bands.

The following makes a nice lattice scatterplot, with a panel for each
age-band, year on the horizontal axis,and frequency on the veritcal
axis.

xyplot(Freq~year | new.age.group, groups=manner,
data=MannerYoung.plot.data, type="l",  auto.key =
list(text=levels(MannerYoung.plot.data$manner), lines=TRUE,
points=FALSE, space = "top"),
layout=c(length(levels(MannerYoung.plot.data$new.age.group)),1))

Except that the year tick labels overlap each other. I thought I would
try to label only every other tick, using scales(), like this:

xyplot(Freq~year | new.age.group, groups=manner,
data=MannerYoung.plot.data, type="l",
scales=list(x=list(at=c(2005,2007,2009),
labels=c("2005","2007","2009"))), auto.key =
list(text=levels(MannerYoung.plot.data$manner), lines=TRUE,
points=FALSE, space = "top"),
layout=c(length(levels(MannerYoung.plot.data$new.age.group)),1))

This draws the panels, vertical axis, horizontal axes, with the
odd-numbered years, everything as desired, except that there are no
lines in the panels--no actual content.

I suppose I could decrease the font size for the year labels, using
cex(), and that cured the overlap, but I wanted to understand
scales().

Where am I going wrong?

Thanks.

--Chris Ryan

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

2011-05-10 Thread Thomas Lumley
On Tue, May 10, 2011 at 2:50 PM, jour4life  wrote:
> Hello all,
>
> I am wondering if there is a way to specify sampling weights for an ols
> model using sample weights.
>
> For instance, right now, my code is:
>
> fit.ex<-lm(y~x1+x2+x3+...xk,data=dataset,weights=weightvariable.)
> summary(fit.ex)
>
> But, there is almost no difference in the coefficients nor standard errors.
> I am skeptical that I am using the option correctly since many posts state
> that the "weights" option is for weighted least squares. Or is this the same
> thing?

You want svyglm() in the survey package.

-thomas

-- 
Thomas Lumley
Professor of Biostatistics
University of Auckland

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

2011-05-10 Thread Thomas Mang

Hi,

Consider the need for a  regression model which can handle an ordered 
multinomial response variable. There are, for example, proportional odds 
/ cumulative logit models, but actually the regression should include 
random effects (a mixed model), and I would not be aware of multinomial 
regression model as part of lme4 (am I wrong here ?). Further, the 
constraint of proportional odd models that predictors have the same 
relative impact across all levels, does most likely not hold for the 
study in question.


I was wondering if an ordinary binomial mixed model can be turned in an 
multinomial one through preparing the input data.frame in a different way:
Consider three response levels, A, B, C, ordered. I can accurately 
describe the occurrence of each of these three realizations using one to 
two Bernoulli random variables:


Let
P(X == A) = a
P(X in {B, C}) = 1 - a
P(X == B | X in {B, C}) = (1 - a) * b
P(X == C | X in {B, C}) = (1 - a) * (1 - b)

so the first comparison checks if A or either of B/C is the case, and 
the second, conditional on it's either B/C, checks which of these two 
holds. Sort of traversing sequentially the hierarchy of the ordered levels.
In terms of the likelihood of the desired model, the probabilities on 
the right hand side would be exactly achieved if I use one input row in 
case the random variable takes on the value A and assign the response 
variable the value 0, while in the other cases the probabilities are 
achieved by using two input table rows, with the first one having value 
1 for the response variable so the random variate is either B/C) and a 
second row with response equal to 0 if B is the case, and 1 otherwise, 
that is C is the case.


Certainly, degrees of freedom must be manually adjusted in inferences, 
as every measured response should provide only a single degree of freedom.


Question: Do I overlook here something, or is above outlined way a valid 
method to yield an ordered multinomial mixed model by tweaking the input 
table in such manners ?


many thanks and best,
Thomas

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


Re: [R] for loop populate matrix

2011-05-10 Thread Chris English

Indeed! Your example works perfectly.
I think I was using the for-loop to try to address cases where thesize of the 
matrix was unknown. But your suggestion works forarbitrary sizes, just have to 
know the range of values to rep (repeat) and have"each" and "ncol" equal.
Thank you.ChrisMy apologies to Dennis for direct reply, I'm not used to using 
'reply all': tales of havoc wrought

> Date: Tue, 10 May 2011 11:19:51 -0700
> Subject: Re: [R] for loop populate matrix
> From: djmu...@gmail.com
> To: sgl...@hotmail.com
> CC: r-help@r-project.org
>
> Why loop?
>
> matrix(rep(3:0, each = 5), ncol = 5, byrow = TRUE)
> [,1] [,2] [,3] [,4] [,5]
> [1,] 3 3 3 3 3
> [2,] 2 2 2 2 2
> [3,] 1 1 1 1 1
> [4,] 0 0 0 0 0
>
> HTH,
> Dennis
>
> On Tue, May 10, 2011 at 11:10 AM, Chris English  wrote:
> >
> > Dear List,
> > I am trying to populate a matrix using a for loop.This works.
> >> four<- matrix(nrow=4,ncol=5)> for (j in 1:ncol(four)){+ for (i in 
> >> 1:nrow(four)){+ four[i,j]<-i-1}+ }> print(four)á á á[,1] [,2] [,3] [,4] 
> >> [,5][1,] á á0 á á0 á á0 á á0 á á0[2,] á á1 á á1 á á1 á á1 á á1[3,] á á2 á 
> >> á2 á á2 á á2 á á2[4,] á á3 á á3 á á3 á á3 á á3
> > What I would like is the inverse ofthe above, i.e.: á á [,1] [,2] [,3] [,4] 
> > [,5][1,] á á 3 á á3 á á3 á á3 á á3[2,] á á 2 á á2 á á2 á á2 á á2[3,] á á 1 
> > á á1 á á1 á á1 á á1[4,] á á 0 á á0 á á0 á á0 á á0
> > and can't seem to get there. áThanks,Chris
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/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] sample weights in ols

2011-05-10 Thread jour4life
On 5/10/2011 3:12 PM, Thomas Lumley-2 [via R] wrote:
> On Tue, May 10, 2011 at 2:50 PM, jour4life <[hidden email] 
> > wrote:
>
> > Hello all,
> >
> > I am wondering if there is a way to specify sampling weights for an ols
> > model using sample weights.
> >
> > For instance, right now, my code is:
> >
> > fit.ex<-lm(y~x1+x2+x3+...xk,data=dataset,weights=weightvariable.)
> > summary(fit.ex)
> >
> > But, there is almost no difference in the coefficients nor standard 
> errors.
> > I am skeptical that I am using the option correctly since many posts 
> state
> > that the "weights" option is for weighted least squares. Or is this 
> the same
> > thing?
>
> You want svyglm() in the survey package.
>
> -thomas
>
> -- 
> Thomas Lumley
> Professor of Biostatistics
> University of Auckland
>
> __
> [hidden email] 
>  mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>
> 
> If you reply to this email, your message will be added to the 
> discussion below:
> http://r.789695.n4.nabble.com/sample-weights-in-ols-tp3510865p3512894.html 
>
> To unsubscribe from sample weights in ols, click here 
> .
>  
>
Thanks very much!!

Carlos


--
View this message in context: 
http://r.789695.n4.nabble.com/sample-weights-in-ols-tp3510865p3512975.html
Sent from the R help mailing list archive at Nabble.com.
[[alternative HTML version deleted]]

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


Re: [R] ANOVA 1 too few degrees of freedom

2011-05-10 Thread Rovinpiper
So what is the difference between a colon and an asterisk in this code? For
that matter what does the slash mean?

--
View this message in context: 
http://r.789695.n4.nabble.com/ANOVA-1-too-few-degrees-of-freedom-tp3493349p3512977.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] Power Spectrum from STFT (e1071)?

2011-05-10 Thread ajf226
Hello.

Does anyone know how to generate a power spectrum from the STFT function in
package e1071? The output for this function supplies the Fourier
coefficients, but I do not know how to relate these back to the power
spectrum.

Alternatively, if anyone is aware of other packages/functions that perform
short time Fourier transforms, that would also be helpful.

Thanks.

--
View this message in context: 
http://r.789695.n4.nabble.com/Power-Spectrum-from-STFT-e1071-tp3512709p3512709.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] metaMDS and envfit: Help reading output

2011-05-10 Thread Songer, Katherine B - DNR
Hello R experts,

I've used metaMDS to run NMDS on some fish abundance data, and am also working 
on correlating environmental data to the NMDS coordinates. I'm fairly new to 
metaMDS and NMDS in general, so I have what are probably some very basic 
questions. My fish abundance data consists of 66 sites for which up to 20 
species of fish were identified and counted. I ran metaMDS on this data in 3 
dimensions (after using a scree plot to check for stress levels in the 
different dimensions). I then used envfit to correlate a predictor dataset of 
environmental variables with the NMDS results, using the following code.

>Fish<-as.data.frame(read.csv("Fish.csv",header=TRUE, sep = ","))
>Fish.mds<-metaMDS(Fish,zerodist = "add",k=3,trymax=20)

>Predictors<-as.data.frame(read.csv("Predictors.csv",header=TRUE, sep = ","))
>Fish.fit <- envfit(Fish.mds$points, Predictors, k=3, 1000, na.rm = TRUE)
>Fish.fit

The output of Fish.fit was as follows (table truncated):

Dim1Dim2r2  Pr(>r)
DrainArea   -0.5923233  -0.8057004  0.7674  0.000999 ***
Flow-0.5283236  -0.8490431  0.7847  0.000999 ***
StrmWidth   -0.6993457  -0.7147836  0.6759  0.000999 ***
Gradient0.4541225   0.8909392   0.2085  0.010989 *

I'd like to better understand how to read this table. I understand that Dim1 
and Dim2 refer to the dimensions of the vectors produced by envfit, and r2 is 
the r-squared of those vectors. But how do I visualize these vectors in a 3-d 
plot? To which of the 3 NMDS dimensions are these vectors being correlated? Is 
there code to produce the x, y, and z coordinates of each of the sites in 
Fish.mds?

Thanks very much.
Katie


[[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] specifying scales in lattice xyplot makes the lines disappear?

2011-05-10 Thread Peter Ehlers

On 2011-05-10 13:10, Christopher W Ryan wrote:

I have a dataframe concerning manner of death from death certificates,
from 2005 to 2009 inclusive, with the following structure:


str(MannerYoung.plot.data)

'data.frame':   245 obs. of  4 variables:
  $ year : Factor w/ 5 levels "2005","2006",..: 1 2 3 4 5 1 2 3 4 5 ...
  $ manner   : Factor w/ 7 levels "Accident","Homicide",..: 1 1 1 1
1 2 2 2 2 2 ...
  $ new.age.group: Factor w/ 7 levels "1 to 4","5 to 9",..: 1 1 1 1 1 1
1 1 1 1 ...
  $ Freq : int  [omitted from this post]...




new.age.group categorizes age of death in 5-year age bands.

The following makes a nice lattice scatterplot, with a panel for each
age-band, year on the horizontal axis,and frequency on the veritcal
axis.

xyplot(Freq~year | new.age.group, groups=manner,
data=MannerYoung.plot.data, type="l",  auto.key =
list(text=levels(MannerYoung.plot.data$manner), lines=TRUE,
points=FALSE, space = "top"),
layout=c(length(levels(MannerYoung.plot.data$new.age.group)),1))

Except that the year tick labels overlap each other. I thought I would
try to label only every other tick, using scales(), like this:

xyplot(Freq~year | new.age.group, groups=manner,
data=MannerYoung.plot.data, type="l",
scales=list(x=list(at=c(2005,2007,2009),
labels=c("2005","2007","2009"))), auto.key =
list(text=levels(MannerYoung.plot.data$manner), lines=TRUE,
points=FALSE, space = "top"),
layout=c(length(levels(MannerYoung.plot.data$new.age.group)),1))

This draws the panels, vertical axis, horizontal axes, with the
odd-numbered years, everything as desired, except that there are no
lines in the panels--no actual content.

I suppose I could decrease the font size for the year labels, using
cex(), and that cured the overlap, but I wanted to understand
scales().

Where am I going wrong?


You're using 'year' as a factor, but thinking about it as a numeric.
Use at=c(1,3,5) in your scales list.

But I would have expected the lines to plot and the axis labels
not to be visible. Are you telling the whole story?
_Minimal reproducible_ code is always a good thing.

Peter Ehlers



Thanks.

--Chris Ryan

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] ANOVA 1 too few degrees of freedom

2011-05-10 Thread Richard M. Heiberger
David,

Please read the first two paragraphs of Details from ?formula

You can try out the formulas with the notation

> expand.formula <- function(f) colnames(attr(terms(f), "factors"))
> expand.formula(~a+b)
[1] "a" "b"
> expand.formula(~a:b)
[1] "a:b"
> expand.formula(~a*b)
[1] "a"   "b"   "a:b"
> expand.formula(~a/b)
[1] "a"   "a:b"
>

See a design text for more information.

Rich


On Tue, May 10, 2011 at 4:51 PM, Rovinpiper wrote:

> So what is the difference between a colon and an asterisk in this code? For
> that matter what does the slash mean?
>
> --
> View this message in context:
> http://r.789695.n4.nabble.com/ANOVA-1-too-few-degrees-of-freedom-tp3493349p3512977.html
>  Sent from the R help mailing list archive at Nabble.com.
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[alternative HTML version deleted]]

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


[R] DCC-GARCH model and AR(1)-GARCH(1, 1) regression model - help needed..

2011-05-10 Thread Marcin P?�ciennik
Hello,
I have a rather complex problem... I will have to explain everything in
detail because I cannot solve it by myself...i just ran out of ideas. So
here is what I want to do:
I take quotes of two indices - S&P500 and DJ. And my first aim is to
estimate coefficients of the DCC-GARCH model for them. This is how I do it:


library(tseries)
p1 = get.hist.quote(instrument = "^gspc",start = "2005-01-07",end =
"2009-09-04",compression = "w", quote="AdjClose")
p2 = get.hist.quote(instrument = "^dji",start = "2005-01-07",end =
"2009-09-04",compression = "w", quote="AdjClose")
p = cbind(p1,p2)
y = diff(log(p))*100
y[,1] = y[,1]-mean(y[,1])
y[,2] = y[,2]-mean(y[,2])
T = length(y[,1])

library(ccgarch)
library(fGarch)

f1 = garchFit(~ garch(1,1), data=y[,1],include.mean=FALSE)
f1 = f1@fit$coef
f2 = garchFit(~ garch(1,1), data=y[,2],include.mean=FALSE)
f2 = f2@fit$coef

a = c(f1[1], f2[1])
A = diag(c(f1[2],f2[2]))
B = diag(c(f1[3], f2[3]))
dccpara = c(0.2,0.6)
dccresults = dcc.estimation(inia=a, iniA=A, iniB=B, ini.dcc=dccpara,dvar=y,
model="diagonal")

dccresults$out
DCCrho = dccresults$DCC[,2]
matplot(DCCrho, type='l')

dccresults$out deliver me the estimated coefficients of the DCC-GARCH model.
And here is my first question:
How can I check if these coefficients are significant or not? How can I test
them for significance?

second question would be:
Is this true that matplot(DCCrho, type='l') shows conditional correlation
between the two indices in question?


Ok. This would be it when it comes to DCC-GARCH.

Now, using conditional correlation obtained from the DCC-GARCH model, I want
to test for structural shifts in conditional correlations. To be precise, I
want to test whether the conditional correlations significantly increase in
the turmoil period / during the Subprime crisis.
The regression model is AR(1)-GARCH(1,1), using a dummy variable specified
as

where the first equation is the conditional correlation among the two
indices during the Subprime crisis, Dt is a dummy variable for the turmoil
period, and the second equation (hij,t) is the conditional variance of eij,t

The aim is, of course, to find the estimates of the regression model on
structural shifts in the conditional correlations obtained in the DCC-GARCH
model.

I found an information that there is no function for AR(1)-GARCH(1,1)
regression model. That's why it has to be done in two steps:
1) estimate the AR parameters
2) estimate the GARCH part of the model on the residuals from the AR model

And this would be my rather poor idea of how to do it...


library(timeSeries)
library(fSeries)
step1 = arma(DCCrho, order = c(1,0), include.intercept = TRUE)
step1$res
step11 = na.remove(step1$res)
step2 = garch (step11, order = c(1,1), include.intercept = TRUE)


To be honest I have no clue how to do it. I don't even now why do I get a
missing value as a result of step1 (step1$res[1]) and how to account for it?
Above, I just removed it but then I have a smaller number of
observations...and this is probably wrong.
And then these GARCH estimates on the residuals...does that make sense at
all?


I know the mail is quite looong, but hopefully, someone will find time to
give me a hand because I have to solve the problem and I reached the point
where I cannot move forward without someone's help. There is not much
information on how to apply DCC-GARCH model and AR(1)-GARCH(1,1) regression
model in the Internet. Hopefully, some of you are familiar with it.

Thank you very much in advance, people of good will, for looking at what I
wrote and helping me.

Best regards
Marcin

[[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] Converting ordinal dates and time into sensible formats

2011-05-10 Thread Steven Kennedy
How about something like:

> data$DateTime<- strptime(paste(data$day.hour,data$min.sec,sep=" 
> "),format="%j%H %M%S")
> data$Date <- strftime(data$DateTime,format="%Y-%m-%d")
> data$Time <- strftime(data$DateTime,format="%H:%M:%S")
> data
  day.hour min.secDateTime   Date Time
1125172034 2011-05-05 17:20:34 2011-05-05 17:20:34
2125202004 2011-05-05 20:20:04 2011-05-05 20:20:04
3125072045 2011-05-05 07:20:45 2011-05-05 07:20:45
4125022000 2011-05-05 02:20:00 2011-05-05 02:20:00
5125002004 2011-05-05 00:20:04 2011-05-05 00:20:04

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


Re: [R] Converting ordinal dates and time into sensible formats

2011-05-10 Thread Sam Albers
Yep that did the trick! I was sort stumbling along to this solution albeit
less efficiently. Thanks for your help!

Sam

On Tue, May 10, 2011 at 3:01 PM, Steven Kennedy  wrote:

> How about something like:
>
> > data$DateTime<- strptime(paste(data$day.hour,data$min.sec,sep="
> "),format="%j%H %M%S")
> > data$Date <- strftime(data$DateTime,format="%Y-%m-%d")
> > data$Time <- strftime(data$DateTime,format="%H:%M:%S")
> > data
>  day.hour min.secDateTime   Date Time
> 1125172034 2011-05-05 17:20:34 2011-05-05 17:20:34
> 2125202004 2011-05-05 20:20:04 2011-05-05 20:20:04
> 3125072045 2011-05-05 07:20:45 2011-05-05 07:20:45
> 4125022000 2011-05-05 02:20:00 2011-05-05 02:20:00
> 5125002004 2011-05-05 00:20:04 2011-05-05 00:20:04
>

[[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] metaMDS and envfit: Help reading output

2011-05-10 Thread Scott Chamberlain
Here is a thread I found on this 
issue:http://www.mail-archive.com/r-sig-ecology@r-project.org/msg00755.html 
It's not a simple correlation apparently, but I can't figure out how they get 
the output from envfit using lm(Y ~ scores1 + scores2). 

You can plot the data like:
library(MASS); require(vegan)
data(varespec) 
data(varechem) 
ord <- metaMDS(varespec) 
fit <- envfit(ord, varechem)
plot(ord)
plot(fit)


On Tuesday, May 10, 2011 at 3:51 PM, Songer, Katherine B - DNR wrote: 
> Hello R experts,
> 
> I've used metaMDS to run NMDS on some fish abundance data, and am also 
> working on correlating environmental data to the NMDS coordinates. I'm fairly 
> new to metaMDS and NMDS in general, so I have what are probably some very 
> basic questions. My fish abundance data consists of 66 sites for which up to 
> 20 species of fish were identified and counted. I ran metaMDS on this data in 
> 3 dimensions (after using a scree plot to check for stress levels in the 
> different dimensions). I then used envfit to correlate a predictor dataset of 
> environmental variables with the NMDS results, using the following code.
> 
> > Fish<-as.data.frame(read.csv("Fish.csv",header=TRUE, sep = ","))
> > Fish.mds<-metaMDS(Fish,zerodist = "add",k=3,trymax=20)
> 
> > Predictors<-as.data.frame(read.csv("Predictors.csv",header=TRUE, sep = ","))
> > Fish.fit <- envfit(Fish.mds$points, Predictors, k=3, 1000, na.rm = TRUE)
> > Fish.fit
> 
> The output of Fish.fit was as follows (table truncated):
> 
>  Dim1 Dim2 r2 Pr(>r)
> DrainArea -0.5923233 -0.8057004 0.7674 0.000999 ***
> Flow -0.5283236 -0.8490431 0.7847 0.000999 ***
> StrmWidth -0.6993457 -0.7147836 0.6759 0.000999 ***
> Gradient 0.4541225 0.8909392 0.2085 0.010989 *
> 
> I'd like to better understand how to read this table. I understand that Dim1 
> and Dim2 refer to the dimensions of the vectors produced by envfit, and r2 is 
> the r-squared of those vectors. But how do I visualize these vectors in a 3-d 
> plot? To which of the 3 NMDS dimensions are these vectors being correlated? 
> Is there code to produce the x, y, and z coordinates of each of the sites in 
> Fish.mds?
> 
> Thanks very much.
> Katie
> 
> 
>  [[alternative HTML version deleted]]
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 

[[alternative HTML version deleted]]

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


[R] Help documentation in extractAIC

2011-05-10 Thread Adler, Avraham
Hello.

The sentence in extractAIC's help 
 which 
discusses AIC's estimate of -2logL from RSS reads: "AIC only handles unknown 
scale and uses the formula n log (RSS/n) - n + n log 2pi - sum(log w) where w 
are the weights. Further AIC counts the scale estimation as a parameter in the 
edf and extractAIC does not." Am I misunderstanding something, or should it 
have a +n as the second term and not a -n? That's what seems to be happening in 
gaussian()$aic (not counting an extra addition of +2). It also could be that 
I've confused something somewhere, so I appreciate any elucidation.

Thank you,

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


[R] Problem rgl

2011-05-10 Thread tarnus
Dear List,

I'm using R on a Fedora 14 32 bit machine and I'm trying to use the
rgl-package, but get the following error-messages:

> library(rgl)
Warnmeldungen:
1: In rgl.init(initValue) : RGL: no suitable visual available
2: In fun(...) : error in rgl_init

> plot3d()
Fehler in rgl.cur() : rgl_dev_getcurrent

Does anyone know how to solve this problem?

Thanks in advance,

Marco

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


Re: [R] Problem rgl

2011-05-10 Thread Duncan Murdoch

On 11-05-10 6:23 PM, tarnus wrote:

Dear List,

I'm using R on a Fedora 14 32 bit machine and I'm trying to use the
rgl-package, but get the following error-messages:


library(rgl)

Warnmeldungen:
1: In rgl.init(initValue) : RGL: no suitable visual available
2: In fun(...) : error in rgl_init


plot3d()

Fehler in rgl.cur() : rgl_dev_getcurrent

Does anyone know how to solve this problem?


I don't know how to solve it, but I can tell you what the message means:

On your platform rgl needs to use X11, and it is apparently not able to 
get the type of "visual" it needs from X11 on your host.  Do other X11 
programs work?


Duncan Murdoch

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


Re: [R] Problem rgl

2011-05-10 Thread Prof Brian Ripley

On Tue, 10 May 2011, Duncan Murdoch wrote:


On 11-05-10 6:23 PM, tarnus wrote:

Dear List,

I'm using R on a Fedora 14 32 bit machine and I'm trying to use the
rgl-package, but get the following error-messages:


library(rgl)

Warnmeldungen:
1: In rgl.init(initValue) : RGL: no suitable visual available
2: In fun(...) : error in rgl_init


plot3d()

Fehler in rgl.cur() : rgl_dev_getcurrent

Does anyone know how to solve this problem?


I don't know how to solve it, but I can tell you what the message means:

On your platform rgl needs to use X11, and it is apparently not able to get 
the type of "visual" it needs from X11 on your host.  Do other X11 programs 
work?


I expect some do.  The requirements are in x11rgl.cpp, namely

GLX_RGBA,
GLX_DOUBLEBUFFER,
GLX_RED_SIZE, 1,
GLX_GREEN_SIZE, 1,
GLX_BLUE_SIZE, 1,
GLX_ALPHA_SIZE, 0,
GLX_DEPTH_SIZE, 1,

and 'man glXChooseVisual' explains them.  I suspect it is one of the 
first two which is the problem.


We are way off topic here (and Fedora 14 on our systems does provide 
suitable visuals), but that should be enough for you to get local 
Linux help on setting up your X server suitably.


--
Brian D. Ripley,  rip...@stats.ox.ac.uk
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] ggplot2 and add circle

2011-05-10 Thread Alaios
Thanks a lot this worked nice. 
It it possible also in ggplot2 to add a small figure? Let's say that I want to 
have somewhere in the plot the value of r printed.

How can I do that with ggplot 2?

Regards
Alex

--- On Tue, 5/10/11, Dennis Murphy  wrote:

> From: Dennis Murphy 
> Subject: Re: [R] ggplot2 and add circle
> To: "Alaios" 
> Cc: R-help@r-project.org
> Date: Tuesday, May 10, 2011, 7:12 PM
> Hi:
> 
> Here's one way:
> 
> plot_shad <- function(d, r) {
>    require('ggplot2')
>    plotdata <- melt(d)
>    names(plotdata)<-c('x','y','z')
>    xc <- mean(range(plotdata$x))
>    yc <- mean(range(plotdata$y))
>    theta <- seq(-pi, pi, length = 200)
>    circ <- data.frame(xv = xc + r *
> cos(theta),
>                
>       yv = yc + r * sin(theta))
>    v <- ggplot(plotdata)
>    print(v + geom_tile(aes(x = x, y = y,
> fill = z)) +
>          geom_path(data =
> circ, aes(x = xv, y = yv), color = 'white',
> size = 1) +
>          coord_equal()
>         )
>  }
> 
> plot_shad(f, 10)
> 
> HTH,
> Dennis
> 
> On Tue, May 10, 2011 at 10:15 AM, Alaios 
> wrote:
> > Here you are :)
> >
> >
> > plot_shad_f<-function(f){
> >   library(ggplot2)
> >   dev.new()
> >   plotdata<-melt(f)
> >   names(plotdata)<-c('x','y','z')
> >   v<-ggplot(plotdata, aes(x, y, z = z))
> >   print(v + geom_tile(aes(fill=z)))
> >
> > }
> >
> >
> > f<-matrix(data=seq(1:1),nrow=100,ncol=100)
> > plot_shad_f(f)
> >
> >
> > I would like to add a circle at the middle of this
> region with a range of 10.
> >
> > Best Regards
> >
> > Alex
> >
> > --- On Tue, 5/10/11, Scott Chamberlain 
> wrote:
> >
> > From: Scott Chamberlain 
> > Subject: Re: [R] ggplot2 and add circle
> > To: "Alaios" 
> > Cc: R-help@r-project.org
> > Date: Tuesday, May 10, 2011, 5:59 PM
> >
> >
> >
> >                You should provide reproducible
> data in addition to your code.
> > S
> >
> >
> >
> >
> >
> >                On Tuesday, May 10, 2011 at
> 11:54 AM, Alaios wrote:
> >
> >                    Dear all,
> > today I have writted the following code,
> > to plot the contents of some matrices I have
> >
> > plot_shad_f
> > function(f){
> >  library(ggplot2)
> >  dev.new()
> >  plotdata<-melt(f)
> >  names(plotdata)<-c('x','y','z')
> >  v<-ggplot(plotdata, aes(x, y, z = z))
> >  print(v + geom_tile(aes(fill=z)))
> > }
> >
> > I would like to ask your help add a small circle in
> this plotting. What would be the easiest way to do that in
> ggplot2?
> >
> > Best Regards
> > Alex
> >
> > __
> > R-help@r-project.org
> mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained,
> reproducible code.
> >
> >
> >
> >
> >
> >
> >
> >
> >
> >
> >
> >
> >
> >        [[alternative HTML version deleted]]
> >
> >
> > __
> > R-help@r-project.org
> mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained,
> reproducible code.
> >
> >
>

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