Re: [R] [Rcpp-devel] Help with integrating R and c/c++

2011-02-19 Thread Rohit Pandey
Hi Christopher/ Dirk,

Thank you very much for your replys. I think the idea of using inline as you
suggest is the best way to start off with using c++ with R. I went through
your examples and also plenty I found on the net. I have been trying to run
some of these over the past few days, but have consistently been getting
this error message as soon as I run the 'cfunction' or 'cxxfunction'. For
example, in the syntax that Dirk sent,
R> library(inline) #Runs perfectly
 R> src <- 'std::cout << "Hello C++_From_R World" << std::endl;
return(Rcpp::wrap(42));' #Runs perfectly
 R> rohit <- cxxfunction(signature(), src, plugin="Rcpp")
Now, as soon as I run this line, R spills a whole lot of text out and gives
an error and a warning:
ERROR(s) during compilation: source code errors or compiler configuration
errors!
Program source:
  1:
  2: // includes from the plugin
  3:
  4: #include 
  5:
  6:
  7: #ifndef BEGIN_RCPP
  8: #define BEGIN_RCPP
  9: #endif
 10:
 11: #ifndef END_RCPP
 12: #define END_RCPP
 13: #endif
 14:
 15: using namespace Rcpp;
 16:
 17:
 18: // user includes
 19:
 20:
 21: // declarations
 22: extern "C" {
 23: SEXP file59046688( ) ;
 24: }
 25:
 26: // definition
 27:
 28: SEXP file59046688(  ){
 29: BEGIN_RCPP
 30: std::cout << "Hello C++_From_R World" << std::endl;
return(Rcpp::wrap(42));
 31: END_RCPP
 32: }
 33:
 34:
Error in compileCode(f, code, language = language, verbose = verbose) :
  Compilation ERROR, function(s)/method(s) not created!
In addition: Warning message:
running command 'C:\PROGRA~1\R\R-212~1.1/bin/i386/R CMD SHLIB
file59046688.cpp 2> file59046688.cpp.err.txt' had status 1

The "file59046688.cpp 2"  changes every time I run a different function, but
the problem seems to be the same.

I installed and loaded the inline package (0.3.8) and then the Rcpp package
(0.9.0). I also tried reversing the order in which I load these, but still
no luck. I think if I can get just one of these programs to work, I will be
on my way. Can any of you tell me what I might be doing wrong?

For your question on what exacly I require, Christopher - I just need to use
a while loop. I have always been able to substitute for loops with some of
the apply functions in R, but can't seem to be able to replace the while
with a more efficient function. But the things that are required inside the
while loop, I have already implemented in R efficiently. So, I thought of
transfering just the while loop to a language that is faster with loops.

Thanks in advance.
On Mon, Feb 7, 2011 at 6:21 AM, Wray, Christopher <
christopher.wray...@ucl.ac.uk> wrote:

> As Dirk says, using "inline" makes it real simple to start and to prototype
> code.
>
> You mention you have R functions you wish to "call" via Rcpp. Im not
> certain I fully understand what you require here, but it is pretty simple to
> pass R-side functions to C++ via Rcpp, and similarly its simple to send
> compiled functions back to R-side as external pointers (and reuse them
> elsewhere). Here is a toy example using a simple "user" function defined on
> the R-side (user_F).
>
> This is evaluated in R, passed as a "function" parameter to compiled C++
> function and evaluated there, and then a compiled version of the function is
> passed back to R as an external pointer, which you can send back to the C
> side and evaluate:
>
>  user_F=function(v){sum((1-v)*exp(-0.5*v))}
>
> cpp <- '
> NumericVector numvec(xvec);
> NumericVector RetVec;
> Function userR_f(fun);
> List ATR;
> typedef SEXP (*h_ptr)(SEXP);
> RetVec = userR_f(numvec);
> ATR["Fn_ptr"]=XPtr (new h_ptr(&fme));
> ATR["Fn_VAL"]=RetVec;
> return ATR;
> '
> inc<-'
> using namespace Rcpp;
> SEXP fme(SEXP x){
> NumericVector xx(x);
> NumericVector yy=(1-xx)*exp(-0.5*xx);
> double Ans=sum(yy);
> return wrap(Ans);
> }
> '
>  FN_CSide <- cxxfunction(signature(xvec = "numeric", fun = "function"
> ),cpp, , plugin = "Rcpp", includes = inc)
>
> cppP <- '
> typedef SEXP (*h_ptr)(SEXP);
> XPtrp(x);
> NumericVector yy(y);
> yy=(*p)(yy);
> return yy;
> '
>  px<-cxxfunction(signature(x = "externalptr", y = "numeric" ),cppP, ,
> plugin = "Rcpp", includes = "using namespace Rcpp; ")
>
>  R_side=FN_CSide(seq(1,5),user_F)
>  user_F(seq(1,5))
> [1] -1.548486
>  R_side
> $Fn_ptr
> 
> $Fn_VAL
> [1] -1.548486
>  px(R_side$Fn_ptr,seq(1,5))
> [1] -1.548486
>
> I've broken it out to make the logic explicit. The above is sloppy - and
> makes no attempt to treat memory de/allocation/gc/protection issues (which
> you should consider) especially if you pass objects around, but it does try
> to show you how the mechanics work.
>
> Aside from looking up the Rcpp::sugar functionality used in places above,
> you might also want to look at package "RcppArmadillo", if functions and
> calculations involve linear algebra calcs/matrix operations. Hopefully this
> gives you some food-for-thought to delve further into the
> archives/documentation to find what you need.
> chris
>
>
>
> 
> From: rcpp-

Re: [R] censoring symbols on survfit plot

2011-02-19 Thread Patrick Connolly
On Fri, 18-Feb-2011 at 06:01AM -0600, Terry Therneau wrote:

|> --begin included message 
|> Hi, when ploting Kaplan-Meier estimate curves as below, the censoring
|> symbols
|> (crosses) to not change thickness along the lines 
|> plot(survfit(surv ~ I(x>=cut.off) ),lty=c(1,2), lwd=2)
|> 
|> is there any strightforward way to make it happen? thanks
|> 
|> -- End inclusion ---
|> 
|>  The symbols on the plot are characters, and are controled by the cex
|> parameter.  I don't know how to make a character "thicker" but not

The character will be thicker in a call to points() if the lwd
argument is set to something larger than 1.  I see in plot.survfit
there are quite a few calls to points().  Some or all of those calls
could have an lwd argument that was set somewhere in the call to
plot.survfit, but it would not be straightforward without that
modification.  

The OP could make a local version of plot.survfit with such a change
easily enough.  That's one of the great things about R: you can change
it if it's not quite what you need.  Whether it's a good idea to do
that to the package's own function is a decision for the maintainer.
There could well be occasions where it's not a good idea.

HTH

-- 
~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.   
   ___Patrick Connolly   
 {~._.~}   Great minds discuss ideas
 _( Y )_ Average minds discuss events 
(:_~*~_:)  Small minds discuss people  
 (_)-(_)  . Eleanor Roosevelt
  
~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.~.

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


[R] help about loop

2011-02-19 Thread 韩文衡


Hi all,
I am try to do a loop  on a data .such as:


X<-data.frame(name=c(1:9),SN=c(1,1,1,2,2,2,3,3,3),EW=rep(1:3,3))
X

  name SN EW
11  1  1
22  1  2
33  1  3
44  2  1
55  2  2
66  2  3
77  3  1
88  3  2
99  3  3



out<-list()
for ( i in 1:3) {out[[i]]<-subset( X,X$SN==1&EW==i) }  ##this expression  
can have the subset of EW=1 (and EW=2 EW=3) when SN=1

out

[[1]]
  name SN EW
11  1  1

[[2]]
  name SN EW
22  1  2

[[3]]
  name SN EW
33  1  3





Question:  What can I do if I want to have the all subset at once  when  
SN=1:3 ? How to write the loop?


thank you very much!

hanwenheng

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


[R] How to detect the spatial point pattern with Thomas process based on pcf

2011-02-19 Thread Jeff Fang
Hi all,

I am using "spatstat" to investigate the spatial structure of some plant
populations, but I have no idea about detecting the spatial point pattern
with Thomas process based on pcf. Additionally, generating simulation
envelope using this null model is another problem for me. I am not very
familiar with R, so I hope anyone have any experience with this work can
give me any assistance.

Many thanks


Jeff

[[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] help about loop

2011-02-19 Thread David Winsemius


On Feb 18, 2011, at 11:06 PM, 韩文衡 wrote:



Hi all,
I am try to do a loop  on a data .such as:


X<-data.frame(name=c(1:9),SN=c(1,1,1,2,2,2,3,3,3),EW=rep(1:3,3))
X

 name SN EW
11  1  1
22  1  2
33  1  3
44  2  1
55  2  2
66  2  3
77  3  1
88  3  2
99  3  3



out<-list()
for ( i in 1:3) {out[[i]]<-subset( X,X$SN==1&EW==i) }  ##this  
expression can have the subset of EW=1 (and EW=2 EW=3) when SN=1

out

[[1]]
 name SN EW
11  1  1

[[2]]
 name SN EW
22  1  2

[[3]]
 name SN EW
33  1  3





Question:  What can I do if I want to have the all subset at once   
when SN=1:3 ? How to write the loop?


?split



thank you very much!

hanwenheng

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

2011-02-19 Thread Dieter Menne


Elliot Joel Bernstein wrote:
> 
> I'm trying to use lattice graphics to produce some small plots for
> inclusion in a LaTeX file. I want the LaTeX fonts to be used in the plots,
> but to be scaled down to match the size of the plot. 
> 
> 

I agree with Deepayan. With tikzdevice, tust try to get the overall look
good, and use \resizebox or \scalebox in latex.

Dieter


-- 
View this message in context: 
http://r.789695.n4.nabble.com/Scaling-Lattice-Graphics-for-tikzDevice-tp3313399p3314045.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] clustered bar chart help

2011-02-19 Thread Jim Lemon
On 02/19/2011 06:46 AM, Li, Qinghong, ST. LOUIS, NRC St. Louis Petcare 
wrote:

Sorry the data were all messed up. Let me try again.
Can anyone help me to plot a chart graph? I have a data set like this. I
would like a bar chart graph which mouse1's two treatments are clustered
together, so on.

I tried with barplot, but couldn't get it right.

Mouse1  Mouse2  Mouse3  
Trial   C   T   C   T   C   T

1   2.4 4.3 2   4.3 2   7
2   2.8 5.4 2.8 6   2   8
3   3.1 4.5 3   4.5 3   7
4   3.1 4.4 3.1 4   3.1 7
5   2   4.1 2   4.1 2   7
6   2.2 4   2.2 3   4   6.9
7   2.7 5   4   5   2.7 6.5
8   4.5 5.6 4.5 6   4   5.6
9   2.3 4.5 2.3 4.5 2.3 4
10  3   4   4   6   2   4

-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org]
On Behalf Of Li,Qinghong,ST. LOUIS,NRC St. Louis Petcare
Sent: Friday, February 18, 2011 12:06 PM
To: r-help@r-project.org
Subject: [R] clustered bar chart help

Hi,



Can anyone help me to plot a chart graph? I have a data set like this. I
would like a bar chart graph which mouse1's two treatments are clustered
together, so on.

I tried with barplot, but couldn't get it right.


Hi Johnny,
Try this (assume your data frame is named "mp"):

# string the scores out in a vector
mpv<-as.vector(as.matrix(mp))
# make another two vectors specifying which score
# belongs to which mouse and which condition
# and combine these into a data frame
mp2<-data.frame(mpv,mouse=rep(c("M1","M2","M3"),each=20),
 cond=rep(c("C","T","C","T","C","T"),each=10))
library(plotrix)
# make a list of colors for the bars
mc_colors<-list("lightgray",2:4,c("orange","brown"))
# display a clustered bar chart
barNest(mpv~mouse+cond,mp2,col=mc_colors,main="Mouse performance")

Jim

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


Re: [R] breaks/bins question

2011-02-19 Thread Jim Lemon

On 02/19/2011 08:50 AM, Julie McWhorter wrote:


Thank you Joshua and Jim--I got the tick marks on finally!
One last question: I am supposed to make bins in 1, 2 and 5 mm increments for 
fish up to 170 mm.  I assume these are the breaks?  So breaks = 170 for 1 mm 
and 340 for 2 mm and 1550 for 5 mm?  I couldn't find any specific examples in 
the R manual.  The problem is I don't think that is correct, as the histogram 
disappears at breaks = 1550.


Hi Julie,
Not quite. I would guess:

breaks=0:170
breaks=seq(0,170,by=2)
breaks=seq(0,170,by=5)

I'm not surprised that something went wrong with 1550 breaks.

Jim

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


Re: [R] segfault during example(svm)

2011-02-19 Thread rose

Hi Peter,

Quoting Peter Ehlers :


On 2011-02-18 12:32, Juergen Rose wrote:

Am Freitag, den 18.02.2011, 11:53 -0800 schrieb Peter Ehlers:

On 2011-02-18 11:16, Juergen Rose wrote:

If do:

library("e1071")
example(svm)


I get:


svm>   data(iris)

svm>   attach(iris)

svm>   ## classification mode
svm>   # default with factor response:
svm>   model<- svm(Species ~ ., data = iris)

svm>   # alternatively the traditional interface:
svm>   x<- subset(iris, select = -Species)

svm>   y<- Species

svm>   model<- svm(x, y)

svm>   print(model)

Call:
svm.default(x = x, y = y)


Parameters:
SVM-Type:  C-classification
  SVM-Kernel:  radial
cost:  1
   gamma:  0.25

Number of Support Vectors:  51


svm>   summary(model)

Call:
svm.default(x = x, y = y)


Parameters:
SVM-Type:  C-classification
  SVM-Kernel:  radial
cost:  1
   gamma:  0.25

Number of Support Vectors:  51

  ( 8 22 21 )


Number of Classes:  3

Levels:
  setosa versicolor virginica




svm>   # test with train data
svm>   pred<- predict(model, x)

svm>   # (same as:)
svm>   pred<- fitted(model)

svm>   # Check accuracy:
svm>   table(pred, y)
 y
pred setosa versicolor virginica
   setosa 50  0 0
   versicolor  0 48 2
   virginica   0  248

svm>   # compute decision values and probabilities:
svm>   pred<- predict(model, x, decision.values = TRUE)

svm>   attr(pred, "decision.values")[1:4,]
   setosa/versicolor setosa/virginica versicolor/virginica
1  1.196152 1.0914600.6705626
2  1.064621 1.0563320.8479934
3  1.180842 1.0745340.6436474
4  1.110699 1.0531430.6778595

svm>   # visualize (classes by color, SV by crosses):
svm>   plot(cmdscale(dist(iris[,-5])),
svm+  col = as.integer(iris[,5]),
svm+  pch = c("o","+")[1:150 %in% model$index + 1])

  *** caught segfault ***
address (nil), cause 'unknown'

Traceback:
  1: .Call("La_rs", x, only.values, PACKAGE = "base")
  2: eigen(-x/2, symmetric = TRUE)
  3: cmdscale(dist(iris[, -5]))
  4: plot(cmdscale(dist(iris[, -5])), col = as.integer(iris[, 5]),
pch = c("o", "+")[1:150 %in% model$index + 1])
  5: eval.with.vis(expr, envir, enclos)
  6: eval.with.vis(ei, envir)
  7: source(tf, local, echo = echo, prompt.echo = paste(prompt.prefix,
getOption("prompt"), sep = ""), continue.echo = paste(prompt.prefix,
getOption("continue"), sep = ""), verbose = verbose, max.deparse.length
= Inf, encoding = "UTF-8", skip.echo = skips, keep.source = TRUE)
  8: example(svm)

Possible actions:
1: abort (with core dump, if enabled)
..

I did already "update.packages(), what can I still do.


Works just fine for me. What's your sessionInfo()?
Here's mine:
 >  sessionInfo()
R version 2.12.1 Patched (2010-12-27 r53883)
Platform: i386-pc-mingw32/i386 (32-bit)

locale:
[1] LC_COLLATE=English_Canada.1252  LC_CTYPE=English_Canada.1252
[3] LC_MONETARY=English_Canada.1252 LC_NUMERIC=C
[5] LC_TIME=English_Canada.1252

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

other attached packages:
[1] e1071_1.5-24 class_7.3-3

loaded via a namespace (and not attached):
[1] tools_2.12.1




sessionInfo()

R version 2.12.1 (2010-12-16)
Platform: x86_64-pc-linux-gnu (64-bit)

locale:
 [1] LC_CTYPE=en_US.UTF-8   LC_NUMERIC=C
 [3] LC_TIME=en_US.UTF-8LC_COLLATE=en_US.UTF-8
 [5] LC_MONETARY=C  LC_MESSAGES=en_US.UTF-8
 [7] LC_PAPER=en_US.UTF-8   LC_NAME=C
 [9] LC_ADDRESS=C   LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

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

It is working at some of my systems and is failing at the most.



It would be good to know what version of e1071 you're using;
presumably it's 1.5-24. I've run the example on both 32-bit
and 64-bit R 12.2.1pat and 2.13.0dev with no problems.
Could you be having a problem with your graphics device?

Peter Ehlers


I now have not acces to the computer, where the error happens first  
time On an other computer where the same issue occurs sessionInfo says:



sessionInfo()

R version 2.12.1 (2010-12-16)
Platform: x86_64-pc-linux-gnu (64-bit)
...
other attached packages:
[1] e1071_1.5-24 class_7.3-3

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

2011-02-19 Thread danielepippo

Hi,
I have two vector with the marginal distribution like this:
> a
[1] -0.419 -0.364 -0.159 -0.046 -0.010 -0.002  0.000  0.000  0.000
> b
[1] 0.125 0.260 0.270 0.187 0.097 0.041 0.014 0.004 0.001

How can I calculate the joint distribution with R?

Thank you to all
Dan
-- 
View this message in context: 
http://r.789695.n4.nabble.com/Joiny-probability-with-R-tp3314059p3314059.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] lmer, MCMCsamp and ranef samples?

2011-02-19 Thread Jukka Koskela

I really hope sombody could help me with the following,

I'm having problems accessing the random effect samples following the  
example on MCMCsamp:


(fm1 <- lmer(Reaction ~ Days + (1|Subject) + (0+Days|Subject), sleepstudy))

set.seed(101); samp0 <- mcmcsamp(fm1, n = 1000, saveb=TRUE)

str(samp0)

Formal class 'merMCMC' [package "lme4"] with 9 slots
  ..@ Gp  : int [1:3] 0 18 36
  ..@ ST  : num [1:2, 1:1000] 0.98 0.234 1.097 0.258 0.915 ...
  ..@ call: language lmer(formula = Reaction ~ Days + (1 |  
Subject) + (0 + Days |  Subject), data = sleepstudy)

  ..@ deviance: num [1:1000] 1752 1752 1753 1752 1752 ...
  ..@ dims: Named int [1:18] 2 180 2 36 1 2 0 1 2 5 ...
  .. ..- attr(*, "names")= chr [1:18] "nt" "n" "p" "q" ...
  ..@ fixef   : num [1:2, 1:1000] 251.41 10.47 254.56 9.71 262.67 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "(Intercept)" "Days"
  .. .. ..$ : NULL
  ..@ nc  : int [1:2] 1 1
  ..@ ranef   : num [1:36, 1:1000] 1.51 -40.37 -39.18 24.52 22.91 ...
  ..@ sigma   : num [1, 1:1000] 25.6 23.6 22.9 22.1 26.2 ...



I assume that there is random effect samples for the intercept  
(1|Subject) and the slope (0+Days|Subject).


The trouble is that I would like to get the samp0@ranef for the  
(0+Days|Subject) term, but I don't know how to separate the samples..


If I try for example:

samp0@ranef[1,],

I get n=1000 samples for what? I think it is for the intercept term..


Jukka

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] CHANGE OF SUMISSION ADDRESS, ASA John M. Chambers Statistical Software Award - 2011

2011-02-19 Thread Fei Chen

A reminder that submission deadline for John Chambers Award is Monday 2/21.

Also please note my email address has changed. Please submit your material to

fch...@its.jnj.com

Thank you.

From: f...@live.com
To: r-help@r-project.org
Subject: ASA John M. Chambers Statistical Software Award - 2011
Date: Tue, 14 Sep 2010 11:56:40 -0400








John M. Chambers Statistical Software Award - 2011
Statistical Computing Section
American Statistical Association

The Statistical Computing Section of the American Statistical
Association announces the competition for the John M.  Chambers
Statistical Software Award. In 1998 the Association for Computing
Machinery presented its Software System Award to John Chambers for the
design and development of S. Dr. Chambers generously donated his award
to the Statistical Computing Section to endow an annual prize for
statistical software written by an undergraduate or graduate student.
The prize carries with it a cash award of $1000, plus a substantial
allowance for travel to the annual Joint Statistical Meetings where
the award will be presented.

Teams of up to 3 people can participate in the competition, with the
cash award being split among team members. The travel allowance will
be given to just one individual in the team, who will be presented the
award at JSM.  To be eligible, the team must have designed and
implemented a piece of statistical software.  The individual within
the team indicated to receive the travel allowance must have begun the
development while a student, and must either currently be a student,
or have completed all requirements for her/his last degree after
January 1, 2009.  To apply for the award, teams must provide the
following materials:

   Current CV's of all team members.

   A letter from a faculty mentor at the academic institution of the
   individual indicated to receive the travel award.  The letter
   should confirm that the individual had substantial participation in
   the development of the software, certify her/his student status
   when the software began to be developed (and either the current
   student status or the date of degree completion), and briefly
   discuss the importance of the software to statistical practice.

   A brief, one to two page description of the software, summarizing
   what it does, how it does it, and why it is an important
   contribution.  If the team member competing for the travel
   allowance has continued developing the software after finishing
   her/his studies, the description should indicate what was developed
   when the individual was a student and what has been added since.

   An installable software package with its source code for use by the
   award committee. It should be accompanied by enough information to allow
   the judges to effectively use and evaluate the software (including
   its design considerations.)  This information can be provided in a
   variety of ways, including but not limited to a user manual (paper
   or electronic), a paper, a URL, and online help to the system.

All materials must be in English.  We prefer that electronic text be
submitted in Postscript or PDF.  The entries will be judged on a
variety of dimensions, including the importance and relevance for
statistical practice of the tasks performed by the software, ease of
use, clarity of description, elegance and availability for use by the
statistical community. Preference will be given to those entries that
are grounded in software design rather than calculation.  The decision
of the award committee is final.

All application materials must be received by 5:00pm EST, Monday,
February 21, 2011 at the address below.  The winner will be announced
in May and the award will be given at the 2011 Joint Statistical
Meetings.

Information on the competition can also be accessed on the website of
the Statistical Computing Section (www.statcomputing.org or see the
ASA website, www.amstat.org for a pointer), including the names and
contributions of previous winners.  Inquiries and application
materials should be emailed or mailed to:

Chambers Software Award
fch...@its.jnj.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] Hartley's table

2011-02-19 Thread Peter Ehlers

On 2011-02-16 06:08, Silvano wrote:

Hi,

I used the commands below to make Hartley's table,
but some values are NA.


You might want to contact the maintainer of SuppDists (cc'd).
Presumably, the C code uses inappropriate starting values
for those cases.

Meanwhile, here's a simpler way to generate your table:

 for(i in seq_len(ncol(har)))
   har[, i] <- qmaxFratio(.95, 2:(nrow(har) + 1), i+1)


Peter Ehlers



  require(SuppDists)
  trat = seq(2, 15, 1)
  gl = seq(2, 40, 1)

har = matrix(0, nr=length(gl), nc=length(trat))

for(i in 1:length(gl))
  for(j in 1:length(trat))
  har[i,j]<- qmaxFratio(.95, df=gl[i], k=trat[j])

rownames(har)<- gl
colnames(har)<- trat

head(har)


The output (head):


[... snip ...]

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

2011-02-19 Thread Ben Ward
I've just realised the couple of graphs I put on here have been stripped 
off. If anyone has to see them and can't see my problem from code, I can 
send them directly to anyone who thinks they can help but wants to see them.


Thanks,
Ben W.

On 18/02/2011 23:29, Ben Ward wrote:

Hi, I wonder if anyone could advise me with this:

I've been trying to make a standard curve in R with lm() of some 
standards from a spectrophotometer, so as I can express the curve as a 
formula, and so obtain values from my treated samples by plugging in 
readings into the formula, instead of trying to judge things by eye, 
with a curve drawn by hand.


It is a curve and so I used the following formula:

model <- lm(Approximate.Counts~X..Light.Transmission + 
I(Approximate.Counts^2), data=Standards)


It gives me a pretty decent graph:
xyplot(Approximate.Counts + fitted(model) ~ X..Light.Transmission, 
data=Standards)


I'm pretty happy with it, and looking at the model summary, to my 
inexperienced eyes it seems pretty good:


lm(formula = Approximate.Counts ~ X..Light.Transmission + 
I(Approximate.Counts^2),

data = Standards)

Residuals:
   Min 1Q Median 3QMax
-91.75 -51.04  27.33  37.28  49.72

Coefficients:
  Estimate Std. Error t value Pr(>|t|)
(Intercept)  9.868e+02  2.614e+01   37.75 <2e-16 ***
X..Light.Transmission   -1.539e+01  8.116e-01  -18.96 <2e-16 ***
I(Approximate.Counts^2)  2.580e-04  6.182e-06   41.73 <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 48.06 on 37 degrees of freedom
Multiple R-squared: 0.9956,Adjusted R-squared: 0.9954
F-statistic:  4190 on 2 and 37 DF,  p-value: < 2.2e-16

I tried to put some 95% confidence interval lines on a plot, as 
advised by my tutor, to see how they looked, and I used a function I 
found in "The R Book":


se.lines <- function(model){
b1<-coef(model)[2]+ summary(model)[[4]][4]
b2<-coef(model)[2]- summary(model)[[4]][4]
xm<-mean(model[[12]][2])
ym<-mean(model[[12]][1])
a1<-ym-b1*xm
a2<-ym-b2*xm
abline(a1,b1,lty=2)
abline(a2,b2,lty=2)
}
se.lines(model)

but when I do this on a plot I get an odd result:


They looks to me, to lie in the same kind of area, that my regression 
line did, before I used polynomial regression, by squaring 
"Approximate.Counts":


lm(formula = Approximate.Counts ~ X..Light.Transmission + 
I(Approximate.Counts^2), data = Standards)


Is there something else I should be doing? I've seen several ways of 
dealing with non-linear relationships, from log's of certain 
variables, and quadratic regression, and using sin and other 
mathematical devices. I'm not completely sure if I'm "allowed" to 
square the y variable, the book only squared the x variable in 
quadratic regression, which I did first, and it fit quite well, but 
not as good squaring Approximate Counts does:


model <- lm(Approximate.Counts~X..Light.Transmission + 
I(X..Light.Transmission^2), data=Standards)



Any advice is greatly appreciated, it's the first time I've really had 
to look at regression with data in my coursework that isn't a straight 
line.


Thanks,
Ben Ward.
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Building an array from matrix blocks

2011-02-19 Thread Eduardo de Oliveira Horta
Hello,

I've googled for a while and couldn't find anything on this topic: say
I have a matrix A and want to build matrices B1, B2,... using blocks
from A (or equivalently an array B with B[,,i] being a block from A),
and that I must sum the B[,,i]'s.

I've come up with this rather non-elegant code:

> n = 6
> p = 3
>
> A <- matrix(1:(n^2), n, n, byrow=TRUE)
>
> B <- array(0, c(n-p, n-p, p+1))
> for (i in 1:(p+1)) B[,,i] <- A[i:(n-p-1+i), i:(n-p-1+i)]
>
> X <- matrix(0, n-p, n-p)
> for (i in 1:(p+1)) X <- X + B[,,i]
> A
 [,1] [,2] [,3] [,4] [,5] [,6]
[1,]123456
[2,]789   10   11   12
[3,]   13   14   15   16   17   18
[4,]   19   20   21   22   23   24
[5,]   25   26   27   28   29   30
[6,]   31   32   33   34   35   36
> B
, , 1

 [,1] [,2] [,3]
[1,]123
[2,]789
[3,]   13   14   15

, , 2

 [,1] [,2] [,3]
[1,]89   10
[2,]   14   15   16
[3,]   20   21   22

, , 3

 [,1] [,2] [,3]
[1,]   15   16   17
[2,]   21   22   23
[3,]   27   28   29

, , 4

 [,1] [,2] [,3]
[1,]   22   23   24
[2,]   28   29   30
[3,]   34   35   36

> X
 [,1] [,2] [,3]
[1,]   46   50   54
[2,]   70   74   78
[3,]   94   98  102

Note that the blocks B[,,i] are obtained by sweeping the diagonal of
A. I wonder if there is a better and faster way to achieve this using
block matrix operations for instance. Actually what matters most for
me is getting to the matrix X, so if it is possible to do this without
having to construct the array B it would be ok as well...

Interesting observation:

> system.time(for (j in 1:1) {X <- matrix(0, n-p, n-p); for (i in 1:(p+1)) 
> X <- X + B[,,i]})
   user  system elapsed
   0.270.000.26
> system.time(for (j in 1:1) {X <- apply(B,c(1,2),sum)})
   user  system elapsed
   1.820.021.86

Thanks in advance, and best regards,

Eduardo Horta

> sessionInfo()
R version 2.11.1 (2010-05-31)
x86_64-pc-mingw32

locale:
[1] LC_COLLATE=Portuguese_Brazil.1252  LC_CTYPE=Portuguese_Brazil.1252
[3] LC_MONETARY=Portuguese_Brazil.1252 LC_NUMERIC=C
[5] LC_TIME=Portuguese_Brazil.1252

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

other attached packages:
[1] Revobase_4.2.0   RevoScaleR_1.1-1 lattice_0.19-13

loaded via a namespace (and not attached):
[1] grid_2.11.1   pkgXMLBuilder_1.0 revoIpe_1.0   tools_2.11.1
[5] XML_3.1-0

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

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 5:47 AM, danielepippo wrote:



Hi,
   I have two vector with the marginal distribution like this:

a

[1] -0.419 -0.364 -0.159 -0.046 -0.010 -0.002  0.000  0.000  0.000

b

[1] 0.125 0.260 0.270 0.187 0.097 0.041 0.014 0.004 0.001

How can I calculate the joint distribution with R?


Marginal distributions do not uniquely determine the joint  
distribution, Furthermore, the first one doesn't even look like a  
distribution or even a density. Densities are positive and CDFs range  
from 0 to 1. (The second on could be a discrete density for some set  
of values.) So you need to explain where those numbers come from and  
why you think we should apply sort of "distributional" assumptions  
about them.


--

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] When is *interactive* data visualization useful to use?

2011-02-19 Thread Liviu Andronic
On Fri, Feb 18, 2011 at 8:00 PM, Tom Hopper  wrote:
> Tal,
>
> One interactive capability that I have repeatedly wished for (but
> never taken the time to develop with the existing R tools) is the
> ability to interactively zoom in on and out of a data set,
>
I believe that you can do this with playwith. See this [1]. Regards
Liviu

[1] 
http://code.google.com/p/playwith/wiki/Screenshots#Time_series_plot_(Lattice)


and to
> interactively create "call-outs of sections of the data. Much of the
> data that I deal with takes the form of time series where both the
> full data and small section carry meaningful information.
>
> Some of the capabilities of Deducer approach interactive graphing,
> such as adjusting alpha values or smoothers, though the updates don't
> happen in quite real-time.
>
> - Tom
>
> On Friday, February 11, 2011, Tal Galili  wrote:
>> Hello all,
>>
>> Before getting to my question, I would like to apologize for asking this
>> question here.  My question is not directly an R question, however, I still
>> find the topic relevant to R community of users  - especially due to only *
>> partial* (current) support for interactive data visualization (see here:
>> http://cran.r-project.org/web/views/Graphics.html  were with iplots we are
>> waiting for iplots extreme, and with rggobi, it currently can not run with R
>> 2.12 and windows 7 OS).
>>
>> And now for my question:
>>
>> While preparing for a talk I will give soon, I recently started digging into
>> two major (Free) tools for interactive data visualization:
>> GGobi
>>  and mondrian  - both offer a great range of
>> capabilities (even if they're a bit buggy).
>>
>> I wish to ask for your help in articulating (both to myself, and for my
>> future audience) *When is it helpful to use interactive plots? Either for
>> data exploration (for ourselves) and data presentation (for a "client")?*
>>
>> For when explaining the data to a client, I can see the value of animation
>> for:
>>
>>    - Using "identify/linking/brushing" for seeing which data point in the
>>    graph is what.
>>    - Presenting a sensitivity analysis of the data (e.g: "if we remove this
>>    point, here is what we will get)
>>    - Showing the effect of different groups in the data (e.g: "let's look at
>>    our graphs for males and now for the females")
>>    - Showing the effect of time (or age, or in general, offering another
>>    dimension to the presentation)
>>
>> For when exploring the data ourselves, I can see the value of
>> identify/linking/brushing when exploring an outlier in a dataset we are
>> working on.
>>
>> But other then these two examples, I am not sure what other practical use
>> these techniques offer. Especially for our own data exploration!
>>
>> It could be argued that the interactive part is good for exploring (For
>> example) a different behavior of different groups/clusters in the data. But
>> when (in practice) I approached such situation, what I tended to do was to
>> run the relevant statistical procedures (and post-hoc tests) - and what I
>> found to be significant I would then plot with colors clearly dividing the
>> data to the relevant groups. From what I've seen, this is a safer approach
>> then "wondering around" the data (which could easily lead to data dredging
>> (were the scope of the multiple comparison needed for correction is not even
>> clear).
>>
>> I'd be very happy to read your experience/thoughts on this matter.
>>
>>
>> Thanks in advance,
>> Tal
>>
>>
>> Contact
>> Details:---
>> Contact me: tal.gal...@gmail.com |  972-52-7275845
>> Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) |
>> www.r-statistics.com (English)
>> --
>>
>>         [[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.
>



-- 
Do you know how to read?
http://www.alienetworks.com/srtest.cfm
http://goodies.xfce.org/projects/applications/xfce4-dict#speed-reader
Do you know how to write?
http://garbl.home.comcast.net/~garbl/stylemanual/e.htm#e-mail

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

[R] simple recoding problem, but a trouble !

2011-02-19 Thread Umesh Rosyara
Just a correction. My expected outdata frame was somehow distorted to a
single, one column. So correct one is:
 
marker1a markerb marker2amarker2b   
11   1   1  
13   1   3  
33   3   3  
33   3   3  
13   1   3  
13   1   3  
 
Thanks;
 
Umesh R 
 
  _  

From: Umesh Rosyara [mailto:rosyar...@gmail.com] 
Sent: Friday, February 18, 2011 10:09 PM
To: 'Joshua Wiley'
Cc: 'r-help@r-project.org'
Subject: RE: [R] recoding a data in different way: please help


Hi Josh and R community members 
 
Thank you for quick response. I am impressed with the help. 
 
To solve my problems, I tried recode options and I had the following problem
and which motivated me to leave it. Thank you for remind me the option
again, might help to solve my problem in different way. 
 
marker1 <- c("AA", "AC", "CC", "CC", "AC", "AC")

marker2 <- c("AA", "AC", "CC", "CC", "AC", "AC")

dfr <- data.frame(cbind(marker1, marker2))

Objective: replace A with 1, C with 3, and split AA into 1 1 (two columns
numeric). So the intended output for the above dataframe is:   



marker1a
markerb
marker2a
marker2b

1
1
1
1

1
3
1
3

3
3
3
3

3
3
3
3

1
3
1
3

1
3
1
3

I tried the following: 

 for(i in 1:length(dfr)) 
   {
 dfr[[i]]=recode (dfr[[i]],"c('AA')= '1,1'; c('AC')= '1,3'; c('CA')=
'1,3';  c('CC')= '3,3' ")
}

write.table(dfr,"dfr.out", sep=" ,", col.names = T) 
dfn=read.table("dfr.out",header=T, sep="," ) 

# just trying to cheat R, unfortunately the marker1 and marker columns
remained non-numeric, even when opened in excel !!


Unfortunately I got the following result ! 

   marker1 marker2
1 1,1  1,1
2 1,2  1,2
3 2,2  2,2
4 2,2  2,2
5 1,2  1,2
6 1,2  1,2

 
Sorry to bother all of you, but simple things are being complicated these
days to me. 
 
Thank you so much
Umesh R 

 
  _  

From: Joshua Wiley [mailto:jwiley.ps...@gmail.com] 
Sent: Friday, February 18, 2011 12:15 AM
Cc: r-help@r-project.org
Subject: Re: [R] recoding a data in different way: please help



Dear Umesh,

I could not figure out exactly what your recoding scheme was, so I do
not have a specific solution for you.  That said, the following
functions may help you get started.

?ifelse # vectorized and different from using if () statements
?if #
?Logic ## logical operators for your tests
## if you install and load the "car" package by John Fox
?recode # a function for recoding in package "car"

I am sure it is possible to string together some massive series of if
statements and then use a for loop, but that is probably the messiest
and slowest possible way.  I suspect there will be faster, neater
options, but I cannot say for certain without having a better feel for
how all the conditions work.

Best regards,

Josh

On Thu, Feb 17, 2011 at 6:21 PM, Umesh Rosyara  wrote:
> Dear R users
>
> The following question looks simple but I have spend alot of time to solve
> it. I would highly appeciate your help.
>
> I have following dataset from family dataset :
>
> Here we have individuals and their two parents and their marker scores
> (marker1, marker2,and so on). 0 means that their parent information
not
> available.
>
>
> Individual  Parent1  Parent2 mark1   mark2
> 10   0   12  11
> 20   0   11  22
> 30   0   13  22
> 40   0   13  11
> 51   2   11  12
> 61   2   12  12
> 73   4   11  12
> 83   4   13  12
> 91   4   11  12
> 10   1   4   11  12
>
> I want to recode mark1 and other mark2.and so on column by looking
> indvidual parent (Parent1 and Parent2).
>
> For example
>
> Take case of Individual 5, who's Parent 1 is 1 (has mark1 score 12) and
> Parent 2 is 2 (has mark1 score 11). Individual 5 has mark1 score 11.
Suppose
> I have following condition to recode Individual 5's mark1 score:
>
> For mark1 variable, If Parent1 score "11" and Parent2 score "22" and
recode
> indvidual 5's score, "12"=1, else 0
>If Parent1 score "12" and Parent2 score
> "22" and recode individual 5's score, "22"=1, "12"= 0.5, else 0
>.more
conditions
>
> Similarly the pointer should move from individual 5 to n individuals at
the
> end of the file.
>
>  Thank you in advance
>
> Umesh R
>
>
>
>
>
>[[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.
>



--
Joshua Wiley
Ph.D. Student, Health Psychology
University of California,

[R] contrasting Somer's D from Design package

2011-02-19 Thread vikkiyft

Dear R help,

I am having a problem with the Design package and my problem is detailed
here.

I fit a cox model to my data and validate the Somer's Dxy using the Design
package.
(Because of computation time problem, i only try 10 bootstrap samples for
the time being)

This is the model without stratification:
> library(Design)
> cox1a<-cph(surv.obj~factor(ecog)+factor(grade)+factor(tumor)+factor(extra),x=T,y=T)
> coef1a<-coef(cox1a)
> coef1a
ecog=1  ecog=2   grade=2   grade=3  tumor=2  tumor=3
extra=1 
  0.3578954   0.8993140   0.4834090   0.5716166   0.7600330   1.5974558  
0.8112942
> validate(cox1a,dxy=T,method="b",B=10)
 index.orig  training  test  optimism
index.corrected  n
Dxy   -5.450122e-01 -5.434214e-01 -5.437194e-01  2.979857e-04  
-5.453102e-01 10
R2 4.470271e-01  4.470172e-01  4.458043e-01  1.212920e-03   
4.458142e-01 10
Slope  1.00e+00  1.00e+00  9.910835e-01  8.916521e-03   
9.910835e-01 10
D  5.315241e-02  5.346057e-02  5.295422e-02  5.063561e-04   
5.264606e-02 10
U -4.651872e-05 -4.677365e-05  2.708456e-05 -7.385821e-05   
2.733949e-05 10
Q  5.319893e-02  5.350735e-02  5.292713e-02  5.802143e-04   
5.261872e-02 10

But if I fit a stratified cox model to the same data, the result becomes:
> cox1b<-cph(surv.obj~factor(ecog)+factor(grade)+factor(tumor)+factor(extra)+strat(ft),x=T,y=T,surv=T)
> coef1b<-coef(cox1b)
> coef1b
ecog=1  ecog=2   grade=2   grade=3  tumor=2  tumor=3
extra=1 
   0.1058620   0.5074692   0.3190904   0.4728515   0.5631035   1.2190005  
0.3500670 
> validate(cox1a,dxy=T,method="b",B=10,u=12) #for whatever u values I try,
> similar result is returned
 index.orig  training   test optimism index.corrected  n
Dxy6.004062e-01  5.469688e-01 0.54184724  0.005121538 0.595284689 10
R2 1.702137e-01  4.463073e-01 0.16343334  0.282873951-0.112660217 10
Slope  1.00e+00  1.00e+00 0.64516984  0.354830161 0.645169839 10
D  2.110517e-02  6.405669e-02 0.02018354  0.043873150-0.022767980 10
U -5.878953e-05 -5.625082e-05 0.00600501 -0.006061261 0.006002471 10
Q  2.116396e-02  6.411294e-02 0.01417853  0.049934410-0.028770451 10

The coefficients are different between the models, but the relative
magnitudes among the coefficients in each model are in fact quite stable. I
expect the results would be somewhat similar, yet the two models give
totally contrasting Dxy, one very negative and the other very positive. Also
the corrected R2 of the stratified model is negative, which is impossible.
Is there something wrong or is there a sensible interpretation for this?

Many thanks in advance.


Best regards,
Vikki
-- 
View this message in context: 
http://r.789695.n4.nabble.com/contrasting-Somer-s-D-from-Design-package-tp3314217p3314217.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] RCurl HTTP Post ?

2011-02-19 Thread Duncan Temple Lang


On 2/17/11 3:54 PM, Hasan Diwan wrote:
> According to [1] and [2], using RCurl to post a form with basic
> authentication is done using the postForm method. I'm trying to post
> generated interpolation data from R onto an HTTP form. The call I'm using is
> page <- postForm('http://our.server.com/dbInt/new', opts =
> curlOptions=(userpwd="test:test", verbose=T), profileid = "-1",
> value="1.801", type="history"). The page instance shows the HTTP response
> 500 screen and I get a nullpointerexception in the server logs. 

Do you mean that the R variable page gives information about the request error
and contains the 500 error code? Not sure what you mean by "screen" here.

Client-server interactions are hard to debug as the problems can be on either
side or in the communication.  The error can be in your request, in RCurl,
on the server side receiving the request or in the script processing the request
on the server.
So it is imperative to try to get diagnostic information.

You used verbose = T  (TRUE).  What did that display?

postForm() has a style parameter. It controls how the POST request is
submitted, either application/x-www-form-urlencoded or multipart/form-data.
Your server script might be expecting the data in a different format
than is being sent. postForm() defaults to the www-form-urlencoded.

But we will need more information to help you if these are not the
cause of the problem.

  D.

> The line it
> points to is dealing with getting an integer out of "profileid". Help?
> Many thanks in advance...

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] simple recoding problem, but a trouble !

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 8:40 AM, Umesh Rosyara wrote:

Just a correction. My expected outdata frame was somehow distorted  
to a

single, one column. So correct one is:

marker1a markerb marker2amarker2b   
11   1   1  
13   1   3  
33   3   3  
33   3   3  
13   1   3  
13   1   3  



func <- function(x) {sapply( strsplit(x, ""),
match, table= c("A", NA, "C"))}
t( apply(dfr, 1, func) )

 [,1] [,2] [,3] [,4]
[1,]1111
[2,]1313
[3,]3333
[4,]3333
[5,]1313
[6,]1313


It's amatrix rather than a dataframe and doesn't have colnames but  
that should be trivial to fix.




Thanks;

Umesh R

 _

From: Umesh Rosyara [mailto:rosyar...@gmail.com]
Sent: Friday, February 18, 2011 10:09 PM
To: 'Joshua Wiley'
Cc: 'r-help@r-project.org'
Subject: RE: [R] recoding a data in different way: please help


Hi Josh and R community members

Thank you for quick response. I am impressed with the help.

To solve my problems, I tried recode options and I had the following  
problem

and which motivated me to leave it. Thank you for remind me the option
again, might help to solve my problem in different way.

marker1 <- c("AA", "AC", "CC", "CC", "AC", "AC")

marker2 <- c("AA", "AC", "CC", "CC", "AC", "AC")

dfr <- data.frame(cbind(marker1, marker2))

Objective: replace A with 1, C with 3, and split AA into 1 1 (two  
columns

numeric). So the intended output for the above dataframe is:



marker1a
markerb
marker2a
marker2b

1
1
1
1

1
3
1
3

3
3
3
3

3
3
3
3

1
3
1
3

1
3
1
3

I tried the following:

for(i in 1:length(dfr))
  {
dfr[[i]]=recode (dfr[[i]],"c('AA')= '1,1'; c('AC')= '1,3';  
c('CA')=

'1,3';  c('CC')= '3,3' ")
}

write.table(dfr,"dfr.out", sep=" ,", col.names = T)
dfn=read.table("dfr.out",header=T, sep="," )

# just trying to cheat R, unfortunately the marker1 and marker columns
remained non-numeric, even when opened in excel !!


Unfortunately I got the following result !

  marker1 marker2
1 1,1  1,1
2 1,2  1,2
3 2,2  2,2
4 2,2  2,2
5 1,2  1,2
6 1,2  1,2


Sorry to bother all of you, but simple things are being complicated  
these

days to me.

Thank you so much
Umesh R


 _

From: Joshua Wiley [mailto:jwiley.ps...@gmail.com]
Sent: Friday, February 18, 2011 12:15 AM
Cc: r-help@r-project.org
Subject: Re: [R] recoding a data in different way: please help



Dear Umesh,

I could not figure out exactly what your recoding scheme was, so I do
not have a specific solution for you.  That said, the following
functions may help you get started.

?ifelse # vectorized and different from using if () statements
?if #
?Logic ## logical operators for your tests
## if you install and load the "car" package by John Fox
?recode # a function for recoding in package "car"

I am sure it is possible to string together some massive series of if
statements and then use a for loop, but that is probably the messiest
and slowest possible way.  I suspect there will be faster, neater
options, but I cannot say for certain without having a better feel for
how all the conditions work.

Best regards,

Josh

On Thu, Feb 17, 2011 at 6:21 PM, Umesh Rosyara   
wrote:

Dear R users

The following question looks simple but I have spend alot of time  
to solve

it. I would highly appeciate your help.

I have following dataset from family dataset :

Here we have individuals and their two parents and their marker  
scores
(marker1, marker2,and so on). 0 means that their parent  
information

not

available.


Individual  Parent1  Parent2 mark1   mark2
10   0   12  11
20   0   11  22
30   0   13  22
40   0   13  11
51   2   11  12
61   2   12  12
73   4   11  12
83   4   13  12
91   4   11  12
10   1   4   11  12

I want to recode mark1 and other mark2.and so on column by  
looking

indvidual parent (Parent1 and Parent2).

For example

Take case of Individual 5, who's Parent 1 is 1 (has mark1 score 12)  
and

Parent 2 is 2 (has mark1 score 11). Individual 5 has mark1 score 11.

Suppose

I have following condition to recode Individual 5's mark1 score:

For mark1 variable, If Parent1 score "11" and Parent2 score "22" and

recode

indvidual 5's score, "12"=1, else 0
  If Parent1 score "12" and Parent2  
score

"22" and recode individual 5's score, "22"=1, "12"= 0.5, else 0
  .more

conditions


Similarly the pointer should move from individual 5 to n  
individuals at

the

end of the file.

Thank you in advance

Umesh R





  [

[R] reading simulations

2011-02-19 Thread garciap

Hi to all the people (again),

I'm doing some simulations with the memisc package of an own function, but
I've the problem that I didn't know how to read the result of such
simulations. My function is:

> Torre<-function(a1,N1,a2,N2)
+ {Etorre<-(a1*N1)/(1+a1*N1)
+ Efuera<-(a2*N2)/(1+a2*N2)
+ if(Etorre>Efuera)Subir=TRUE
+ if (Etorre Torre<-Simulate(Torre(a1,N1,a2,N2),expand.grid(a1=3,N1=0.5,a2=c(0.01,0.02,0.05,0.1),N2=0.1),nsim=1000,seed=1,trace=50,keep.data=TRUE)

---
  a1  N1   a2  N2
1  3 0.5 0.01 0.1
Replication  50 
Replication  100 
Replication  150 
Replication  200 
Replication  250 
Replication  300 
Replication  350 
Replication  400 
Replication  450 
Replication  500 
Replication  550 
Replication  600 
Replication  650 
Replication  700 
Replication  750 
Replication  800 
Replication  850 
Replication  900 
Replication  950 
Replication  1000 

---
  a1  N1   a2  N2
2  3 0.5 0.02 0.1
Replication  50 
Replication  100 
Replication  150 
Replication  200 
Replication  250 
Replication  300 
Replication  350 
Replication  400 
Replication  450 
Replication  500 
Replication  550 
Replication  600 
Replication  650 
Replication  700 
Replication  750 
Replication  800 
Replication  850 
Replication  900 
Replication  950 
Replication  1000 

---
  a1  N1   a2  N2
3  3 0.5 0.05 0.1
Replication  50 
Replication  100 
Replication  150 
Replication  200 
Replication  250 
Replication  300 
Replication  350 
Replication  400 
Replication  450 
Replication  500 
Replication  550 
Replication  600 
Replication  650 
Replication  700 
Replication  750 
Replication  800 
Replication  850 
Replication  900 
Replication  950 
Replication  1000 

---
  a1  N1  a2  N2
4  3 0.5 0.1 0.1
Replication  50 
Replication  100 
Replication  150 
Replication  200 
Replication  250 
Replication  300 
Replication  350 
Replication  400 
Replication  450 
Replication  500 
Replication  550 
Replication  600 
Replication  650 
Replication  700 
Replication  750 
Replication  800 
Replication  850 
Replication  900 
Replication  950 
Replication  1000 
> print(Torre)
A 'default_bucket' object with 5 variables and 4000 observations

When I try print the results I obtain this, but how can I read it? It will
seem stupid to people, but I can't! even though function summary(Torre):

> summary(Torre)
Error in object[[i]] : wrong arguments for subsetting an environment
> 


Best regards (again) and many thanks for your help

Pablo
-- 
View this message in context: 
http://r.789695.n4.nabble.com/reading-simulations-tp3314281p3314281.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] [Rcpp-devel] Help with integrating R and c/c++

2011-02-19 Thread Douglas Bates
On Sat, Feb 19, 2011 at 12:56 AM, Rohit Pandey  wrote:
> Hi Christopher/ Dirk,

> Thank you very much for your replys. I think the idea of using inline as you
> suggest is the best way to start off with using c++ with R. I went through
> your examples and also plenty I found on the net. I have been trying to run
> some of these over the past few days, but have consistently been getting
> this error message as soon as I run the 'cfunction' or 'cxxfunction'. For
> example, in the syntax that Dirk sent,
> R> library(inline) #Runs perfectly
>  R> src <- 'std::cout << "Hello C++_From_R World" << std::endl;
> return(Rcpp::wrap(42));' #Runs perfectly
>  R> rohit <- cxxfunction(signature(), src, plugin="Rcpp")
> Now, as soon as I run this line, R spills a whole lot of text out and gives
> an error and a warning:
> ERROR(s) during compilation: source code errors or compiler configuration
> errors!
> Program source:
>  1:
>  2: // includes from the plugin
>  3:
>  4: #include 
>  5:
>  6:
>  7: #ifndef BEGIN_RCPP
>  8: #define BEGIN_RCPP
>  9: #endif
>  10:
>  11: #ifndef END_RCPP
>  12: #define END_RCPP
>  13: #endif
>  14:
>  15: using namespace Rcpp;
>  16:
>  17:
>  18: // user includes
>  19:
>  20:
>  21: // declarations
>  22: extern "C" {
>  23: SEXP file59046688( ) ;
>  24: }
>  25:
>  26: // definition
>  27:
>  28: SEXP file59046688(  ){
>  29: BEGIN_RCPP
>  30: std::cout << "Hello C++_From_R World" << std::endl;
> return(Rcpp::wrap(42));
>  31: END_RCPP
>  32: }
>  33:
>  34:
> Error in compileCode(f, code, language = language, verbose = verbose) :
>  Compilation ERROR, function(s)/method(s) not created!
> In addition: Warning message:
> running command 'C:\PROGRA~1\R\R-212~1.1/bin/i386/R CMD SHLIB
> file59046688.cpp 2> file59046688.cpp.err.txt' had status 1
>
> The "file59046688.cpp 2"  changes every time I run a different function, but
> the problem seems to be the same.
>
> I installed and loaded the inline package (0.3.8) and then the Rcpp package
> (0.9.0). I also tried reversing the order in which I load these, but still
> no luck. I think if I can get just one of these programs to work, I will be
> on my way. Can any of you tell me what I might be doing wrong?
>
> For your question on what exacly I require, Christopher - I just need to use
> a while loop. I have always been able to substitute for loops with some of
> the apply functions in R, but can't seem to be able to replace the while
> with a more efficient function. But the things that are required inside the
> while loop, I have already implemented in R efficiently. So, I thought of
> transfering just the while loop to a language that is faster with loops.

What operating system and version of R are you using?  It would help
if you included the results of executing

sessionInfo()

in R.  More importantly, do you have the compiler tools installed and
configured?  You need to have certain tools installed on Windows or
Mac OS X before you can compile packages (as opposed to installing
pre-compiled binary packages).  Most Linux distributions assume that
their users are adults and provide them with the tools to compile
programs.

>
> Thanks in advance.
> On Mon, Feb 7, 2011 at 6:21 AM, Wray, Christopher <
> christopher.wray...@ucl.ac.uk> wrote:
>
>> As Dirk says, using "inline" makes it real simple to start and to prototype
>> code.
>>
>> You mention you have R functions you wish to "call" via Rcpp. Im not
>> certain I fully understand what you require here, but it is pretty simple to
>> pass R-side functions to C++ via Rcpp, and similarly its simple to send
>> compiled functions back to R-side as external pointers (and reuse them
>> elsewhere). Here is a toy example using a simple "user" function defined on
>> the R-side (user_F).
>>
>> This is evaluated in R, passed as a "function" parameter to compiled C++
>> function and evaluated there, and then a compiled version of the function is
>> passed back to R as an external pointer, which you can send back to the C
>> side and evaluate:
>>
>>  user_F=function(v){sum((1-v)*exp(-0.5*v))}
>>
>> cpp <- '
>> NumericVector numvec(xvec);
>> NumericVector RetVec;
>> Function userR_f(fun);
>> List ATR;
>> typedef SEXP (*h_ptr)(SEXP);
>> RetVec = userR_f(numvec);
>> ATR["Fn_ptr"]=XPtr (new h_ptr(&fme));
>> ATR["Fn_VAL"]=RetVec;
>> return ATR;
>> '
>> inc<-'
>> using namespace Rcpp;
>> SEXP fme(SEXP x){
>> NumericVector xx(x);
>> NumericVector yy=(1-xx)*exp(-0.5*xx);
>> double Ans=sum(yy);
>> return wrap(Ans);
>> }
>> '
>>  FN_CSide <- cxxfunction(signature(xvec = "numeric", fun = "function"
>> ),cpp, , plugin = "Rcpp", includes = inc)
>>
>> cppP <- '
>> typedef SEXP (*h_ptr)(SEXP);
>> XPtrp(x);
>> NumericVector yy(y);
>> yy=(*p)(yy);
>> return yy;
>> '
>>  px<-cxxfunction(signature(x = "externalptr", y = "numeric" ),cppP, ,
>> plugin = "Rcpp", includes = "using namespace Rcpp; ")
>>
>>  R_side=FN_CSide(seq(1,5),user_F)
>>  user_F(seq(1,5))
>> [1] -1.548486
>>  R_side
>> $Fn_ptr
>> 
>> $Fn_VAL
>> [1] -1.5484

Re: [R] reading simulations

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 10:35 AM, garciap wrote:



Hi to all the people (again),

I'm doing some simulations with the memisc package of an own  
function, but

I've the problem that I didn't know how to read the result of such
simulations. My function is:


Torre<-function(a1,N1,a2,N2)

+ {Etorre<-(a1*N1)/(1+a1*N1)
+ Efuera<-(a2*N2)/(1+a2*N2)
+ if(Etorre>Efuera)Subir=TRUE
+ if (EtorreTorre<- 
Simulate 
(Torre 
(a1 
,N1 
,a2 
,N2 
),expand 
.grid 
(a1 
= 
3 
,N1 
= 
0.5 
,a2 
= 
c 
(0.01,0.02,0.05,0.1 
),N2=0.1),nsim=1000,seed=1,trace=50,keep.data=TRUE)




?Simulate

Try something like this, since I thought your example was far from  
"minimal". After reading that the value was a dataframe, I just used  
standard indexing to access individual rows. I used instead the first  
example from the help page:

require(memisc)
Normal.example <- function(mean=0,sd=1,n=10){
 x <- rnorm(n=n,mean=mean,sd=sd)
 c( Mean=mean(x), Median=median(x), Var=var(x) ) }
Normal.simres <- Simulate(
 Normal.example(mean,sd,n),
 expand.grid( mean=0, sd=c(1,10), n=c(10,100) ), nsim=200, trace=50)
> Normal.simres[1, ]
  mean sd  n   Mean Median  Var
10  1 10 0.01784074 -0.3662563 1.381853
> Normal.simres[200, ]
mean sd  n  MeanMedian  Var
2000  1 10 0.4004172 0.3280587 1.207104

That object does pose some surprises:
> str(Normal.simres[, 1])  # regular extraction doesn't return what I  
expect

'data.frame':   800 obs. of  6 variables:
 $ mean  : num  0 0 0 0 0 0 0 0 0 0 ...
 $ sd: num  1 1 1 1 1 1 1 1 1 1 ...
 $ n : num  10 10 10 10 10 10 10 10 10 10 ...
 $ Mean  : num  0.0178 -0.0803 0.4951 0.2081 -0.3496 ...
 $ Median: num  -0.36626 -0.29166 0.27018 -0.00887 -0.5638 ...
 $ Var   : num  1.382 0.798 1.18 0.48 0.571 ...
> is.data.frame(Normal.simres)
[1] FALSE
> is.data.frame(Normal.simres[1])
[1] TRUE
> str(Normal.simres[, 1][2])
'data.frame':   800 obs. of  1 variable:
 $ sd: num  1 1 1 1 1 1 1 1 1 1 ...
> str(Normal.simres[, 1]["Mean"])  # but it is possible to extract  
the "good stuff"

'data.frame':   800 obs. of  1 variable:
 $ Mean: num  0.0178 -0.0803 0.4951 0.2081 -0.3496 ...

So it appears to be something not quite a dataframe.

> quantile(Normal.simres[]$Mean , c(.25,.5,.75))  # works
25% 50% 75%
-0.22652591  0.02842814  0.42513166

> str(Normal.simres$Mean ) # doesn't work
 NULL


---
 a1  N1   a2  N2
1  3 0.5 0.01 0.1
Replication  50
Replication  100
Replication  150
Replication  200
Replication  250
Replication  300
Replication  350
Replication  400
Replication  450
Replication  500
Replication  550
Replication  600
Replication  650
Replication  700
Replication  750
Replication  800
Replication  850
Replication  900
Replication  950
Replication  1000

---
 a1  N1   a2  N2
2  3 0.5 0.02 0.1
Replication  50
Replication  100
Replication  150
Replication  200
Replication  250
Replication  300
Replication  350
Replication  400
Replication  450
Replication  500
Replication  550
Replication  600
Replication  650
Replication  700
Replication  750
Replication  800
Replication  850
Replication  900
Replication  950
Replication  1000

---
 a1  N1   a2  N2
3  3 0.5 0.05 0.1
Replication  50
Replication  100
Replication  150
Replication  200
Replication  250
Replication  300
Replication  350
Replication  400
Replication  450
Replication  500
Replication  550
Replication  600
Replication  650
Replication  700
Replication  750
Replication  800
Replication  850
Replication  900
Replication  950
Replication  1000

---
 a1  N1  a2  N2
4  3 0.5 0.1 0.1
Replication  50
Replication  100
Replication  150
Replication  200
Replication  250
Replication  300
Replication  350
Replication  400
Replication  450
Replication  500
Replication  550
Replication  600
Replication  650
Replication  700
Replication  750
Replication  800
Replication  850
Replication  900
Replication  950
Replication  1000

print(Torre)

A 'default_bucket' object with 5 variables and 4000 observations

When I try print the results I obtain this, but how can I read it?  
It will
seem stupid to people, but I can't! even though function  
summary(Torre):



summary(Torre)

Error in object[[i]] : wrong arguments for subsetting an environment





Best regards (again) and many thanks for your help

Pablo
--
View this message in context: 
http://r.789695.n4.nabble.com/reading-simulations-tp3314281p3314281.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.


David Winsemius, MD
West Hartford, CT

__
R-help@

Re: [R] Joiny probability with R

2011-02-19 Thread danielepippo

I'm sorry it was my mistake. The two vectors are found as:
dpois(num,0.87)
dpois(num,2.08)
and represent the discrete density  of 8 intervals (num=0:8).


-- 
View this message in context: 
http://r.789695.n4.nabble.com/Joiny-probability-with-R-tp3314059p3314301.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] Joiny probability with R

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 11:05 AM, danielepippo wrote:



I'm sorry it was my mistake. The two vectors are found as:
dpois(num,0.87)
dpois(num,2.08)
and represent the discrete density  of 8 intervals (num=0:8).



Is this homework?



--
View this message in context: 
http://r.789695.n4.nabble.com/Joiny-probability-with-R-tp3314059p3314301.html
Sent from the R help mailing list archive at Nabble.com.


--

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] little help with Nonlinear regression needed

2011-02-19 Thread Paka yag

  Hello

as the subject says I need a little help with following nonlinear regression

y(t) = a + b*x(t-1)+c*y_(t-1)+G*(a+b*x_(t-1)) 

where G=[1 + exp(-s(x_(t-1) - k)]^(-1)

(In reality there are more variables)

Please could anybody give me a hint how I can estimate this?? 
Should I use nls()? what parameter settings would I use then?

And one more: how would I estimate it if the independant variables werent 
lagged, would it make a difference?

Your help is much appreciated!!



Thank you

___
Schon gehört? WEB.DE hat einen genialen Phishing-Filter in die

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

2011-02-19 Thread Ted Harding
On 19-Feb-11 14:48:53, David Winsemius wrote:
> On Feb 19, 2011, at 5:47 AM, danielepippo wrote:
>> Hi,
>>I have two vector with the marginal distribution like this:
>>> a
>> [1] -0.419 -0.364 -0.159 -0.046 -0.010 -0.002  0.000  0.000  0.000
>>> b
>> [1] 0.125 0.260 0.270 0.187 0.097 0.041 0.014 0.004 0.001
>>
>> How can I calculate the joint distribution with R?
> 
> Marginal distributions do not uniquely determine the joint  
> distribution, Furthermore, the first one doesn't even look like a  
> distribution or even a density. Densities are positive and CDFs range  
> from 0 to 1. (The second on could be a discrete density for some set  
> of values.) So you need to explain where those numbers come from and  
> why you think we should apply sort of "distributional" assumptions  
> about them.
> --
> David Winsemius, MD
> West Hartford, CT

Following up on David's reply (I was in the midst of thinking about
this when he responded):

I strongly suspect that the values of 'a' have had their signs
changed. If we simply make them all positive:

  a <- c(0.419,0.364,0.159,0.046,0.010,0.002,0.000,0.000,0.000)

then now sum(a) = 1 and these can be probabilities. Also, with
the values given, sum(b) = 0.999; so (in the relatively least
disturbing way) I shall increase the largest value by 0.001:

  b <- c(0.125,0.260,0.271,0.187,0.097,0.041,0.014,0.004,0.001)

(0.270 -> 0.271). Now also sum(b) = 1.

Now, as a somewhat trivial example to illustrate the non-uniqueness
that David points out, I shall invent a very arbitrary set of
probabilities P = {P[i,j]} such that they have the given marginal
distributions.

Since we have the marginals given to 3 decimal places, now consider
allocating 1000 pairs of random numbers (x,y). Then we pick out the
smallest 419 of the x's to go into the bottom class of 'a'
(proportion = 419/1000 = 0.419), then the 364 next smallest for the
second class of 'a', and so on. And similarly for the numbers 'y'
to go into the classes of 'b'.

Then the inequalities which define these marginal divisions can be
applied jointly to produce the joint divisions. In the first instance
the results will be computed as counts

Therefore:

  X0 <- runif(1000) ; X<- sort(X0)
  Y0 <- runif(1000) ; Y<- sort(Y0)

  A <- cumsum(c(419,364,159, 46, 10,  2,  0,  0,  0))
  B <- cumsum(c(125,260,271,187, 97, 41, 14,  4,  1))

  Xdiv <- numeric(8)
  Ydiv <- numeric(8)

  for(i in (1:8)){
Xdiv[i] <- 0.5*(max(X[1:A[i]])+min(X[(A[i]+1):A[i+1]]))
Ydiv[i] <- 0.5*(max(Y[1:B[i]])+min(Y[(B[i]+1):B[i+1]]))
Xdiv[is.na(Xdiv)] <- 1.0
Ydiv[is.na(Ydiv)] <- 1.0
  }
  Xdiv <-c(0,Xdiv,1)
  Ydiv <- c(0,Ydiv,1)

  N <- matrix(0,nrow=9,ncol=9)

  for(i in (1:9)){
for(j in (1:9)){
  N[i,j] <- sum((Xdiv[i]
Fax-to-email: +44 (0)870 094 0861
Date: 19-Feb-11   Time: 16:30:08
-- XFMail --

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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 in plotting numeric x by POSIXt class with lattice

2011-02-19 Thread Aviad Klein
# hi all,

# I'm trying to plot temperatures by date in a trellis plot by their
stations

# I'm plotting the following data.frame

library(lattice)

h <- structure(list(station_name = structure(c(3L, 4L, 2L, 10L, 11L,
12L, 6L, 7L, 5L, 8L, 9L, 3L, 4L, 2L, 10L, 11L, 12L, 6L, 7L, 5L,
8L, 9L, 3L, 4L, 2L, 10L, 11L, 12L, 6L, 7L), .Label = c("Ashqelon",
"Beer Sheva", "Beit Dagan", "Eden", "Ein Carmel", "Jerusalem",
"Negba", "Nir HaEmek", "Ramat HaNegev", "Ramon", "Tel Mond",
"Zemah"), class = "factor"), max_dry = c(23.4, 24.3, 25.5, 21.8,
23.1, 22, 23, 24, 22.7, 24.9, 23.3, 20.2, 20.8, 22.5, 20, 20.6,
20.5, 20.9, 21, 19.2, 21.8, 20.7, 21, 21.8, 20.2, 18.4, 20.5,
20.3, 17.1, 21.4), measure_date = structure(list(sec = c(0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0), min = c(0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L), hour = c(0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L), mday = c(1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), mon = c(0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L), year = c(100, 100, 100,
100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100,
100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100,
100), wday = c(3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L,
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L,
5L, 5L), yday = c(0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L), isdst = c(0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L)), .Names = c("sec", "min", "hour", "mday", "mon",
"year", "wday", "yday", "isdst"), class = c("POSIXlt", "POSIXt"
))), .Names = c("station_name", "max_dry", "measure_date"), row.names =
c(NA,
30L), class = "data.frame")

# using the following :

xyplot(max_dry~measure_date|station_name,data=h)

# I receive the following error
# Error in structure(xx, class = c("POSIXct", "POSIXt"), tzone = tz) :
invalid 'x' argument

# does anyone know why?

# thanks

[[alternative HTML version deleted]]

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


Re: [R] problem in plotting numeric x by POSIXt class with lattice

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 11:44 AM, Aviad Klein wrote:


# hi all,

# I'm trying to plot temperatures by date in a trellis plot by their
stations

# I'm plotting the following data.frame

library(lattice)

h <- structure(list(station_name = structure(c(3L, 4L, 2L, 10L, 11L,
12L, 6L, 7L, 5L, 8L, 9L, 3L, 4L, 2L, 10L, 11L, 12L, 6L, 7L, 5L,
8L, 9L, 3L, 4L, 2L, 10L, 11L, 12L, 6L, 7L), .Label = c("Ashqelon",
"Beer Sheva", "Beit Dagan", "Eden", "Ein Carmel", "Jerusalem",
"Negba", "Nir HaEmek", "Ramat HaNegev", "Ramon", "Tel Mond",
"Zemah"), class = "factor"), max_dry = c(23.4, 24.3, 25.5, 21.8,
23.1, 22, 23, 24, 22.7, 24.9, 23.3, 20.2, 20.8, 22.5, 20, 20.6,
20.5, 20.9, 21, 19.2, 21.8, 20.7, 21, 21.8, 20.2, 18.4, 20.5,
20.3, 17.1, 21.4), measure_date = structure(list(sec = c(0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0), min = c(0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L), hour = c(0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L), mday = c(1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), mon = c(0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L), year = c(100, 100, 100,
100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100,
100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100,
100), wday = c(3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L,
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L,
5L, 5L), yday = c(0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L), isdst = c(0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, 0L, 0L)), .Names = c("sec", "min", "hour", "mday", "mon",
"year", "wday", "yday", "isdst"), class = c("POSIXlt", "POSIXt"
))), .Names = c("station_name", "max_dry", "measure_date"),  
row.names =

c(NA,
30L), class = "data.frame")

# using the following :

xyplot(max_dry~measure_date|station_name,data=h)

# I receive the following error
# Error in structure(xx, class = c("POSIXct", "POSIXt"), tzone = tz) :
invalid 'x' argument

# does anyone know why?


POSIXlt dates are not always handled well by plotting routines. Try:

> h$measure_date <- as.POSIXct(h$measure_date)
> xyplot(max_dry ~ measure_date | station_name, data=h)

--
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] Joiny probability with R

2011-02-19 Thread Gabor Grothendieck
On Sat, Feb 19, 2011 at 11:30 AM, Ted Harding  wrote:
> On 19-Feb-11 14:48:53, David Winsemius wrote:
>> On Feb 19, 2011, at 5:47 AM, danielepippo wrote:
>>> Hi,
>>>    I have two vector with the marginal distribution like this:
 a
>>> [1] -0.419 -0.364 -0.159 -0.046 -0.010 -0.002  0.000  0.000  0.000
 b
>>> [1] 0.125 0.260 0.270 0.187 0.097 0.041 0.014 0.004 0.001
>>>
>>> How can I calculate the joint distribution with R?
>>
>> Marginal distributions do not uniquely determine the joint
>> distribution, Furthermore, the first one doesn't even look like a
>> distribution or even a density. Densities are positive and CDFs range
>> from 0 to 1. (The second on could be a discrete density for some set
>> of values.) So you need to explain where those numbers come from and
>> why you think we should apply sort of "distributional" assumptions
>> about them.
>> --
>> David Winsemius, MD
>> West Hartford, CT
>
> Following up on David's reply (I was in the midst of thinking about
> this when he responded):
>
> I strongly suspect that the values of 'a' have had their signs
> changed. If we simply make them all positive:
>
>  a <- c(0.419,0.364,0.159,0.046,0.010,0.002,0.000,0.000,0.000)
>
> then now sum(a) = 1 and these can be probabilities. Also, with
> the values given, sum(b) = 0.999; so (in the relatively least
> disturbing way) I shall increase the largest value by 0.001:
>
>  b <- c(0.125,0.260,0.271,0.187,0.097,0.041,0.014,0.004,0.001)
>
> (0.270 -> 0.271). Now also sum(b) = 1.
>
> Now, as a somewhat trivial example to illustrate the non-uniqueness
> that David points out, I shall invent a very arbitrary set of
> probabilities P = {P[i,j]} such that they have the given marginal
> distributions.
>
> Since we have the marginals given to 3 decimal places, now consider
> allocating 1000 pairs of random numbers (x,y). Then we pick out the
> smallest 419 of the x's to go into the bottom class of 'a'
> (proportion = 419/1000 = 0.419), then the 364 next smallest for the
> second class of 'a', and so on. And similarly for the numbers 'y'
> to go into the classes of 'b'.
>
> Then the inequalities which define these marginal divisions can be
> applied jointly to produce the joint divisions. In the first instance
> the results will be computed as counts
>
> Therefore:
>
>  X0 <- runif(1000) ; X<- sort(X0)
>  Y0 <- runif(1000) ; Y<- sort(Y0)
>
>  A <- cumsum(c(419,364,159, 46, 10,  2,  0,  0,  0))
>  B <- cumsum(c(125,260,271,187, 97, 41, 14,  4,  1))
>
>  Xdiv <- numeric(8)
>  Ydiv <- numeric(8)
>
>  for(i in (1:8)){
>    Xdiv[i] <- 0.5*(max(X[1:A[i]])+min(X[(A[i]+1):A[i+1]]))
>    Ydiv[i] <- 0.5*(max(Y[1:B[i]])+min(Y[(B[i]+1):B[i+1]]))
>    Xdiv[is.na(Xdiv)] <- 1.0
>    Ydiv[is.na(Ydiv)] <- 1.0
>  }
>  Xdiv <-c(0,Xdiv,1)
>  Ydiv <- c(0,Ydiv,1)
>
>  N <- matrix(0,nrow=9,ncol=9)
>
>  for(i in (1:9)){
>    for(j in (1:9)){
>      N[i,j] <- sum((Xdiv[i]                    (Ydiv[j]    }
>  }
>
> N
> rowSums(N)
> colSums(N)
>
> Here is one typical result (depending on the random numbers) for N:
>
>  N
>  #       [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
>  #  [1,]   54  113  112   76   36   20    7    1    0
>  #  [2,]   44   90  101   75   35   13    3    2    1
>  #  [3,]   20   40   42   28   21    6    2    0    0
>  #  [4,]    5   13   13    6    5    2    1    1    0
>  #  [5,]    2    3    3    1    0    0    1    0    0
>  #  [6,]    0    1    0    1    0    0    0    0    0
>  #  [7,]    0    0    0    0    0    0    0    0    0
>  #  [8,]    0    0    0    0    0    0    0    0    0
>  #  [9,]    0    0    0    0    0    0    0    0    0
>  rowSums(N)
>  # [1] 419 364 159  46  10   2   0   0   0 ## Agreeing with (A)
>  colSums(N)
>  # [1] 125 260 271 187  97  41  14   4   1 ## Agreeing with (B)
>
> And here is another:
>
>  #       [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
>  #  [1,]   52  109  113   89   35   15    3    3    0
>  #  [2,]   45   95  108   50   40   16    9    0    1
>  #  [3,]   22   44   33   33   16    9    1    1    0
>  #  [4,]    4    9   14   13    4    1    1    0    0
>  #  [5,]    2    3    3    1    1    0    0    0    0
>  #  [6,]    0    0    0    1    1    0    0    0    0
>  #  [7,]    0    0    0    0    0    0    0    0    0
>  #  [8,]    0    0    0    0    0    0    0    0    0
>  #  [9,]    0    0    0    0    0    0    0    0    0
>  rowSums(N)
>  # [1] 419 364 159  46  10   2   0   0   0
>  colSums(N)
>  # [1] 125 260 271 187  97  41  14   4   1
>
> Then N/1000 gives the joint probability distribution, and
> rowSums(N) and colSums(N) give the marginal distributions.

Here is another way to demonstrate two different joint probability
distributions with the same marginals using R:

set.seed(123)

a <- c(-0.419, -0.364, -0.159, -0.046, -0.01, -0.002, 0, 0, 0)
b <- c(0.125, 0.26, 0.27, 0.187, 0.097, 0.041, 0.014, 0.004, 0.001)

# tabs will be a list of two tables with same marginals
tabs <- lapply(r2dtable(2, 1000 * a / sum(a), 1000 * b / sum(b)), "/",
1000); tabs

# row

Re: [R] help files

2011-02-19 Thread Uwe Ligges



On 19.02.2011 06:08, Erin Hodgess wrote:

Dear R People:

I downloaded the binary for R-2.12.1 on Windows 7.

When checking help files, I got the following:


?mean

starting httpd help server ... done
Error in shell.exec(url) :
   problem in displaying 'http://127.0.0.1:16945/library/base/html/mean.html'




I did check the URL in my IE browser and all was well.

Any suggestions would be much appreciated.



Perhaps some policy applies to your machine that prevents from accessing 
the http help server (on your machine)?


Uwe







Sincerely,
Erin




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

2011-02-19 Thread Uwe Ligges

On 18.02.2011 11:02, Rosario Garcia Gil wrote:

Hello

I have a multiple linear regression with two cofactors, I would like to 
represent a plane but I could not find any help which worked out.

Any suggestions.



One way is  explained in:

library("scatterplot3d")
?scatterplot3d

Now see example 5.

Uwe Ligges





Regards and thanks in advance.
Rosario
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Integration with an Indicator Function in R

2011-02-19 Thread Uwe Ligges



On 17.02.2011 22:10, li li wrote:

Hi all,
I have some some problem with regard to finding the integral of a
function containing an indicator function.
please see the code below:

func1<- function(x, mu){
   (mu^2)*dnorm(x, mean = mu, sd = 1)*dgamma(x, shape=2)}
m1star<- function(x){
integrate(func1, lower = 0, upper = Inf,x)$val}
T<- function(x){
0.3*dnorm(x)/(0.3*dnorm(x)+0.7*m1star(x))}

func2<- function(x,c){(T(x)<=c)*0.3*dnorm(x)}
func3<- function(x,c){(T(x)<= c)*(0.3*dnorm(x)+0.7*m1star(x))}
numer<- function(c){
integrate(func2, -Inf, Inf, c)$val}

denom<- function(c){
integrate(func3, lower, Inf,c)$val}



The error message is as below :


numer(0.5)

Error in integrate(func1, lower = 0, upper = Inf, x) :
   the integral is probably divergent



Several problems:

1. Your m1star does not work vectorized which is required here.
This can be fixed by changing the function as follows:

m1star <- function(x){
sapply(x, function(x) integrate(func1, lower = 0, upper = Inf, x = 
x)$val)

}

2. You will find that func2 calculates T(x). Note that T(x) is not 
undefined for x=0.


I have not investigated further, but my guess is that your problem needs 
some rethinking ...


Best,
Uwe Ligges






Thank you very much!
Hannah

[[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] [Rcpp-devel] Help with integrating R and c/c++

2011-02-19 Thread Rohit Pandey
Hi Douglas,

Sorry for leaving that information out earlier.

I am running windows XP. I'm not sure about the tools you mention. I thought
installing the packages in R was enough. Are these tools like a program you
install?

The results of running sessionInfo() are:

R version 2.12.1 (2010-12-16)
Platform: i386-pc-mingw32/i386 (32-bit)

locale:
[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United
States.1252LC_MONETARY=English_United States.1252
LC_NUMERIC=C
[5] LC_TIME=English_United States.1252

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



On Sat, Feb 19, 2011 at 9:21 PM, Douglas Bates  wrote:

> On Sat, Feb 19, 2011 at 12:56 AM, Rohit Pandey 
> wrote:
> > Hi Christopher/ Dirk,
>
> > Thank you very much for your replys. I think the idea of using inline as
> you
> > suggest is the best way to start off with using c++ with R. I went
> through
> > your examples and also plenty I found on the net. I have been trying to
> run
> > some of these over the past few days, but have consistently been getting
> > this error message as soon as I run the 'cfunction' or 'cxxfunction'. For
> > example, in the syntax that Dirk sent,
> > R> library(inline) #Runs perfectly
> >  R> src <- 'std::cout << "Hello C++_From_R World" << std::endl;
> > return(Rcpp::wrap(42));' #Runs perfectly
> >  R> rohit <- cxxfunction(signature(), src, plugin="Rcpp")
> > Now, as soon as I run this line, R spills a whole lot of text out and
> gives
> > an error and a warning:
> > ERROR(s) during compilation: source code errors or compiler configuration
> > errors!
> > Program source:
> >  1:
> >  2: // includes from the plugin
> >  3:
> >  4: #include 
> >  5:
> >  6:
> >  7: #ifndef BEGIN_RCPP
> >  8: #define BEGIN_RCPP
> >  9: #endif
> >  10:
> >  11: #ifndef END_RCPP
> >  12: #define END_RCPP
> >  13: #endif
> >  14:
> >  15: using namespace Rcpp;
> >  16:
> >  17:
> >  18: // user includes
> >  19:
> >  20:
> >  21: // declarations
> >  22: extern "C" {
> >  23: SEXP file59046688( ) ;
> >  24: }
> >  25:
> >  26: // definition
> >  27:
> >  28: SEXP file59046688(  ){
> >  29: BEGIN_RCPP
> >  30: std::cout << "Hello C++_From_R World" << std::endl;
> > return(Rcpp::wrap(42));
> >  31: END_RCPP
> >  32: }
> >  33:
> >  34:
> > Error in compileCode(f, code, language = language, verbose = verbose) :
> >  Compilation ERROR, function(s)/method(s) not created!
> > In addition: Warning message:
> > running command 'C:\PROGRA~1\R\R-212~1.1/bin/i386/R CMD SHLIB
> > file59046688.cpp 2> file59046688.cpp.err.txt' had status 1
> >
> > The "file59046688.cpp 2"  changes every time I run a different function,
> but
> > the problem seems to be the same.
> >
> > I installed and loaded the inline package (0.3.8) and then the Rcpp
> package
> > (0.9.0). I also tried reversing the order in which I load these, but
> still
> > no luck. I think if I can get just one of these programs to work, I will
> be
> > on my way. Can any of you tell me what I might be doing wrong?
> >
> > For your question on what exacly I require, Christopher - I just need to
> use
> > a while loop. I have always been able to substitute for loops with some
> of
> > the apply functions in R, but can't seem to be able to replace the while
> > with a more efficient function. But the things that are required inside
> the
> > while loop, I have already implemented in R efficiently. So, I thought of
> > transfering just the while loop to a language that is faster with loops.
>
> What operating system and version of R are you using?  It would help
> if you included the results of executing
>
> sessionInfo()
>
> in R.  More importantly, do you have the compiler tools installed and
> configured?  You need to have certain tools installed on Windows or
> Mac OS X before you can compile packages (as opposed to installing
> pre-compiled binary packages).  Most Linux distributions assume that
> their users are adults and provide them with the tools to compile
> programs.
>
> >
> > Thanks in advance.
> > On Mon, Feb 7, 2011 at 6:21 AM, Wray, Christopher <
> > christopher.wray...@ucl.ac.uk> wrote:
> >
> >> As Dirk says, using "inline" makes it real simple to start and to
> prototype
> >> code.
> >>
> >> You mention you have R functions you wish to "call" via Rcpp. Im not
> >> certain I fully understand what you require here, but it is pretty
> simple to
> >> pass R-side functions to C++ via Rcpp, and similarly its simple to send
> >> compiled functions back to R-side as external pointers (and reuse them
> >> elsewhere). Here is a toy example using a simple "user" function defined
> on
> >> the R-side (user_F).
> >>
> >> This is evaluated in R, passed as a "function" parameter to compiled C++
> >> function and evaluated there, and then a compiled version of the
> function is
> >> passed back to R as an external pointer, which you can send back to the
> C
> >> side and evaluate:
> >>
> >>  user_F=function(v){sum((1-v)*exp(-0.5*v))}
> >>

Re: [R] Confidence Intervals on Standard Curve

2011-02-19 Thread Ben Ward
Hi Graham,

Thanks, that does explain lots. I've been playing with making log's of 
data in models to make the relationship linear, which it does, which 
suggests to me that lm() is the right way to go, however, after if  try 
to predict after y values after about 60% on the x axis for light 
transmission, the y value, for bacterial numbers, crosses the axis and 
gives me negative values for y, which on a practical level isn't 
possible, as one can't have less than no bacteria in a culture. On a 
practical level, when I include the cirve in my appendix I could say 
anything above around 60% is 0, and mention the negative results from 
the standard curve's prediction capabilities are not literal, and to say 
turn any negative bacterial count obtained as a result of the curve to 
0. I've not had to deal with such pleatauing curves before. The values I 
have for the curve don't go above 50%, so anything above it is 
prediction, and my experiment probably won't result in x values above 
50% as the death of the culture proceeds slowly, but that depending on 
relative amounts of culture and antimicrobial I use the rate 'could' go 
faster or slower, so could go above 50%. I was wondering if non-linear 
regression is better for such a thing, but I'm hesitant to go into it in 
more detail for now because of the danger of drastically increacing 
complexity, if on a practical level, what I currenly have, works and is 
very accurate, within the range I will most likely be using it.

Thanks,
Ben.


On 19/02/2011 15:39, Graham Smith wrote:
> Ben,
>
> Does this help.
>
> http://r-eco-evo.blogspot.com/2011/01/confidence-intervals-for-regression.html
>
> Not sure if it will work with your particular model, but may be worth 
> a try.
>
>
> Graham
>
> On 18 February 2011 23:29, Ben Ward  > wrote:
>
> Hi, I wonder if anyone could advise me with this:
>
> I've been trying to make a standard curve in R with lm() of some
> standards from a spectrophotometer, so as I can express the curve
> as a formula, and so obtain values from my treated samples by
> plugging in readings into the formula, instead of trying to judge
> things by eye, with a curve drawn by hand.
>
> It is a curve and so I used the following formula:
>
> model <- lm(Approximate.Counts~X..Light.Transmission +
> I(Approximate.Counts^2), data=Standards)
>
> It gives me a pretty decent graph:
> xyplot(Approximate.Counts + fitted(model) ~ X..Light.Transmission,
> data=Standards)
>
> I'm pretty happy with it, and looking at the model summary, to my
> inexperienced eyes it seems pretty good:
>
> lm(formula = Approximate.Counts ~ X..Light.Transmission +
> I(Approximate.Counts^2),
>data = Standards)
>
> Residuals:
>   Min 1Q Median 3QMax
> -91.75 -51.04  27.33  37.28  49.72
>
> Coefficients:
>  Estimate Std. Error t value Pr(>|t|)
> (Intercept)  9.868e+02  2.614e+01   37.75 <2e-16 ***
> X..Light.Transmission   -1.539e+01  8.116e-01  -18.96 <2e-16 ***
> I(Approximate.Counts^2)  2.580e-04  6.182e-06   41.73 <2e-16 ***
> ---
> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>
> Residual standard error: 48.06 on 37 degrees of freedom
> Multiple R-squared: 0.9956,Adjusted R-squared: 0.9954
> F-statistic:  4190 on 2 and 37 DF,  p-value: < 2.2e-16
>
> I tried to put some 95% confidence interval lines on a plot, as
> advised by my tutor, to see how they looked, and I used a function
> I found in "The R Book":
>
> se.lines <- function(model){
> b1<-coef(model)[2]+ summary(model)[[4]][4]
> b2<-coef(model)[2]- summary(model)[[4]][4]
> xm<-mean(model[[12]][2])
> ym<-mean(model[[12]][1])
> a1<-ym-b1*xm
> a2<-ym-b2*xm
> abline(a1,b1,lty=2)
> abline(a2,b2,lty=2)
> }
> se.lines(model)
>
> but when I do this on a plot I get an odd result:
>
>
> They looks to me, to lie in the same kind of area, that my
> regression line did, before I used polynomial regression, by
> squaring "Approximate.Counts":
>
> lm(formula = Approximate.Counts ~ X..Light.Transmission +
> I(Approximate.Counts^2), data = Standards)
>
> Is there something else I should be doing? I've seen several ways
> of dealing with non-linear relationships, from log's of certain
> variables, and quadratic regression, and using sin and other
> mathematical devices. I'm not completely sure if I'm "allowed" to
> square the y variable, the book only squared the x variable in
> quadratic regression, which I did first, and it fit quite well,
> but not as good squaring Approximate Counts does:
>
> model <- lm(Approximate.Counts~X..Light.Transmission +
> I(X..Light.Transmission^2), data=Standards)
>
>
> Any advice is greatly appreciated, it's the first time I've really
> had to look at regression wi

[R] Sample function in R???

2011-02-19 Thread Hongwei Dong
Hi, R users,

I'm wondering if there a way in R I can select cases based on a probability
vector. if a case is selected, that case is marked as 1, otherwise, 0.

For example:

x<-12:18
y<-1:7
sample(x,2,replace=FALSE,y)

I got:

[1] 15 17

What I want to see is:

[1] 0 0 0 1 0 1 0


Thanks.
Gary

[[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] Confidence Intervals on Standard Curve

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 1:08 PM, Ben Ward wrote:


Hi Graham,

Thanks, that does explain lots. I've been playing with making log's of
data in models to make the relationship linear, which it does, which
suggests to me that lm() is the right way to go, however, after if   
try

to predict after y values after about 60% on the x axis for light
transmission, the y value, for bacterial numbers, crosses the axis and
gives me negative values for y, which on a practical level isn't
possible, as one can't have less than no bacteria in a culture.


Once you have your estimated parameter in the transformed "analysis"  
scale, you need to apply the inverse transformation, in this case  
exp,   to return the estimate to the measured scale. You also need to  
consider that in the process of transforming the values, performing an  
additive estiamtion, and transforming back to the "natural" scale, you  
will have estimated a multiplicative model, since exp(log(x)+log(y)) =  
x*y.


Had you used log(x) inside the lm call and then used predict, the  
predictions would be correct.


You might want to consider glm models with family="poisson".

--
David.


On a
practical level, when I include the cirve in my appendix I could say
anything above around 60% is 0, and mention the negative results from
the standard curve's prediction capabilities are not literal, and to  
say

turn any negative bacterial count obtained as a result of the curve to
0.


I wouldn't. The informed observer would know you were flailing around.


I've not had to deal with such pleatauing curves before. The values I
have for the curve don't go above 50%, so anything above it is
prediction, and my experiment probably won't result in x values above
50% as the death of the culture proceeds slowly, but that depending on
relative amounts of culture and antimicrobial I use the rate 'could'  
go

faster or slower, so could go above 50%. I was wondering if non-linear
regression is better for such a thing, but I'm hesitant to go into  
it in

more detail for now because of the danger of drastically increacing
complexity, if on a practical level, what I currenly have, works and  
is

very accurate, within the range I will most likely be using it.

Thanks,
Ben.


On 19/02/2011 15:39, Graham Smith wrote:

Ben,

Does this help.

http://r-eco-evo.blogspot.com/2011/01/confidence-intervals-for-regression.html

Not sure if it will work with your particular model, but may be worth
a try.


Graham

On 18 February 2011 23:29, Ben Ward mailto:benjamin.w...@bathspa.org>> wrote:

   Hi, I wonder if anyone could advise me with this:

   I've been trying to make a standard curve in R with lm() of some
   standards from a spectrophotometer, so as I can express the curve
   as a formula, and so obtain values from my treated samples by
   plugging in readings into the formula, instead of trying to judge
   things by eye, with a curve drawn by hand.

   It is a curve and so I used the following formula:

   model <- lm(Approximate.Counts~X..Light.Transmission +
   I(Approximate.Counts^2), data=Standards)

   It gives me a pretty decent graph:
   xyplot(Approximate.Counts + fitted(model) ~ X..Light.Transmission,
   data=Standards)

   I'm pretty happy with it, and looking at the model summary, to my
   inexperienced eyes it seems pretty good:

   lm(formula = Approximate.Counts ~ X..Light.Transmission +
   I(Approximate.Counts^2),
  data = Standards)

   Residuals:
 Min 1Q Median 3QMax
   -91.75 -51.04  27.33  37.28  49.72

   Coefficients:
Estimate Std. Error t value Pr(>|t|)
   (Intercept)  9.868e+02  2.614e+01   37.75 <2e-16 ***
   X..Light.Transmission   -1.539e+01  8.116e-01  -18.96 <2e-16 ***
   I(Approximate.Counts^2)  2.580e-04  6.182e-06   41.73 <2e-16 ***
   ---
   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

   Residual standard error: 48.06 on 37 degrees of freedom
   Multiple R-squared: 0.9956,Adjusted R-squared: 0.9954
   F-statistic:  4190 on 2 and 37 DF,  p-value: < 2.2e-16

   I tried to put some 95% confidence interval lines on a plot, as
   advised by my tutor, to see how they looked, and I used a function
   I found in "The R Book":

   se.lines <- function(model){
   b1<-coef(model)[2]+ summary(model)[[4]][4]
   b2<-coef(model)[2]- summary(model)[[4]][4]
   xm<-mean(model[[12]][2])
   ym<-mean(model[[12]][1])
   a1<-ym-b1*xm
   a2<-ym-b2*xm
   abline(a1,b1,lty=2)
   abline(a2,b2,lty=2)
   }
   se.lines(model)

   but when I do this on a plot I get an odd result:


   They looks to me, to lie in the same kind of area, that my
   regression line did, before I used polynomial regression, by
   squaring "Approximate.Counts":

   lm(formula = Approximate.Counts ~ X..Light.Transmission +
   I(Approximate.Counts^2), data = Standards)

   Is there something else I should be doing? I've seen several ways
   of dealing with non-linear relationships, from log's of certain
   variables, and qua

Re: [R] multiple regression plane

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 12:03 PM, Uwe Ligges wrote:


On 18.02.2011 11:02, Rosario Garcia Gil wrote:

Hello

I have a multiple linear regression with two cofactors, I would  
like to represent a plane but I could not find any help which  
worked out.


Any suggestions.



One way is  explained in:

library("scatterplot3d")
?scatterplot3d

Now see example 5.


Multiple examples of this one and other Ligges' work can also be found  
in the R Graphics Gallery:


http://addictedtor.free.fr/graphiques/thumbs.php

A nice 3d intro:

http://www.polisci.ohio-state.edu/faculty/lkeele/3dinR.pdf


Not even indirectly related to your request for planar  
representations, but a very nice collection of R graphics related to  
spatial processes:


http://r-spatial.sourceforge.net/gallery/

And a rather amazing graphics compendium that is apparently able to  
create any package's graphics examples in the fly. When you enter the  
website you get what appears to be a random example, but can then  
navigate to the package of your choice:


http://bm2.genes.nig.ac.jp/RGM2/





Uwe Ligges





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


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] Sample function in R???

2011-02-19 Thread Peter Ehlers

On 2011-02-19 10:32, Hongwei Dong wrote:

Hi, R users,

I'm wondering if there a way in R I can select cases based on a probability
vector. if a case is selected, that case is marked as 1, otherwise, 0.

For example:

x<-12:18
y<-1:7
sample(x,2,replace=FALSE,y)

I got:

[1] 15 17

What I want to see is:

[1] 0 0 0 1 0 1 0


Maybe:

 (x %in% sample(x,2,replace=FALSE,y)) * 1

Peter Ehlers




Thanks.
Gary

[[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] Porting existing Java GUI to R

2011-02-19 Thread Roy Mendelssohn
Hi All:

we have an existing GUI written in Java that allows users to access a slew of 
remote satellite and other environmental data.  The GUI allows the use to 
browse the source, select the dataset, time and spatial extent of the data, and 
then obtain that subset.  There is a standalone version, a version that maps to 
ArcGIS, and in the works a version that maps into Matlab.

We are looking at trying to port this to R, and there appears to be several 
options and I am hoping for suggestions as to which would be the best given 
what we are trying to achieve.  For example, incorporating java code into 
Matlab is straightforward, so we can create a function that calls up the GUI, 
and then on "close", maps the retrieved data  into Matlab data structures.  The 
user is still in their usual Matlab workspace - that is what we prefer - we do 
not want if possible a standalone analysis environment, just an extension to 
the "Open" statement as it were.

One option appears to be using rjava to make the call to open the GUI  (I 
believe this would work but I am not certain) and then use jri to write the 
data structures into the present R workspace.  Another option is to do 
something like Bio7 does, where there is an R "window" inside the java 
application, but if I understand the code of that correctly this is done by 
starting up and calling an RServer.  We would prefer just getting the data into 
the R workspace and then getting out of the way, rather than running RServer or 
something similar.

So I guess my questions are:

1.  Can rjava and jri do the types of things we are contemplating.

2.  Are there any other options besides RServer?

3.  Any existing examples of doing something similar that we might be able to 
look at the code?

Thanks,

-Roy M.




**
"The contents of this message do not reflect any position of the U.S. 
Government or NOAA."
**
Roy Mendelssohn
Supervisory Operations Research Analyst
NOAA/NMFS
Environmental Research Division
Southwest Fisheries Science Center
1352 Lighthouse Avenue
Pacific Grove, CA 93950-2097

e-mail: roy.mendelss...@noaa.gov (Note new e-mail address)
voice: (831)-648-9029
fax: (831)-648-8440
www: http://www.pfeg.noaa.gov/

"Old age and treachery will overcome youth and skill."
"From those who have been given much, much will be expected" 

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

2011-02-19 Thread Ben Ward

Hi David,

I had use log(x)inside the lm call and used predict, although I didn't 
know about logs of data making a multiplacative model  
exp(log(x)+log(y)) = x*y.
I'll have a look at the poisson model. and see what I manage to produce. 
Looking at the internet the Cumulative distribution function looks like 
a sort of S with pleataus at the top and bottom, which I suppose my 
standard curve would be a part/chunk of? And if I predicted outside of 
what I have, it might become like the full S. I apologise for not 
putting it very elequently in mathematical terms, I'm a biology student 
so stats and hypothesis testing is taught, but doing this standard curve 
is unprecedented in my work so far, normally lm() and anova and such is 
enough.


Thanks,
Ben.

On 19/02/2011 19:06, David Winsemius wrote:


On Feb 19, 2011, at 1:08 PM, Ben Ward wrote:


Hi Graham,

Thanks, that does explain lots. I've been playing with making log's of
data in models to make the relationship linear, which it does, which
suggests to me that lm() is the right way to go, however, after if  try
to predict after y values after about 60% on the x axis for light
transmission, the y value, for bacterial numbers, crosses the axis and
gives me negative values for y, which on a practical level isn't
possible, as one can't have less than no bacteria in a culture.


Once you have your estimated parameter in the transformed "analysis" 
scale, you need to apply the inverse transformation, in this case 
exp,   to return the estimate to the measured scale. You also need to 
consider that in the process of transforming the values, performing an 
additive estiamtion, and transforming back to the "natural" scale, you 
will have estimated a multiplicative model, since exp(log(x)+log(y)) = 
x*y.


Had you used log(x) inside the lm call and then used predict, the 
predictions would be correct.


You might want to consider glm models with family="poisson".



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

2011-02-19 Thread Peter Ehlers

On 2011-02-19 07:35, garciap wrote:


Hi to all the people (again),

I'm doing some simulations with the memisc package of an own function, but
I've the problem that I didn't know how to read the result of such
simulations. My function is:


Torre<-function(a1,N1,a2,N2)

+ {Etorre<-(a1*N1)/(1+a1*N1)
+ Efuera<-(a2*N2)/(1+a2*N2)
+ if(Etorre>Efuera)Subir=TRUE
+ if (Etorre
Torre<-Simulate(Torre(a1,N1,a2,N2),expand.grid(a1=3,N1=0.5,a2=c(0.01,0.02,0.05,0.1),N2=0.1),nsim=1000,seed=1,trace=50,keep.data=TRUE)




[... snip ...]


print(Torre)

A 'default_bucket' object with 5 variables and 4000 observations

When I try print the results I obtain this, but how can I read it? It will
seem stupid to people, but I can't! even though function summary(Torre):


summary(Torre)

Error in object[[i]] : wrong arguments for subsetting an environment




You should be able to wrap your Simulate output in as.data.frame:

 dat <- as.data.frame(Simulate())

Peter Ehlers




Best regards (again) and many thanks for your help

Pablo


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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, 'se' variable in geom_errorbar's limits?

2011-02-19 Thread Eric Fail
Can't anybody give me a hint on how to solve this? I even bought the
ggplot2-book, so you could also give a page (or a series of pages).

Thanks,
Eric

On Thu, Feb 17, 2011 at 10:19 AM, Eric Fail  wrote:
>
> Dear R-list
>
> I'm working with with geom_errorbar; specifically I'm trying to
> reproduce the example Hadley Wickham have on
> http://had.co.nz/ggplot2/geom_errorbar.html (all in the button of the
> page) where he makes an nice plot with errorbars and then draw lines
> between the points.
>
> What confuses me is the 'limits' he defines for the errorbars from the
> se variable.
>
> First he creates a dataset,
>
> df <- data.frame(
>  trt = factor(c(1, 1, 2, 2)),
>  resp = c(1, 5, 3, 4),
>  group = factor(c(1, 2, 1, 2)),
>  se = c(0.1, 0.3, 0.3, 0.2)
> )
>
> # library(ggplot2)
>
> and then he creates some limits from the se variables.
>
> limits <- aes(ymax = resp + se, ymin=resp - se)
>
> [elements omitted]
>
> # and then he creates the plot (I'm interested in).
>
> p <- ggplot(df, aes(colour=group, y=resp, x=trt))
> p + geom_line(aes(group=group)) + geom_errorbar(limits, width=0.2)
>
> I can (of course) get Hadley's example to run, but I can't do it on my
> data as I don't have a 'se' variable/don't know how to create it. I
> have a group variable, a treatment variable, and a response variable,
> but no se variable.
>
> Could anyone out there explain how I create a 'se' variable in my data?
>
> I'm sure my reasoning is the one that is off, and not ggplot2 (I'm a big fan).
>
> Your help is appreciated!
>
> Thanks,
> Eric

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

2011-02-19 Thread Krishnan Viswanathan
I am trying to recode a variable into another variable and while the package
'car' works well when it is only recoding A into B, I am not sure I can do
the same with recoding (A or C) into B. If i can use recode please advise on
how to. So i am using an if/else if conditions.

My sample dataset is below along with the code and the warning and results i
get.

TIA
Krishnan

#Code
attach(airport_survey)
if(AccessMode == 1 || AccessMode == 8 || AccessMode == 9) {
   airport_survey$amodecat <- 1
} else if (AccessMode == 2) {
   airport_survey$amodecat <- 2
} else if (AccessMode == 3 & rentalcat == 1 || AccessMode == 4 & rentalcat
== 1) {
airport_survey$amodecat <- 3
} else if (AccessMode == 5) {
   airport_survey$amodecat <- 8
} else if (AccessMode == 6 || AccessMode == 14) {
   airport_survey$amodecat <- 9
} else if (AccessMode == 7) {
   airport_survey$amodecat <- 5
} else if (AccessMode == 10) {
   airport_survey$amodecat <- 7
} else if (AccessMode == 11 || AccessMode == 12) {
   airport_survey$amodecat <- 4
} else if (AccessMode == 13) {
   airport_survey$amodecat <- 6
} else if(AccessMode == 3 & rentalcat == 2 ||AccessMode == 4 & rentalcat ==
2) {
   airport_survey$amodecat <- 10
} else { # else statement here just to close out the condition
   airport_survey$amodecat <- 0
}


#*dataset

AccessMode RentalCat
4 2
4 1
4 2
3 1
3 2
14 1
1 1
2 1
8 1
9 1
10 1
11 1
12 1
13 1
6 1
5 1
7 2

# Warning + Results

Warning messages:
1: In if (AccessMode == 2) { :
  the condition has length > 1 and only the first element will be used
2: In if (AccessMode == 5) { :
  the condition has length > 1 and only the first element will be used
3: In if (AccessMode == 7) { :
  the condition has length > 1 and only the first element will be used
4: In if (AccessMode == 10) { :
  the condition has length > 1 and only the first element will be used
5: In if (AccessMode == 13) { :
  the condition has length > 1 and only the first element will be used

# After Running summary on the variable (it seems to only take the last if
condition)

   amodecat
 Min.   :10
 1st Qu.:10
 Median :10
 Mean   :10
 3rd Qu.:10
 Max.   :10

-- 
Krishnan Viswanathan
1101 High Meadow Dr
Tallahassee FL 32311

[[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] Seeking help in Package development

2011-02-19 Thread Nipesh Bajaj
Dear all, I am a new user of R and currently trying hard to develop my
own package. Here I am following this tutorial
'http://www.mathfinance.cn/how-to-create-an-R-package-in-windows/'

Here it says that (step 8): "open a “command prompt” window, change
the directory to where your package is, type the command “R CMD build
MonteCarloPi” to build the package, this will generate a file called
MonteCarloPi_1.0.tar.gz. "

According to that, I have opened the Windows command prompt window (a
black screen window) and then changed the directory, where my new
package (a folder in current working directory in R, as created by
'package.skeleton') is there. Then typed 'R CMD build MyPackage' (I
named my package as 'MyPackage'). However doing so I got following
error in that command prompt:

'R' is not recognized as an internal or external command, operable
program or batch file.

Can somebody please guide me what to do in this situation?

Best thanks

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

2011-02-19 Thread Gabor Grothendieck
On Sat, Feb 19, 2011 at 4:39 PM, Nipesh Bajaj  wrote:
> Dear all, I am a new user of R and currently trying hard to develop my
> own package. Here I am following this tutorial
> 'http://www.mathfinance.cn/how-to-create-an-R-package-in-windows/'
>
> Here it says that (step 8): "open a “command prompt” window, change
> the directory to where your package is, type the command “R CMD build
> MonteCarloPi” to build the package, this will generate a file called
> MonteCarloPi_1.0.tar.gz. "
>
> According to that, I have opened the Windows command prompt window (a
> black screen window) and then changed the directory, where my new
> package (a folder in current working directory in R, as created by
> 'package.skeleton') is there. Then typed 'R CMD build MyPackage' (I
> named my package as 'MyPackage'). However doing so I got following
> error in that command prompt:
>
> 'R' is not recognized as an internal or external command, operable
> program or batch file.
>
> Can somebody please guide me what to do in this situation?
>

You have to modify you path, adding the directory containing R.exe so
that Windows can find R.

Alternately you can download R.bat from batchfiles.googlecode.com and
place it anywhere on your path to avoid having to modify your path
itself.  It will find R from the registry and run it with the same
arguments.


-- 
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] barplot, different color for shading lines and bar

2011-02-19 Thread Markus Loecher
Dear all,
might there be a modified barplot function out there which allows the user
to specify a fill color for the bars and independent parameters for the
overlaid shading lines ?
Currently, when I specify density and col, the fill color for the bars is
white.

Thanks!

Markus

[[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] Seeking help in Package development

2011-02-19 Thread Joshua Wiley
Dear Nipesh,

On Sat, Feb 19, 2011 at 1:39 PM, Nipesh Bajaj  wrote:
> Dear all, I am a new user of R and currently trying hard to develop my
> own package. Here I am following this tutorial
> 'http://www.mathfinance.cn/how-to-create-an-R-package-in-windows/'
>
> Here it says that (step 8): "open a “command prompt” window, change
> the directory to where your package is, type the command “R CMD build
> MonteCarloPi” to build the package, this will generate a file called
> MonteCarloPi_1.0.tar.gz. "
>
> According to that, I have opened the Windows command prompt window (a
> black screen window) and then changed the directory, where my new
> package (a folder in current working directory in R, as created by
> 'package.skeleton') is there. Then typed 'R CMD build MyPackage' (I
> named my package as 'MyPackage'). However doing so I got following
> error in that command prompt:
>
> 'R' is not recognized as an internal or external command, operable
> program or batch file.

This suggests that either your current directory does not contain R
and you have not added the directory containing R to the PATH
environment variable in Windows.  For that command to work, Windows
needs to know where to find the program to execute.  You might find
this site useful for becoming more familiar with the Windows command
prompt (geared towards XP, but cmd.exe has changed little from XP to
Vista, to 7---though increasing MS is encouraging users to switch over
to using the Windows powershell):

http://www.microsoft.com/resources/documentation/windows/xp/all/proddocs/en-us/ntcmds.mspx?mfr=true

>
> Can somebody please guide me what to do in this situation?

You should also read the official R manual on installing and building
R from source.  Particularly pay attention to the sections for Windows
users.  You may want to get Rtools and certainly follow the
instructions to add all the necessary directories to your PATH
variable.  The steps will be something like:

Right click My Computer -> Click Properties -> Click Advanced -> Click
Environment Variables -> edit the variable "PATH" to include relevant
R directories.  Please note that this is not an exact step-by-step
process as it varies slightly by different versions of Windows.

Official R Installation Manual (very relevant for building packages):
http://cran.r-project.org/doc/manuals/R-admin.html

Good luck,

Josh

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



-- 
Joshua Wiley
Ph.D. Student, Health Psychology
University of California, Los Angeles
http://www.joshuawiley.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] Accessing DF index

2011-02-19 Thread eric

I have a dataframe called x2. It seems to have a date column but I can't
access it or give it a name or convert it to a date. How do I refer to that
first column and make it a date ? When I try x2[1,] I get the second column.

head(x2)
  FAIRXSP500delta
2000-08-31  0.010101096  0.007426964  0.002674132
2000-09-29  0.096679730 -0.054966292  0.151646023
2000-10-31 -0.008245580 -0.004961785 -0.003283795
2000-11-30  0.037024545 -0.083456134  0.120480679
2000-12-29  0.080042708  0.004045193  0.075997514
2001-01-31 -0.009042396  0.034050246 -0.043092643

x[1,1]
FAIRX
2000-08-31 0.01010110

x[1,2]
 SP500
2000-08-31 0.007426964

str(x2)
'data.frame':   127 obs. of  3 variables:
 $ FAIRX: num  0.0101 0.09668 -0.00825 0.03702 0.08004 ...
 $ SP500: num  0.00743 -0.05497 -0.00496 -0.08346 0.00405 ...
 $ delta: num  0.00267 0.15165 -0.00328 0.12048 0.076 ..
-- 
View this message in context: 
http://r.789695.n4.nabble.com/Accessing-DF-index-tp3314649p3314649.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] Conditional recoding

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 4:36 PM, Krishnan Viswanathan wrote:

I am trying to recode a variable into another variable and while the  
package
'car' works well when it is only recoding A into B, I am not sure I  
can do
the same with recoding (A or C) into B. If i can use recode please  
advise on

how to. So i am using an if/else if conditions.


If you are using if-else logic to transform vectors, you should use  
the ifelse function which will work on vectors rather than the if() 
{}else{} control structure which only accepts a single element. You  
should also learn to distinguish the "|" operator that retruns vectors  
from vectors, versus the "||" operator which is used with if(). Notice  
that it would be improper to pair "||" operations with "&"


(All this is explained in ?ifelse .)


My sample dataset is below along with the code and the warning and  
results i

get.

TIA
Krishnan

#Code
attach(airport_survey)


If airport_survey is an available dataset then you should specify  
which package contains it.



if(AccessMode == 1 || AccessMode == 8 || AccessMode == 9) {
  airport_survey$amodecat <- 1
} else if (AccessMode == 2) {
  airport_survey$amodecat <- 2
} else if (AccessMode == 3 & rentalcat == 1 || AccessMode == 4 &  
rentalcat

== 1) {
   airport_survey$amodecat <- 3
} else if (AccessMode == 5) {
  airport_survey$amodecat <- 8
} else if (AccessMode == 6 || AccessMode == 14) {
  airport_survey$amodecat <- 9
} else if (AccessMode == 7) {
  airport_survey$amodecat <- 5
} else if (AccessMode == 10) {
  airport_survey$amodecat <- 7
} else if (AccessMode == 11 || AccessMode == 12) {
  airport_survey$amodecat <- 4
} else if (AccessMode == 13) {
  airport_survey$amodecat <- 6
} else if(AccessMode == 3 & rentalcat == 2 ||AccessMode == 4 &  
rentalcat ==

2) {
  airport_survey$amodecat <- 10
} else { # else statement here just to close out the condition
  airport_survey$amodecat <- 0
}


#*dataset


And if this is what you are calling airport_survey then you should  
offer code that will read it in:


airport_survey <- read.table(textConnection("AccessMode RentalCat
4 2
4 1
4 2
3 1
3 2
14 1
1 1
2 1
8 1
9 1
10 1
11 1
12 1
13 1
6 1
5 1
7 2"), header=TRUE)

May want to do this in stages although you can nest up to 7 ifelse's.  
I will illustrate the first couple:

(Note: I NEVER USE attach(), use with() instead.

airport_survey$amodecat <- NA
airport_survey$amodecat <- with(airport_survey, ifelse( AccessMode %in 
% c(1,8,9) ,  1, amodecat)
airport_survey$amodecat <- with(airport_survey, ifelse( AccessMode ==  
2 ,  2, amodecat)
airport_survey$amodecat <- with(airport_survey, ifelse( (AccessMode ==  
3 & rentalcat == 1) | (AccessMode == 4 & rentalcat == 1) ,  3, amodecat)

... and so on.




# Warning + Results

Warning messages:
1: In if (AccessMode == 2) { :
 the condition has length > 1 and only the first element will be used
2: In if (AccessMode == 5) { :
 the condition has length > 1 and only the first element will be used
3: In if (AccessMode == 7) { :
 the condition has length > 1 and only the first element will be used
4: In if (AccessMode == 10) { :
 the condition has length > 1 and only the first element will be used
5: In if (AccessMode == 13) { :
 the condition has length > 1 and only the first element will be used

# After Running summary on the variable (it seems to only take the  
last if

condition)


Right. You made them all equal with the incorrect use of if(){}else{}.


  amodecat
Min.   :10
1st Qu.:10
Median :10
Mean   :10
3rd Qu.:10
Max.   :10

--
Krishnan Viswanathan
1101 High Meadow Dr



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] Accessing DF index

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 4:50 PM, eric wrote:



I have a dataframe called x2. It seems to have a date column but I  
can't
access it or give it a name or convert it to a date. How do I refer  
to that
first column and make it a date ? When I try x2[1,] I get the second  
column.


It is not actually a column but is rather a set of row names:

?rownames

rownames(x2)




head(x2)
 FAIRXSP500delta
2000-08-31  0.010101096  0.007426964  0.002674132
2000-09-29  0.096679730 -0.054966292  0.151646023
2000-10-31 -0.008245580 -0.004961785 -0.003283795
2000-11-30  0.037024545 -0.083456134  0.120480679
2000-12-29  0.080042708  0.004045193  0.075997514
2001-01-31 -0.009042396  0.034050246 -0.043092643

x[1,1]
   FAIRX
2000-08-31 0.01010110

x[1,2]
SP500
2000-08-31 0.007426964

str(x2)
'data.frame':   127 obs. of  3 variables:
$ FAIRX: num  0.0101 0.09668 -0.00825 0.03702 0.08004 ...
$ SP500: num  0.00743 -0.05497 -0.00496 -0.08346 0.00405 ...
$ delta: num  0.00267 0.15165 -0.00328 0.12048 0.076 ..
--
View this message in context: 
http://r.789695.n4.nabble.com/Accessing-DF-index-tp3314649p3314649.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.


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] Seeking help in Package development

2011-02-19 Thread Nipesh Bajaj
Thanks Gabor for your input. Here what I have done is that:

1. Copy 'MyPackage' folder (developed by package.skeleton) into
'C:\Program Files\R\R-2.12.1\bin' (I found R.exe is there)

2. In the command prompt, I changed the working directory using "CD"
command and run 'R CMD build MyPackage'

3. I have seen that a file named 'MyPackage_1.0.tar' has been created.
Then I pasted that file in ''C:\Program Files\R\R-2.12.1\bin'

4. Again run R CMD INSTALL MyPackage_1.0.tar.gz. However here I got
some error saying:
'Error in Rd_info(db[[i]]) : Rd files must have a non-empty \title'

In my Package, there are 2 user defined functions MyFunction1 &
MyFunction2, 1st function just get the length of it's argument and 2nd
function get the Range. I have modified the manual pages of those 2
functions as (just 1st few lines):

\name{MyFunction1}
\alias{MyFunction1}
%- Also NEED an '\alias' for EACH other topic documented here.
\title{
  ~~(Title) This function try to get the length of it's arguments ~~
}
\description{
  ~~ (Description) This function try to get the length of it's arguments ~~
}
\usage{
MyFunction1(x)
}

\name{MyFunction2}
\alias{MyFunction2}
%- Also NEED an '\alias' for EACH other topic documented here.
\title{
  ~~(Title) This function try to get the range of it's arguments ~~
}
\description{
%%  ~~ A concise (1-5 lines) description of what the function does. ~~
}
\usage{
MyFunction2(x)
}

I understand there is something to do with the 'title' segments of any
(or all) those 2 functions. However could not find what to do. Any
point please?

Actually I wanted to create .zip file as I am working in Vista. For
that I have intalled R_tool as well. However I could not understand
how I can use that. Any help will be highly appreciated.

Thanks,
On Sun, Feb 20, 2011 at 3:22 AM, Gabor Grothendieck
 wrote:
> On Sat, Feb 19, 2011 at 4:39 PM, Nipesh Bajaj  wrote:
>> Dear all, I am a new user of R and currently trying hard to develop my
>> own package. Here I am following this tutorial
>> 'http://www.mathfinance.cn/how-to-create-an-R-package-in-windows/'
>>
>> Here it says that (step 8): "open a “command prompt” window, change
>> the directory to where your package is, type the command “R CMD build
>> MonteCarloPi” to build the package, this will generate a file called
>> MonteCarloPi_1.0.tar.gz. "
>>
>> According to that, I have opened the Windows command prompt window (a
>> black screen window) and then changed the directory, where my new
>> package (a folder in current working directory in R, as created by
>> 'package.skeleton') is there. Then typed 'R CMD build MyPackage' (I
>> named my package as 'MyPackage'). However doing so I got following
>> error in that command prompt:
>>
>> 'R' is not recognized as an internal or external command, operable
>> program or batch file.
>>
>> Can somebody please guide me what to do in this situation?
>>
>
> You have to modify you path, adding the directory containing R.exe so
> that Windows can find R.
>
> Alternately you can download R.bat from batchfiles.googlecode.com and
> place it anywhere on your path to avoid having to modify your path
> itself.  It will find R from the registry and run it with the same
> arguments.
>
>
> --
> 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] Conditional recoding

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 5:14 PM, David Winsemius wrote:



On Feb 19, 2011, at 4:36 PM, Krishnan Viswanathan wrote:

I am trying to recode a variable into another variable and while  
the package
'car' works well when it is only recoding A into B, I am not sure I  
can do
the same with recoding (A or C) into B. If i can use recode please  
advise on

how to. So i am using an if/else if conditions.


If you are using if-else logic to transform vectors, you should use  
the ifelse function which will work on vectors rather than the if() 
{}else{} control structure which only accepts a single element. You  
should also learn to distinguish the "|" operator that retruns  
vectors from vectors, versus the "||" operator which is used with  
if(). Notice that it would be improper to pair "||" operations with  
"&"


(All this is explained in ?ifelse .)


My sample dataset is below along with the code and the warning and  
results i

get.

TIA
Krishnan

#Code
attach(airport_survey)


If airport_survey is an available dataset then you should specify  
which package contains it.



if(AccessMode == 1 || AccessMode == 8 || AccessMode == 9) {
 airport_survey$amodecat <- 1
} else if (AccessMode == 2) {
 airport_survey$amodecat <- 2
} else if (AccessMode == 3 & rentalcat == 1 || AccessMode == 4 &  
rentalcat

== 1) {
  airport_survey$amodecat <- 3
} else if (AccessMode == 5) {
 airport_survey$amodecat <- 8
} else if (AccessMode == 6 || AccessMode == 14) {
 airport_survey$amodecat <- 9
} else if (AccessMode == 7) {
 airport_survey$amodecat <- 5
} else if (AccessMode == 10) {
 airport_survey$amodecat <- 7
} else if (AccessMode == 11 || AccessMode == 12) {
 airport_survey$amodecat <- 4
} else if (AccessMode == 13) {
 airport_survey$amodecat <- 6
} else if(AccessMode == 3 & rentalcat == 2 ||AccessMode == 4 &  
rentalcat ==

2) {
 airport_survey$amodecat <- 10
} else { # else statement here just to close out the condition
 airport_survey$amodecat <- 0
}


#*dataset


And if this is what you are calling airport_survey then you should  
offer code that will read it in:


airport_survey <- read.table(textConnection("AccessMode RentalCat
4 2
4 1
4 2
3 1
3 2
14 1
1 1
2 1
8 1
9 1
10 1
11 1
12 1
13 1
6 1
5 1
7 2"), header=TRUE)

May want to do this in stages although you can nest up to 7  
ifelse's. I will illustrate the first couple:

(Note: I NEVER USE attach(), use with() instead.


Correction:

> airport_survey$amodecat <- NA
> airport_survey$amodecat <- with(airport_survey, ifelse( AccessMode  
%in% c(1,8,9) ,  1, amodecat) )
> airport_survey$amodecat <- with(airport_survey, ifelse( AccessMode  
== 2 ,  2, amodecat) )
> airport_survey$amodecat <- with(airport_survey, ifelse( (AccessMode  
== 3 & RentalCat == 1) | (AccessMode == 4 & RentalCat == 1) ,  3,  
amodecat) )





airport_survey$amodecat <- NA
airport_survey$amodecat <- with(airport_survey, ifelse( AccessMode  
%in% c(1,8,9) ,  1, amodecat)
airport_survey$amodecat <- with(airport_survey, ifelse( AccessMode  
== 2 ,  2, amodecat)
airport_survey$amodecat <- with(airport_survey, ifelse( (AccessMode  
== 3 & rentalcat == 1) | (AccessMode == 4 & rentalcat == 1) ,  3,  
amodecat)

... and so on.




# Warning + Results

Warning messages:
1: In if (AccessMode == 2) { :
the condition has length > 1 and only the first element will be used
2: In if (AccessMode == 5) { :
the condition has length > 1 and only the first element will be used
3: In if (AccessMode == 7) { :
the condition has length > 1 and only the first element will be used
4: In if (AccessMode == 10) { :
the condition has length > 1 and only the first element will be used
5: In if (AccessMode == 13) { :
the condition has length > 1 and only the first element will be used

# After Running summary on the variable (it seems to only take the  
last if

condition)


Right. You made them all equal with the incorrect use of if(){}else{}.


 amodecat
Min.   :10
1st Qu.:10
Median :10
Mean   :10
3rd Qu.:10
Max.   :10

--
Krishnan Viswanathan
1101 High Meadow Dr



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.


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] Accessing DF index

2011-02-19 Thread eric

So how would I convert those row names to dates and give that column the name
"Date" so that I can use subset and other functions on the Date column ?
-- 
View this message in context: 
http://r.789695.n4.nabble.com/Accessing-DF-index-tp3314649p3314689.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] Accessing DF index

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 5:25 PM, eric wrote:



So how would I convert those row names to dates and give that column  
the name
"Date" so that I can use subset and other functions on the Date  
column ?


> x2$dtcol <- as.Date(rownames(x2))
> x2
  FAIRXSP500delta  dtcol
2000-08-31  0.010101096  0.007426964  0.002674132 2000-08-31
2000-09-29  0.096679730 -0.054966292  0.151646023 2000-09-29
2000-10-31 -0.008245580 -0.004961785 -0.003283795 2000-10-31
2000-11-30  0.037024545 -0.083456134  0.120480679 2000-11-30
2000-12-29  0.080042708  0.004045193  0.075997514 2000-12-29
2001-01-31 -0.009042396  0.034050246 -0.043092643 2001-01-31


--
View this message in context: 
http://r.789695.n4.nabble.com/Accessing-DF-index-tp3314649p3314689.html
Sent from the R help mailing list archive at Nabble.com.


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] Accessing DF index

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 5:25 PM, eric wrote:



So how would I convert those row names to dates and give that column  
the name
"Date" so that I can use subset and other functions on the Date  
column ?


You could have used subset without the new column:

> subset(x2, rownames(x2) < "2000-11-30")
 FAIRXSP500delta  dtcol
2000-08-31  0.01010110  0.007426964  0.002674132 2000-08-31
2000-09-29  0.09667973 -0.054966292  0.151646023 2000-09-29
2000-10-31 -0.00824558 -0.004961785 -0.003283795 2000-10-31


--



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] ggplot2, 'se' variable in geom_errorbar's limits?

2011-02-19 Thread Ista Zahn
Hi Eric,
It's not really a ggplot question. If you are plotting means you can
calculate the standard error as sd/sqrt(n). Then add and subtract the
standard error to create ymax and ymin.

Best,
Ista

On Sat, Feb 19, 2011 at 9:12 PM, Eric Fail  wrote:
> Can't anybody give me a hint on how to solve this? I even bought the
> ggplot2-book, so you could also give a page (or a series of pages).
>
> Thanks,
> Eric
>
> On Thu, Feb 17, 2011 at 10:19 AM, Eric Fail  wrote:
>>
>> Dear R-list
>>
>> I'm working with with geom_errorbar; specifically I'm trying to
>> reproduce the example Hadley Wickham have on
>> http://had.co.nz/ggplot2/geom_errorbar.html (all in the button of the
>> page) where he makes an nice plot with errorbars and then draw lines
>> between the points.
>>
>> What confuses me is the 'limits' he defines for the errorbars from the
>> se variable.
>>
>> First he creates a dataset,
>>
>> df <- data.frame(
>>  trt = factor(c(1, 1, 2, 2)),
>>  resp = c(1, 5, 3, 4),
>>  group = factor(c(1, 2, 1, 2)),
>>  se = c(0.1, 0.3, 0.3, 0.2)
>> )
>>
>> # library(ggplot2)
>>
>> and then he creates some limits from the se variables.
>>
>> limits <- aes(ymax = resp + se, ymin=resp - se)
>>
>> [elements omitted]
>>
>> # and then he creates the plot (I'm interested in).
>>
>> p <- ggplot(df, aes(colour=group, y=resp, x=trt))
>> p + geom_line(aes(group=group)) + geom_errorbar(limits, width=0.2)
>>
>> I can (of course) get Hadley's example to run, but I can't do it on my
>> data as I don't have a 'se' variable/don't know how to create it. I
>> have a group variable, a treatment variable, and a response variable,
>> but no se variable.
>>
>> Could anyone out there explain how I create a 'se' variable in my data?
>>
>> I'm sure my reasoning is the one that is off, and not ggplot2 (I'm a big 
>> fan).
>>
>> Your help is appreciated!
>>
>> Thanks,
>> Eric
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
Ista Zahn
Graduate student
University of Rochester
Department of Clinical and Social Psychology
http://yourpsyche.org

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


Re: [R] Seeking help in Package development

2011-02-19 Thread Gabor Grothendieck
On Sat, Feb 19, 2011 at 5:31 PM, Nipesh Bajaj  wrote:
> Thanks Gabor for your input. Here what I have done is that:
>
> 1. Copy 'MyPackage' folder (developed by package.skeleton) into
> 'C:\Program Files\R\R-2.12.1\bin' (I found R.exe is there)
>
> 2. In the command prompt, I changed the working directory using "CD"
> command and run 'R CMD build MyPackage'
>
> 3. I have seen that a file named 'MyPackage_1.0.tar' has been created.
> Then I pasted that file in ''C:\Program Files\R\R-2.12.1\bin'
>
> 4. Again run R CMD INSTALL MyPackage_1.0.tar.gz. However here I got
> some error saying:
> 'Error in Rd_info(db[[i]]) : Rd files must have a non-empty \title'
>

What you have done in #1 is asking for trouble.

You do need to get used to debugging your Rd files and you should
expect to get many warnings and errors.  Read the messages carefully
and keep correcting them and rebuilding and checking until they pass.

-- 
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] Seeking help in Package development

2011-02-19 Thread Nipesh Bajaj
Hi Gabor, can you be more detail on step 01? Is that not the correct
path for R.exe? But I found it there.

When I ran R CMD build MyPackage then actually got these warnings:
cygwin warnings:
MS-DOS style path detected: c:/Program Files/R/R-2.12.1/bin/Mypackage_1.0.tar
prefered POSIX equivalent is /cygdrive/c/Program
Files/R/R-2.12.1/bin/MyPackage_1.0.tar
CYGWIN environment variable option "nodesfilewarning" turns off this
warning.

Is this warning something to do with that error?

Sorry for stretching this thread so long. However I am not really
expert in programming therefore find really hard on understanding
different terminology given in different documentation.

Thanks,


On Sun, Feb 20, 2011 at 4:35 AM, Gabor Grothendieck
 wrote:
> On Sat, Feb 19, 2011 at 5:31 PM, Nipesh Bajaj  wrote:
>> Thanks Gabor for your input. Here what I have done is that:
>>
>> 1. Copy 'MyPackage' folder (developed by package.skeleton) into
>> 'C:\Program Files\R\R-2.12.1\bin' (I found R.exe is there)
>>
>> 2. In the command prompt, I changed the working directory using "CD"
>> command and run 'R CMD build MyPackage'
>>
>> 3. I have seen that a file named 'MyPackage_1.0.tar' has been created.
>> Then I pasted that file in ''C:\Program Files\R\R-2.12.1\bin'
>>
>> 4. Again run R CMD INSTALL MyPackage_1.0.tar.gz. However here I got
>> some error saying:
>> 'Error in Rd_info(db[[i]]) : Rd files must have a non-empty \title'
>>
>
> What you have done in #1 is asking for trouble.
>
> You do need to get used to debugging your Rd files and you should
> expect to get many warnings and errors.  Read the messages carefully
> and keep correcting them and rebuilding and checking until they pass.
>
> --
> 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] ggplot2, 'se' variable in geom_errorbar's limits?

2011-02-19 Thread Eric Fail
Hi Scott,

Thank you for taking the time to look at my problem!

I played around with your example and realized that in solving the
problem with limits by summarizing the data I loose the option to
split the data along some third variable, say the 'color' variable in
the diamonds data.

Any idea on how I can solve the problem directly in ggplot2? Any
ggplot2-expects out there?

Sincerely,
Eric


On Sat, Feb 19, 2011 at 4:51 PM, Scott Chamberlain
 wrote:
> require(ggplot2)
>
> data(diamonds)
>
> diamonds <- diamonds[1:100,c(2,7)]
>
> # use ddply in plyr package (loaded with ggplot2) to get data to plot
>
> diamonds_df <- ddply(diamonds, .(cut), summarise,
>
> mean_price = mean(price),
>
> se_price = sd(price)/sqrt(length(price))
>
> )
>
> limits <- aes(ymax = mean_price + se_price, ymin = mean_price - se_price)
>
> ggplot(diamonds_df, aes(x = cut, y = mean_price)) +
>
> geom_point() +
>
> geom_errorbar(limits, width=0.2)
>
> Sincerely,
>
> Scott Chamberlain
>
> Rice University, EEB Dept.
>
> On Saturday, February 19, 2011 at 3:12 PM, Eric Fail wrote:
>
> Can't anybody give me a hint on how to solve this? I even bought the
> ggplot2-book, so you could also give a page (or a series of pages).
>
> Thanks,
> Eric
>
> On Thu, Feb 17, 2011 at 10:19 AM, Eric Fail  wrote:
>
> Dear R-list
>
> I'm working with with geom_errorbar; specifically I'm trying to
> reproduce the example Hadley Wickham have on
> http://had.co.nz/ggplot2/geom_errorbar.html (all in the button of the
> page) where he makes an nice plot with errorbars and then draw lines
> between the points.
>
> What confuses me is the 'limits' he defines for the errorbars from the
> se variable.
>
> First he creates a dataset,
>
> df <- data.frame(
>  trt = factor(c(1, 1, 2, 2)),
>  resp = c(1, 5, 3, 4),
>  group = factor(c(1, 2, 1, 2)),
>  se = c(0.1, 0.3, 0.3, 0.2)
> )
>
> # library(ggplot2)
>
> and then he creates some limits from the se variables.
>
> limits <- aes(ymax = resp + se, ymin=resp - se)
>
> [elements omitted]
>
> # and then he creates the plot (I'm interested in).
>
> p <- ggplot(df, aes(colour=group, y=resp, x=trt))
> p + geom_line(aes(group=group)) + geom_errorbar(limits, width=0.2)
>
> I can (of course) get Hadley's example to run, but I can't do it on my
> data as I don't have a 'se' variable/don't know how to create it. I
> have a group variable, a treatment variable, and a response variable,
> but no se variable.
>
> Could anyone out there explain how I create a 'se' variable in my data?
>
> I'm sure my reasoning is the one that is off, and not ggplot2 (I'm a big
> fan).
>
> Your help is appreciated!
>
> Thanks,
> Eric
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] Seeking help in Package development

2011-02-19 Thread Gabor Grothendieck
On Sat, Feb 19, 2011 at 6:19 PM, Nipesh Bajaj  wrote:
> Hi Gabor, can you be more detail on step 01? Is that not the correct
> path for R.exe? But I found it there.
>
> When I ran R CMD build MyPackage then actually got these warnings:
> cygwin warnings:
> MS-DOS style path detected: c:/Program Files/R/R-2.12.1/bin/Mypackage_1.0.tar
> prefered POSIX equivalent is /cygdrive/c/Program
> Files/R/R-2.12.1/bin/MyPackage_1.0.tar
> CYGWIN environment variable option "nodesfilewarning" turns off this
> warning.
>
> Is this warning something to do with that error?
>
> Sorry for stretching this thread so long. However I am not really
> expert in programming therefore find really hard on understanding
> different terminology given in different documentation.
>
> Thanks,
>

Don't put your package in the R tree.

You can ignore those warnings.  If you had used R.bat then as in my
original instructions you wouldn't have gotten those warnings in the
first place.

-- 
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] ggplot2, 'se' variable in geom_errorbar's limits?

2011-02-19 Thread Ista Zahn
Hi Eric,
Again, it's not really a ggplot question. If you want to plot
separately by color, you need to calculate means and standard errors
separatly by color also:

data(diamonds)

 diamonds_df <- ddply(diamonds, .(cut, color), summarise,
mean_price = mean(price),
se_price = sd(price)/sqrt(length(price))
)

 limits <- aes(ymax = mean_price + se_price, ymin = mean_price - se_price)

 ggplot(diamonds_df, aes(x = cut, y = mean_price, color=color)) +
 geom_point(position=position_dodge(width=0.9)) +
 geom_errorbar(limits, position=position_dodge(width=0.9), width=0.9)

see

?ddply

For details.

Best,
Ista
On Sat, Feb 19, 2011 at 10:58 PM, Eric Fail  wrote:
> Hi Scott,
>
> Thank you for ta king the time to look at my problem!
>
> I played around with your example and realized that in solving the
> problem with limits by summarizing the data I loose the option to
> split the data along some third variable, say the 'color' variable in
> the diamonds data.
>
> Any idea on how I can solve the problem directly in ggplot2? Any
> ggplot2-expects out there?
>
> Sincerely,
> Eric
>
>
> On Sat, Feb 19, 2011 at 4:51 PM, Scott Chamberlain
>  wrote:
>> require(ggplot2)
>>
>> data(diamonds)
>>
>> diamonds <- diamonds[1:100,c(2,7)]
>>
>> # use ddply in plyr package (loaded with ggplot2) to get data to plot
>>
>> diamonds_df <- ddply(diamonds, .(cut), summarise,
>>
>> mean_price = mean(price),
>>
>> se_price = sd(price)/sqrt(length(price))
>>
>> )
>>
>> limits <- aes(ymax = mean_price + se_price, ymin = mean_price - se_price)
>>
>> ggplot(diamonds_df, aes(x = cut, y = mean_price)) +
>>
>> geom_point() +
>>
>> geom_errorbar(limits, width=0.2)
>>
>> Sincerely,
>>
>> Scott Chamberlain
>>
>> Rice University, EEB Dept.
>>
>> On Saturday, February 19, 2011 at 3:12 PM, Eric Fail wrote:
>>
>> Can't anybody give me a hint on how to solve this? I even bought the
>> ggplot2-book, so you could also give a page (or a series of pages).
>>
>> Thanks,
>> Eric
>>
>> On Thu, Feb 17, 2011 at 10:19 AM, Eric Fail  wrote:
>>
>> Dear R-list
>>
>> I'm working with with geom_errorbar; specifically I'm trying to
>> reproduce the example Hadley Wickham have on
>> http://had.co.nz/ggplot2/geom_errorbar.html (all in the button of the
>> page) where he makes an nice plot with errorbars and then draw lines
>> between the points.
>>
>> What confuses me is the 'limits' he defines for the errorbars from the
>> se variable.
>>
>> First he creates a dataset,
>>
>> df <- data.frame(
>>  trt = factor(c(1, 1, 2, 2)),
>>  resp = c(1, 5, 3, 4),
>>  group = factor(c(1, 2, 1, 2)),
>>  se = c(0.1, 0.3, 0.3, 0.2)
>> )
>>
>> # library(ggplot2)
>>
>> and then he creates some limits from the se variables.
>>
>> limits <- aes(ymax = resp + se, ymin=resp - se)
>>
>> [elements omitted]
>>
>> # and then he creates the plot (I'm interested in).
>>
>> p <- ggplot(df, aes(colour=group, y=resp, x=trt))
>> p + geom_line(aes(group=group)) + geom_errorbar(limits, width=0.2)
>>
>> I can (of course) get Hadley's example to run, but I can't do it on my
>> data as I don't have a 'se' variable/don't know how to create it. I
>> have a group variable, a treatment variable, and a response variable,
>> but no se variable.
>>
>> Could anyone out there explain how I create a 'se' variable in my data?
>>
>> I'm sure my reasoning is the one that is off, and not ggplot2 (I'm a big
>> fan).
>>
>> Your help is appreciated!
>>
>> Thanks,
>> Eric
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/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.
>



-- 
Ista Zahn
Graduate student
University of Rochester
Department of Clinical and Social Psychology
http://yourpsyche.org

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


[R] Extreme Values - Help with GPD function

2011-02-19 Thread aleiche

Hi,

I'm a second year Master's student in Applied Statistics. I am doing a  
project using average weekly U.S. regular gasoline prices (in cents,  
per gallon) from an Excel file (from the years 1990- May 2010). I want  
to find the probability that the average weekly U.S. regular gasoline  
prices (in the long term) goes over 400 cents a gallon (or $4.00 a  
gallon). I am using the extRemes program (in R), and I've already  
uploaded my data, and found the probability using the fitted GEV  
distribution (I used the pgev function from the evd package). I want  
to compare the results from the GEV distribution to those of the GPD  
distribution. I've already fitted the GPD distribution to my data (at  
a threshold of 400 cents, where there are 52 entries a year (for most  
years, at least)) and have obtained the supporting graphs.


My question is about how to find the probability that the average  
weekly U.S. regular gasoline prices (in cents, per gallon) goes over  
400 cents a gallon from the Maximum Likelihood Estimates I found when  
I fit the GPD distribution.


MLE's
MLE Std. Err.
Scale (sigma): 62.63489 1.64e-06
Shape (xi): -11.59905 3.948669e-04
Negative log-likelihood: -12.2703516139708

The results did not give me a location (mu) or a beta value. I'm  
wondering which function I can use to find the probability (like when  
I used the pgev). I've already tried to use the pgpd function from the  
evd library, but it kept giving me an error message (I don't remember  
exactly what), and I've tried to use the ppareto and pgenpareto  
function from the actuar package, but they did not work either (I  
think because the shape parameter has to be positive, which mine is  
not). Since the shape paremeter is less than 0, I believe that the GPD  
is of the beta type.


If any of you could help me answer this quandry (I'm a bit of a  
beginner in R) it would help me a lot.


Thank you.


This message was sent using Illinois State University RedbirdMail

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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, 'se' variable in geom_errorbar's limits?

2011-02-19 Thread Eric Fail
Thank you Scott and Ista,

I really appreciate your help! I solved it with Scott's help on the ddply.

For the record, here is the working example that solves my initial question:

## install.packages(c("ggplot2", "plyr"))
require(ggplot2)
data(diamonds)
diamonds <- diamonds[1:100,c(2,7)]

# use ddply in plyr package (loaded with ggplot2) to get data to plot
diamonds_df <- ddply(diamonds, .(cut, color), summarise,
mean_price = mean(price),
se_price = sd(price)/sqrt(length(price))
)

limits <- aes(ymax = mean_price + se_price, ymin = mean_price - se_price)

ggplot(diamonds_df, aes(colour= color, x = cut, y = mean_price)) +
geom_point() +
geom_line(aes(group= color)) +
geom_errorbar(limits, width=0.2)

Very grateful!

Eric


On Sat, Feb 19, 2011 at 6:17 PM, Scott Chamberlain
 wrote:
> Hi Eric,
> I would just include that third variable in the ddply call, for example:
> ddply(diamonds, .(cut, clarity, etc...), summarise,
> mean = mean(price,
> se = ...
> )
> where you can summarise by multiple variables within the ".(x, y, etc.)"\\
> I think that answers your question. Let me know if not. The example I sent
> earlier was just for simplicity.
> Scott
>
> On Saturday, February 19, 2011 at 4:58 PM, Eric Fail wrote:
>
> Hi Scott,
>
> Thank you for taking the time to look at my problem!
>
> I played around with your example and realized that in solving the
> problem with limits by summarizing the data I loose the option to
> split the data along some third variable, say the 'color' variable in
> the diamonds data.
>
> Any idea on how I can solve the problem directly in ggplot2? Any
> ggplot2-expects out there?
>
> Sincerely,
> Eric
>
>
> On Sat, Feb 19, 2011 at 4:51 PM, Scott Chamberlain
>  wrote:
>
> require(ggplot2)
>
> data(diamonds)
>
> diamonds <- diamonds[1:100,c(2,7)]
>
> # use ddply in plyr package (loaded with ggplot2) to get data to plot
>
> diamonds_df <- ddply(diamonds, .(cut), summarise,
>
> mean_price = mean(price),
>
> se_price = sd(price)/sqrt(length(price))
>
> )
>
> limits <- aes(ymax = mean_price + se_price, ymin = mean_price - se_price)
>
> ggplot(diamonds_df, aes(x = cut, y = mean_price)) +
>
> geom_point() +
>
> geom_errorbar(limits, width=0.2)
>
> Sincerely,
>
> Scott Chamberlain
>
> Rice University, EEB Dept.
>
> On Saturday, February 19, 2011 at 3:12 PM, Eric Fail wrote:
>
> Can't anybody give me a hint on how to solve this? I even bought the
> ggplot2-book, so you could also give a page (or a series of pages).
>
> Thanks,
> Eric
>
> On Thu, Feb 17, 2011 at 10:19 AM, Eric Fail  wrote:
>
> Dear R-list
>
> I'm working with with geom_errorbar; specifically I'm trying to
> reproduce the example Hadley Wickham have on
> http://had.co.nz/ggplot2/geom_errorbar.html (all in the button of the
> page) where he makes an nice plot with errorbars and then draw lines
> between the points.
>
> What confuses me is the 'limits' he defines for the errorbars from the
> se variable.
>
> First he creates a dataset,
>
> df <- data.frame(
>  trt = factor(c(1, 1, 2, 2)),
>  resp = c(1, 5, 3, 4),
>  group = factor(c(1, 2, 1, 2)),
>  se = c(0.1, 0.3, 0.3, 0.2)
> )
>
> # library(ggplot2)
>
> and then he creates some limits from the se variables.
>
> limits <- aes(ymax = resp + se, ymin=resp - se)
>
> [elements omitted]
>
> # and then he creates the plot (I'm interested in).
>
> p <- ggplot(df, aes(colour=group, y=resp, x=trt))
> p + geom_line(aes(group=group)) + geom_errorbar(limits, width=0.2)
>
> I can (of course) get Hadley's example to run, but I can't do it on my
> data as I don't have a 'se' variable/don't know how to create it. I
> have a group variable, a treatment variable, and a response variable,
> but no se variable.
>
> Could anyone out there explain how I create a 'se' variable in my data?
>
> I'm sure my reasoning is the one that is off, and not ggplot2 (I'm a big
> fan).
>
> Your help is appreciated!
>
> Thanks,
> Eric
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] Scaling Lattice Graphics for tikzDevice

2011-02-19 Thread Elliot Joel Bernstein
On Feb 18, 2011 11:02 PM, "Deepayan Sarkar" 
wrote:
>
> On Fri, Feb 18, 2011 at 11:04 PM, Elliot Joel Bernstein
>  wrote:
> > I'm trying to use lattice graphics to produce some small plots for
inclusion in a LaTeX file. I want the LaTeX fonts to be used in the plots,
but to be scaled down to match the size of the plot. I have written the
following code to apply a scaling factor to all the "cex" and "padding"
entries in the trellis parameters, but there is still a large white space
between the key and the top of the plot area. Does anyone know how I can get
rid of that (or if I'm going about this all wrong and there's a much cleaner
way)?
> >
> > ## -- ##
> > ## test.R ##
> > ## -- ##
> >
> > require(tikzDevice)
>
> I have not (yet) used tikzDevice, so don't know what impact it will
> have, but in my experience with pdf(), it's usually more effective to
> create a larger file and then scale it when including it in the .tex
> file.
>
> -Deepayan
>
> > require(lattice)
> >
> > rescale.pars <- function(cex.factor = 0.5,
> > padding.factor = 0.25,
> > pars   = trellis.par.get()) {
> >
> >  result <- NULL
> >  for (i in seq_along(pars)) {
> >
> >if (names(pars)[[i]] == "cex") {
> >
> >  if (is.null(result)) {
> >result <- list()
> >  }
> >  result <- c(result, cex=cex.factor*pars[[i]])
> >
> >} else if (grepl("padding$", names(pars)[[i]])) {
> >
> >  if (is.null(result)) {
> >result <- list()
> >  }
> >  eval(parse(text=sprintf("result <- c(result,
%s=padding.factor*pars[[i]])", names(pars)[[i]])))
> >
> >} else if (inherits(pars[[i]], "list")) {
> >
> >  temp <- rescale.pars(cex.factor, padding.factor, pars[[i]])
> >  if (!is.null(temp)) {
> >result[[names(pars)[[i <- temp
> >  }
> >
> >}
> >  }
> >
> >  return (result)
> > }
> >
> > x <- 1:100
> > y <- cumsum(rnorm(length(x), 1, 10))
> > z <- cumsum(rnorm(length(x), 1, 20))
> >
> > tikz("test-plot.tex", width=2, height=2)
> >
> > xyplot(y + z ~ x,
> >   main = "My Plot",
> >   xlab = "x",
> >   auto.key=list(text=c("Line 1", "Line 2"), points=FALSE,
lines=TRUE),
> >   grid = TRUE,
> >   par.settings = rescale.pars())
> >
> > dev.off()
> >
> > ## -- ##
> > ## END R CODE ##
> > ## -- ##
> >
> > Thanks.
> >
> > - Elliot
> >
> > -
> > Elliot Joel Bernstein, Ph.D. | Research Associate | FDO Partners, LLC
> > 134 Mount Auburn Street | Cambridge, MA | 02138
> > Phone: (617) 503-4619 | Email: elliot.bernst...@fdopartners.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.
> >

Deepayan -

Thank you for your reply. In the past I've generally done what you suggest,
but I'm trying to use tikz so the graphics can inherit the fonts used in the
LaTeX document. Is there a way to do that with the pdf device?

Thanks.

- Elliot

[[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] Random Forest & Cross Validation

2011-02-19 Thread ronzhao

Hi,
I am using randomForest package to do some prediction job on GWAS data. I
firstly split the data into training and testing set (70% vs 30%), then
using training set to grow the trees (ntree=10). It looks that the OOB
error in training set is good (<10%). However, it is not very good for the
test set with a AUC only about 50%. 
Although some people said no cross-validation was necessary for RF, I still
felt unsafe and thought a testing set is important. I felt really frustrated
with the results.


Anyone has some suggestions?

Thanks.

PS: example code I used

RF<-randomForest(PHENOTYPE~.,data=Train,importance=T,ntree=2,do.trace=5000)
rownames(Test)<-Test$IID

Pred<-predict(RF,Test,type="prob")


-- 
View this message in context: 
http://r.789695.n4.nabble.com/Random-Forest-Cross-Validation-tp3314777p3314777.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] Partial italic in graph titles when looping

2011-02-19 Thread Josh B
Dear all,

I have a rather complicated problem. I am trying to loop through making graphs, 
so that the graph-making process is fully automated. For each graph, I'd like 
to 
make sure the corresponding title is formatted properly. The titles will be a 
combination of a gene name and numerical position within the gene. The gene 
name 
should be italic-bold, whereas the gene position should be just bold.

Consider the following:

x <- read.table(textConnection("gene position
FLC 3312
TFL1 687
GA1 1127"), header = TRUE, as.is = TRUE)
closeAllConnections()

Now this, below, is essentially how I am automating the graph-making (imagine 
these graphs contain some sort of real data):

par(mfrow = c(3,1))
for (i in 1:nrow(x)){
plot(z <- sort(rnorm(47)), type = "s", main = "")
points(z, cex = .5, col = "dark red")
title(main = paste(x[i,1], " p", x[i,2], sep = ""))
}

The graphs produced by this method are almost perfect, except that the gene 
names are not italicized (they SHOULD be). 


So, once again, the big question is: how would I italicize the gene names but 
NOT the gene positions, when looping through to make these graphs and graph 
titles? If I WASN'T looping to make my graph titles, I could write:

title(main = expression(paste(bolditalic("FLC"), bold("p3312"), sep = " ")))

...but I can't do that, because I'm looping (or can I?)

Thanks in advance for your help!

---
Josh Banta, Ph.D
Center for Genomics and Systems Biology
New York University
100 Washington Square East
New York, NY 10003
Tel: (212) 998-8465
http://plantevolutionaryecology.org


  
[[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] Partial italic in graph titles when looping

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 7:41 PM, Josh B wrote:


Dear all,

I have a rather complicated problem. I am trying to loop through  
making graphs,
so that the graph-making process is fully automated. For each graph,  
I'd like to
make sure the corresponding title is formatted properly. The titles  
will be a
combination of a gene name and numerical position within the gene.  
The gene name

should be italic-bold, whereas the gene position should be just bold.

Consider the following:

x <- read.table(textConnection("gene position
FLC 3312
TFL1 687
GA1 1127"), header = TRUE, as.is = TRUE)
closeAllConnections()

Now this, below, is essentially how I am automating the graph-making  
(imagine

these graphs contain some sort of real data):

par(mfrow = c(3,1))
for (i in 1:nrow(x)){
   plot(z <- sort(rnorm(47)), type = "s", main = "")
   points(z, cex = .5, col = "dark red")
   title(main = paste(x[i,1], " p", x[i,2], sep = ""))
   }


bquote is your friend (at least if one wants to use expressions rather  
than text):


for (i in 1:nrow(x)){
plot(z <- sort(rnorm(47)), type = "s", main = "")
points(z, cex = .5, col = "dark red")
title(main = bquote(italic(.(x[i,1])*" p"*.(x[i,2]
}

--
David.



The graphs produced by this method are almost perfect, except that  
the gene

names are not italicized (they SHOULD be).


So, once again, the big question is: how would I italicize the gene  
names but
NOT the gene positions, when looping through to make these graphs  
and graph

titles? If I WASN'T looping to make my graph titles, I could write:

title(main = expression(paste(bolditalic("FLC"), bold("p3312"), sep  
= " ")))


...but I can't do that, because I'm looping (or can I?)

Thanks in advance for your help!

---
Josh Banta, Ph.D
Center for Genomics and Systems Biology
New York University
100 Washington Square East
New York, NY 10003
Tel: (212) 998-8465
http://plantevolutionaryecology.org



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


Re: [R] Partial italic in graph titles when looping

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 7:41 PM, Josh B wrote:


Dear all,

I have a rather complicated problem. I am trying to loop through  
making graphs,
so that the graph-making process is fully automated. For each graph,  
I'd like to
make sure the corresponding title is formatted properly. The titles  
will be a
combination of a gene name and numerical position within the gene.  
The gene name

should be italic-bold, whereas the gene position should be just bold.

Consider the following:

x <- read.table(textConnection("gene position
FLC 3312
TFL1 687
GA1 1127"), header = TRUE, as.is = TRUE)
closeAllConnections()

Now this, below, is essentially how I am automating the graph-making  
(imagine

these graphs contain some sort of real data):

par(mfrow = c(3,1))
for (i in 1:nrow(x)){
   plot(z <- sort(rnorm(47)), type = "s", main = "")
   points(z, cex = .5, col = "dark red")
   title(main = paste(x[i,1], " p", x[i,2], sep = ""))
   }


Or perhaps (with a shuffling of the parens):
 for (i in 1:nrow(x)){
plot(z <- sort(rnorm(47)), type = "s", main = "")
points(z, cex = .5, col = "dark red")
title(main = bquote(italic(.(x[i,1]))*" p"*.(x[i,2])))
}

The graphs produced by this method are almost perfect, except that  
the gene

names are not italicized (they SHOULD be).


So, once again, the big question is: how would I italicize the gene  
names but
NOT the gene positions, when looping through to make these graphs  
and graph

titles? If I WASN'T looping to make my graph titles, I could write:

title(main = expression(paste(bolditalic("FLC"), bold("p3312"), sep  
= " ")))


...but I can't do that, because I'm looping (or can I?)

Thanks in advance for your help!

---
Josh Banta, Ph.D
Center for Genomics and Systems Biology
New York University
100 Washington Square East
New York, NY 10003
Tel: (212) 998-8465
http://plantevolutionaryecology.org



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


Re: [R] simple recoding problem, but a trouble !

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 10:19 PM, Umesh Rosyara wrote:


Thank you David

I was able to create dataframe and  restore names with the following:

dfr1 <- data.frame(t( apply(dfr, 1, func) ))
names(dfr1) <- c("marker1a","marker1b", "marker2a",  
"marker2b" ,"marker3a", "marker3b")
Still I wonder if there is easier way to restore the names, in  
situations where there are 1000's of variables making the list as  
above might be tidious.


Well, we wouldn't want life to be tidious, now, would we?

> rep(names(dfr), each=2)
[1] "marker1" "marker1" "marker2" "marker2"
> rep(letters[1:2], each=2)
[1] "a" "a" "b" "b"
> paste(rep(names(dfr), each=2), rep(letters[1:2], each=2), sep="")
[1] "marker1a" "marker1a" "marker2b" "marker2b"

--
David.



Thank you for solving my problem. I appreciate it.
Umesh R
From: David Winsemius [mailto:dwinsem...@comcast.net]
Sent: Saturday, February 19, 2011 10:28 AM
To: Umesh Rosyara
Cc: 'Joshua Wiley'; r-help@r-project.org
Subject: Re: [R] simple recoding problem, but a trouble !


On Feb 19, 2011, at 8:40 AM, Umesh Rosyara wrote:

> Just a correction. My expected outdata frame was somehow distorted
> to a
> single, one column. So correct one is:
>
> marker1a   markerb marker2amarker2b
> 1  1   1   1
> 1  3   1   3
> 3  3   3   3
> 3  3   3   3
> 1  3   1   3
> 1  3   1   3


func <- function(x) {sapply( strsplit(x, ""),
 match, table= c("A", NA, "C"))}
t( apply(dfr, 1, func) )

  [,1] [,2] [,3] [,4]
[1,]1111
[2,]1313
[3,]3333
[4,]3333
[5,]1313
[6,]1313


It's amatrix rather than a dataframe and doesn't have colnames but
that should be trivial to fix.

>
> Thanks;
>
> Umesh R
>
>  _
>
> From: Umesh Rosyara [mailto:rosyar...@gmail.com]
> Sent: Friday, February 18, 2011 10:09 PM
> To: 'Joshua Wiley'
> Cc: 'r-help@r-project.org'
> Subject: RE: [R] recoding a data in different way: please help
>
>
> Hi Josh and R community members
>
> Thank you for quick response. I am impressed with the help.
>
> To solve my problems, I tried recode options and I had the following
> problem
> and which motivated me to leave it. Thank you for remind me the  
option

> again, might help to solve my problem in different way.
>
> marker1 <- c("AA", "AC", "CC", "CC", "AC", "AC")
>
> marker2 <- c("AA", "AC", "CC", "CC", "AC", "AC")
>
> dfr <- data.frame(cbind(marker1, marker2))
>
> Objective: replace A with 1, C with 3, and split AA into 1 1 (two
> columns
> numeric). So the intended output for the above dataframe is:
>
>
>
> marker1a
> markerb
> marker2a
> marker2b
>
> 1
> 1
> 1
> 1
>
> 1
> 3
> 1
> 3
>
> 3
> 3
> 3
> 3
>
> 3
> 3
> 3
> 3
>
> 1
> 3
> 1
> 3
>
> 1
> 3
> 1
> 3
>
> I tried the following:
>
> for(i in 1:length(dfr))
>   {
> dfr[[i]]=recode (dfr[[i]],"c('AA')= '1,1'; c('AC')= '1,3';
> c('CA')=
> '1,3';  c('CC')= '3,3' ")
> }
>
> write.table(dfr,"dfr.out", sep=" ,", col.names = T)
> dfn=read.table("dfr.out",header=T, sep="," )
>
> # just trying to cheat R, unfortunately the marker1 and marker  
columns

> remained non-numeric, even when opened in excel !!
>
>
> Unfortunately I got the following result !
>
>   marker1 marker2
> 1 1,1  1,1
> 2 1,2  1,2
> 3 2,2  2,2
> 4 2,2  2,2
> 5 1,2  1,2
> 6 1,2  1,2
>
>
> Sorry to bother all of you, but simple things are being complicated
> these
> days to me.
>
> Thank you so much
> Umesh R
>
>
>  _
>
> From: Joshua Wiley [mailto:jwiley.ps...@gmail.com]
> Sent: Friday, February 18, 2011 12:15 AM
> Cc: r-help@r-project.org
> Subject: Re: [R] recoding a data in different way: please help
>
>
>
> Dear Umesh,
>
> I could not figure out exactly what your recoding scheme was, so I  
do

> not have a specific solution for you.  That said, the following
> functions may help you get started.
>
> ?ifelse # vectorized and different from using if () statements
> ?if #
> ?Logic ## logical operators for your tests
> ## if you install and load the "car" package by John Fox
> ?recode # a function for recoding in package "car"
>
> I am sure it is possible to string together some massive series of  
if
> statements and then use a for loop, but that is probably the  
messiest

> and slowest possible way.  I suspect there will be faster, neater
> options, but I cannot say for certain without having a better feel  
for

> how all the conditions work.
>
> Best regards,
>
> Josh
>
> On Thu, Feb 17, 2011 at 6:21 PM, Umesh Rosyara 
> wrote:
>> Dear R users
>>
>> The following question looks simple but I have spend alot of time
>> to solve
>> it. I would highly appeciate your help.
>>
>> I have following dataset from family dataset :
>>
>> Here we have individuals and their two parents and their marker
>> scores
>> (marker1, marker2,and so on). 0 means that their parent
>> information
> not
>> available.

[R] loop with model fitting pair of consecutive pair of variable: mailing all of you because it was last option

2011-02-19 Thread John Clark
Dear R-users and experts

I want to create to analyse my data which looks like follows:( I have show
only 8 variables but original variables  much more number >1000)


*sub*

*ydata*

*X1*

*X2*

*X3*

*X4*

*X5*

*X6*

*X7*

*X8*

1

12

1

1

1

2

1

1

1

1

2

13

2

2

1

2

2

1

1

1

3

11

1

1

1

2

1

2

1

2

4

12

1

1

2

1

1

2

2

2

5

14

1

2

2

1

1

2

2

2

6

12

2

1

1

2

2

1

1

1

7

8

2

2

2

2

1

1

2

1

8

19

1

1

1

1

1

2

2

1

9

6

2

2

2

1

1

1

1

2

I want to create look to fit models like the following:

lme(ydata ~ X1+ X2, random= intercept)

lme(ydata ~ X3+X4, random= intercept)

lme(ydata ~ X5+ X6, random= intercept)

lme(ydata ~ X7 + X8, random= intercept)

means that the loop should read all the X variables in files, but I do not
want to avoid litting lme(ydata ~ X2+ X3) or lme(ydata ~ X4+X5) or lme(ydata
~ X6+X7). Also I want only output the P value to a vector so that I can plot
the graph with X versus P value.

I have hard time with picking lme output to a vector. If it is impossible
please just suggest in term of linear model (lm), in this type I had no
option than adopt the lm nor lme.

I am sorry to bother you all but it only option remaining !

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] problem installing R in Ubuntu 10.04 -HELP

2011-02-19 Thread fdp
I am running Ubuntu 10.04 (on a Windows machine with dual boot) and I having
a very hard time trying to recover from what I initially thought was a minor
problem. I was trying to install "rattle()" and it failed, and after that I
cannot get R to run AT ALL! I've tried multiple times to reinstall it with a
clean apt-get removal and install and nothing... I keep get the following
message:

*Error in loadNamespace(name) : there is no package called 'rattle'
Fatal error: unable to restore saved data in .RData
*
I've tried to remove all ".RData" files, but couldn't find many... I even
remove all R files/directories/config-files and still nothing worked.

Has anyone encounter something like this? ANYTHING HELPS!!!

Thanks a lot...
fdp

[[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] simple recoding problem, but a trouble !

2011-02-19 Thread Umesh Rosyara
Thank you David 
 
I was able to create dataframe and  restore names with the following:
 
dfr1 <- data.frame(t( apply(dfr, 1, func) ))

names(dfr1) <- c("marker1a","marker1b", "marker2a", "marker2b" ,"marker3a",
"marker3b")
Still I wonder if there is easier way to restore the names, in situations
where there are 1000's of variables making the list as above might be
tidious. 
Thank you for solving my problem. I appreciate it.
Umesh R 
  _  

From: David Winsemius [mailto:dwinsem...@comcast.net] 
Sent: Saturday, February 19, 2011 10:28 AM
To: Umesh Rosyara
Cc: 'Joshua Wiley'; r-help@r-project.org
Subject: Re: [R] simple recoding problem, but a trouble !




On Feb 19, 2011, at 8:40 AM, Umesh Rosyara wrote:

> Just a correction. My expected outdata frame was somehow distorted 
> to a
> single, one column. So correct one is:
>
> marker1a   markerb marker2amarker2b  
> 1  1   1   1 
> 1  3   1   3 
> 3  3   3   3 
> 3  3   3   3 
> 1  3   1   3 
> 1  3   1   3 


func <- function(x) {sapply( strsplit(x, ""),
 match, table= c("A", NA, "C"))}
t( apply(dfr, 1, func) )

  [,1] [,2] [,3] [,4]
[1,]1111
[2,]1313
[3,]3333
[4,]3333
[5,]1313
[6,]1313


It's amatrix rather than a dataframe and doesn't have colnames but 
that should be trivial to fix.

>
> Thanks;
>
> Umesh R
>
>  _
>
> From: Umesh Rosyara [mailto:rosyar...@gmail.com]
> Sent: Friday, February 18, 2011 10:09 PM
> To: 'Joshua Wiley'
> Cc: 'r-help@r-project.org'
> Subject: RE: [R] recoding a data in different way: please help
>
>
> Hi Josh and R community members
>
> Thank you for quick response. I am impressed with the help.
>
> To solve my problems, I tried recode options and I had the following 
> problem
> and which motivated me to leave it. Thank you for remind me the option
> again, might help to solve my problem in different way.
>
> marker1 <- c("AA", "AC", "CC", "CC", "AC", "AC")
>
> marker2 <- c("AA", "AC", "CC", "CC", "AC", "AC")
>
> dfr <- data.frame(cbind(marker1, marker2))
>
> Objective: replace A with 1, C with 3, and split AA into 1 1 (two 
> columns
> numeric). So the intended output for the above dataframe is:
>
>
>
> marker1a
> markerb
> marker2a
> marker2b
>
> 1
> 1
> 1
> 1
>
> 1
> 3
> 1
> 3
>
> 3
> 3
> 3
> 3
>
> 3
> 3
> 3
> 3
>
> 1
> 3
> 1
> 3
>
> 1
> 3
> 1
> 3
>
> I tried the following:
>
> for(i in 1:length(dfr))
>   {
> dfr[[i]]=recode (dfr[[i]],"c('AA')= '1,1'; c('AC')= '1,3'; 
> c('CA')=
> '1,3';  c('CC')= '3,3' ")
> }
>
> write.table(dfr,"dfr.out", sep=" ,", col.names = T)
> dfn=read.table("dfr.out",header=T, sep="," )
>
> # just trying to cheat R, unfortunately the marker1 and marker columns
> remained non-numeric, even when opened in excel !!
>
>
> Unfortunately I got the following result !
>
>   marker1 marker2
> 1 1,1  1,1
> 2 1,2  1,2
> 3 2,2  2,2
> 4 2,2  2,2
> 5 1,2  1,2
> 6 1,2  1,2
>
>
> Sorry to bother all of you, but simple things are being complicated 
> these
> days to me.
>
> Thank you so much
> Umesh R
>
>
>  _
>
> From: Joshua Wiley [mailto:jwiley.ps...@gmail.com]
> Sent: Friday, February 18, 2011 12:15 AM
> Cc: r-help@r-project.org
> Subject: Re: [R] recoding a data in different way: please help
>
>
>
> Dear Umesh,
>
> I could not figure out exactly what your recoding scheme was, so I do
> not have a specific solution for you.  That said, the following
> functions may help you get started.
>
> ?ifelse # vectorized and different from using if () statements
> ?if #
> ?Logic ## logical operators for your tests
> ## if you install and load the "car" package by John Fox
> ?recode # a function for recoding in package "car"
>
> I am sure it is possible to string together some massive series of if
> statements and then use a for loop, but that is probably the messiest
> and slowest possible way.  I suspect there will be faster, neater
> options, but I cannot say for certain without having a better feel for
> how all the conditions work.
>
> Best regards,
>
> Josh
>
> On Thu, Feb 17, 2011 at 6:21 PM, Umesh Rosyara  
> wrote:
>> Dear R users
>>
>> The following question looks simple but I have spend alot of time 
>> to solve
>> it. I would highly appeciate your help.
>>
>> I have following dataset from family dataset :
>>
>> Here we have individuals and their two parents and their marker 
>> scores
>> (marker1, marker2,and so on). 0 means that their parent 
>> information
> not
>> available.
>>
>>
>> Individual  Parent1  Parent2 mark1   mark2
>> 10   0   12  11
>> 20   0   11  22
>> 30   0   13  22
>> 40   0   13  11
>> 51   2   11  12
>> 61   2   12  12
>> 73   4  

[R] Generating uniformly distributed correlated data.

2011-02-19 Thread Søren Faurby
I wish to generate a vector of uniformly distributed data with a  
defined correlation to another vector


The only function I have been able to find doing something similar is  
corgen from the library ecodist.


The following code generates data with the desired correlation to the  
vector x but the resulting vector y is normal and not uniform  
distributed


library(ecodist)
x <- runif(10^5)
y <- corgen(x=x, r=.5)$y

Do anyone know a similar function generating uniform distributed data  
or a way of transforming y to the desired distribution while keeping  
the correlation between x and y


Kind regards, Soren

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] simple recoding problem, but a trouble !

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 10:28 PM, David Winsemius wrote:



On Feb 19, 2011, at 10:19 PM, Umesh Rosyara wrote:


Thank you David

I was able to create dataframe and  restore names with the following:

dfr1 <- data.frame(t( apply(dfr, 1, func) ))
names(dfr1) <- c("marker1a","marker1b", "marker2a",  
"marker2b" ,"marker3a", "marker3b")
Still I wonder if there is easier way to restore the names, in  
situations where there are 1000's of variables making the list as  
above might be tidious.


Well, we wouldn't want life to be tidious, now, would we?

> rep(names(dfr), each=2)
[1] "marker1" "marker1" "marker2" "marker2"
> rep(letters[1:2], each=2)
[1] "a" "a" "b" "b"
> paste(rep(names(dfr), each=2), rep(letters[1:2], each=2), sep="")
[1] "marker1a" "marker1a" "marker2b" "marker2b"


Thanks, Jorge. Should have been one of these:

paste(rep(names(dfr), each=2), letters[1:2], sep="")
# [[1] "marker1a" "marker1b" "marker2a" "marker2b"

paste(rep(names(dfr), each=2), rep(letters[1:2], 2), sep="")





--
David.



Thank you for solving my problem. I appreciate it.
Umesh R
From: David Winsemius [mailto:dwinsem...@comcast.net]
Sent: Saturday, February 19, 2011 10:28 AM
To: Umesh Rosyara
Cc: 'Joshua Wiley'; r-help@r-project.org
Subject: Re: [R] simple recoding problem, but a trouble !


On Feb 19, 2011, at 8:40 AM, Umesh Rosyara wrote:

> Just a correction. My expected outdata frame was somehow distorted
> to a
> single, one column. So correct one is:
>
> marker1a   markerb marker2amarker2b
> 1  1   1   1
> 1  3   1   3
> 3  3   3   3
> 3  3   3   3
> 1  3   1   3
> 1  3   1   3


func <- function(x) {sapply( strsplit(x, ""),
match, table= c("A", NA, "C"))}
t( apply(dfr, 1, func) )

 [,1] [,2] [,3] [,4]
[1,]1111
[2,]1313
[3,]3333
[4,]3333
[5,]1313
[6,]1313


It's amatrix rather than a dataframe and doesn't have colnames but
that should be trivial to fix.

>
> Thanks;
>
> Umesh R
>
>  _
>
> From: Umesh Rosyara [mailto:rosyar...@gmail.com]
> Sent: Friday, February 18, 2011 10:09 PM
> To: 'Joshua Wiley'
> Cc: 'r-help@r-project.org'
> Subject: RE: [R] recoding a data in different way: please help
>
>
> Hi Josh and R community members
>
> Thank you for quick response. I am impressed with the help.
>
> To solve my problems, I tried recode options and I had the  
following

> problem
> and which motivated me to leave it. Thank you for remind me the  
option

> again, might help to solve my problem in different way.
>
> marker1 <- c("AA", "AC", "CC", "CC", "AC", "AC")
>
> marker2 <- c("AA", "AC", "CC", "CC", "AC", "AC")
>
> dfr <- data.frame(cbind(marker1, marker2))
>
> Objective: replace A with 1, C with 3, and split AA into 1 1 (two
> columns
> numeric). So the intended output for the above dataframe is:
>
>
>
> marker1a
> markerb
> marker2a
> marker2b
>
> 1
> 1
> 1
> 1
>
> 1
> 3
> 1
> 3
>
> 3
> 3
> 3
> 3
>
> 3
> 3
> 3
> 3
>
> 1
> 3
> 1
> 3
>
> 1
> 3
> 1
> 3
>
> I tried the following:
>
> for(i in 1:length(dfr))
>   {
> dfr[[i]]=recode (dfr[[i]],"c('AA')= '1,1'; c('AC')= '1,3';
> c('CA')=
> '1,3';  c('CC')= '3,3' ")
> }
>
> write.table(dfr,"dfr.out", sep=" ,", col.names = T)
> dfn=read.table("dfr.out",header=T, sep="," )
>
> # just trying to cheat R, unfortunately the marker1 and marker  
columns

> remained non-numeric, even when opened in excel !!
>
>
> Unfortunately I got the following result !
>
>   marker1 marker2
> 1 1,1  1,1
> 2 1,2  1,2
> 3 2,2  2,2
> 4 2,2  2,2
> 5 1,2  1,2
> 6 1,2  1,2
>
>
> Sorry to bother all of you, but simple things are being complicated
> these
> days to me.
>
> Thank you so much
> Umesh R
>
>
>  _
>
> From: Joshua Wiley [mailto:jwiley.ps...@gmail.com]
> Sent: Friday, February 18, 2011 12:15 AM
> Cc: r-help@r-project.org
> Subject: Re: [R] recoding a data in different way: please help
>
>
>
> Dear Umesh,
>
> I could not figure out exactly what your recoding scheme was, so  
I do

> not have a specific solution for you.  That said, the following
> functions may help you get started.
>
> ?ifelse # vectorized and different from using if () statements
> ?if #
> ?Logic ## logical operators for your tests
> ## if you install and load the "car" package by John Fox
> ?recode # a function for recoding in package "car"
>
> I am sure it is possible to string together some massive series  
of if
> statements and then use a for loop, but that is probably the  
messiest

> and slowest possible way.  I suspect there will be faster, neater
> options, but I cannot say for certain without having a better  
feel for

> how all the conditions work.
>
> Best regards,
>
> Josh
>
> On Thu, Feb 17, 2011 at 6:21 PM, Umesh Rosyara  


> wrote:
>> Dear R users
>>
>> The following question looks simple but I have spend alot of time
>> to s

Re: [R] Generating uniformly distributed correlated data.

2011-02-19 Thread David Winsemius


On Feb 19, 2011, at 9:17 PM, Søren Faurby wrote:

I wish to generate a vector of uniformly distributed data with a  
defined correlation to another vector


The only function I have been able to find doing something similar  
is corgen from the library ecodist.


The following code generates data with the desired correlation to  
the vector x but the resulting vector y is normal and not uniform  
distributed


library(ecodist)
x <- runif(10^5)
y <- corgen(x=x, r=.5)$y

Do anyone know a similar function generating uniform distributed  
data or a way of transforming y to the desired distribution while  
keeping the correlation between x and y


Package "copula" should support that.

(And.) These citations to the archives identified with an RSiteSearch  
search on terms:

uniform multivariate correlation

http://finzi.psych.upenn.edu/Rhelp10/2010-November/258834.html
http://finzi.psych.upenn.edu/R/Rhelp02/archive/57042.html

(Not technically the Archives.)

--

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] Generating uniformly distributed correlated data.

2011-02-19 Thread Peter Langfelder
On Sat, Feb 19, 2011 at 6:17 PM, Søren Faurby
 wrote:
> I wish to generate a vector of uniformly distributed data with a defined
> correlation to another vector
>
> The only function I have been able to find doing something similar is corgen
> from the library ecodist.
>
> The following code generates data with the desired correlation to the vector
> x but the resulting vector y is normal and not uniform distributed
>
> library(ecodist)
> x <- runif(10^5)
> y <- corgen(x=x, r=.5)$y
>
> Do anyone know a similar function generating uniform distributed data or a
> way of transforming y to the desired distribution while keeping the
> correlation between x and y

Hi Soren,

I'm not aware of such functions, but you can try the following code:

# generate some x
n = 100
x = runif(n)
r = 0.5;

y = r * scale(x) + sqrt(1-r^2) * scale(runif(n));

cor(x,y)

The result is not exactly 0.5 because cor(x, runif(n)) is not exactly
zero, but on average you get 0.5 (try to run it 1000 times and
calculate the mean and standard error of the obtained values).

Note that the correlation will be on average r no matter what the
distribution of x is, but the distribution of y will be uniform only
if the distribution of x is uniform.

Peter

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

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


Re: [R] Generating uniformly distributed correlated data.

2011-02-19 Thread Jorge Ivan Velez
Hi Soren,

Take a look at http://tolstoy.newcastle.edu.au/R/help/05/07/7741.html

HTH,
Jorge


On Sat, Feb 19, 2011 at 9:17 PM, Søren Faurby <> wrote:

> I wish to generate a vector of uniformly distributed data with a defined
> correlation to another vector
>
> The only function I have been able to find doing something similar is
> corgen from the library ecodist.
>
> The following code generates data with the desired correlation to the
> vector x but the resulting vector y is normal and not uniform distributed
>
> library(ecodist)
> x <- runif(10^5)
> y <- corgen(x=x, r=.5)$y
>
> Do anyone know a similar function generating uniform distributed data or a
> way of transforming y to the desired distribution while keeping the
> correlation between x and y
>
> Kind regards, Soren
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] Seeking help in Package development

2011-02-19 Thread Joshua Wiley
Hi Nipesh,

Although "c:/Program Files/R/R-2.12.1/bin/" may be where R.exe is
located, it is not where you want to locate your package, and, as
Gabor pointed out, it is best to keep your package outside of the main
R tree.  Your package does not need to be in the same location as
R.exe, you just need the command prompt to know where R.exe is.  If
R.exe is located in location A and your package is in location B, at
the command prompt you can change directories to A, and then tell it
to build the package found in B.  Alternately (again as I mentioned),
you can setyour PATH environment variable to include R's directory,
then it will be found automatically no matter what your current
directory is.

The cygwin warnings are okay in this case.

Unless you urgently need to develop your own package, it might be
easier to learn more about programming in R first and then try
creating a package.  The references I sent in my earlier email would
be good; I also found Software for Data Analysis: Programming with R
by John Chambers to be a very useful resource---I have read most of it
twice now and each time I feel like I have picked up on details or
nuances I missed the first time through.  It is also relatively
inexpensive.  Another helpful book is S Programming by Venables and
Ripley.

You might also consider searching the archives of the R-devel mailing
list.  Although much of it will not be relevant to general package
development, there is a lot you can learn from reading what the people
there discuss.

Cheers,

Josh

On Sat, Feb 19, 2011 at 3:19 PM, Nipesh Bajaj  wrote:
> Hi Gabor, can you be more detail on step 01? Is that not the correct
> path for R.exe? But I found it there.
>
> When I ran R CMD build MyPackage then actually got these warnings:
> cygwin warnings:
> MS-DOS style path detected: c:/Program Files/R/R-2.12.1/bin/Mypackage_1.0.tar
> prefered POSIX equivalent is /cygdrive/c/Program
> Files/R/R-2.12.1/bin/MyPackage_1.0.tar
> CYGWIN environment variable option "nodesfilewarning" turns off this
> warning.
>
> Is this warning something to do with that error?
>
> Sorry for stretching this thread so long. However I am not really
> expert in programming therefore find really hard on understanding
> different terminology given in different documentation.
>
> Thanks,
>
>
> On Sun, Feb 20, 2011 at 4:35 AM, Gabor Grothendieck
>  wrote:
>> On Sat, Feb 19, 2011 at 5:31 PM, Nipesh Bajaj  wrote:
>>> Thanks Gabor for your input. Here what I have done is that:
>>>
>>> 1. Copy 'MyPackage' folder (developed by package.skeleton) into
>>> 'C:\Program Files\R\R-2.12.1\bin' (I found R.exe is there)
>>>
>>> 2. In the command prompt, I changed the working directory using "CD"
>>> command and run 'R CMD build MyPackage'
>>>
>>> 3. I have seen that a file named 'MyPackage_1.0.tar' has been created.
>>> Then I pasted that file in ''C:\Program Files\R\R-2.12.1\bin'
>>>
>>> 4. Again run R CMD INSTALL MyPackage_1.0.tar.gz. However here I got
>>> some error saying:
>>> 'Error in Rd_info(db[[i]]) : Rd files must have a non-empty \title'
>>>
>>
>> What you have done in #1 is asking for trouble.
>>
>> You do need to get used to debugging your Rd files and you should
>> expect to get many warnings and errors.  Read the messages carefully
>> and keep correcting them and rebuilding and checking until they pass.
>>
>> --
>> 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.
>



-- 
Joshua Wiley
Ph.D. Student, Health Psychology
University of California, Los Angeles
http://www.joshuawiley.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] Accessing DF index

2011-02-19 Thread Prof Brian Ripley

On Sat, 19 Feb 2011, David Winsemius wrote:



On Feb 19, 2011, at 4:50 PM, eric wrote:



I have a dataframe called x2. It seems to have a date column but I can't
access it or give it a name or convert it to a date. How do I refer to that
first column and make it a date ? When I try x2[1,] I get the second 
column.


It is not actually a column but is rather a set of row names:

?rownames

rownames(x2)


More precisely, ?row.names for a DF: rownames is for a matrix.  As 
the help you quoted says


 For a data frame, ‘rownames’
 and ‘colnames’ eventually call ‘row.names’ and ‘names’
 respectively, but the latter are preferred.







head(x2)
FAIRXSP500delta
2000-08-31  0.010101096  0.007426964  0.002674132
2000-09-29  0.096679730 -0.054966292  0.151646023
2000-10-31 -0.008245580 -0.004961785 -0.003283795
2000-11-30  0.037024545 -0.083456134  0.120480679
2000-12-29  0.080042708  0.004045193  0.075997514
2001-01-31 -0.009042396  0.034050246 -0.043092643

x[1,1]
  FAIRX
2000-08-31 0.01010110

x[1,2]
   SP500
2000-08-31 0.007426964

str(x2)
'data.frame':   127 obs. of  3 variables:
$ FAIRX: num  0.0101 0.09668 -0.00825 0.03702 0.08004 ...
$ SP500: num  0.00743 -0.05497 -0.00496 -0.08346 0.00405 ...
$ delta: num  0.00267 0.15165 -0.00328 0.12048 0.076 ..
--
View this message in context: 
http://r.789695.n4.nabble.com/Accessing-DF-index-tp3314649p3314649.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.


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.


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