Re: [R] get mean from cdf

2011-09-09 Thread peter dalgaard

On Sep 9, 2011, at 03:51 , Jeff Newmiller wrote:

 The diff function would be helpful.


...or not:

 sum(ppois(0:1000, lambda=2.345, lower=F))
[1] 2.345

(if this was homework: the hard bit is to figure out _why_ this works.)

 annie Zhang annie.zhang2...@gmail.com wrote:
 
 Hi All,
 
 How can I get the expected value from a discrete cdf? Is there any R
 function that can do this?
 
 Thanks,
 Annie

-- 
Peter Dalgaard, Professor,
Center for Statistics, Copenhagen Business School
Solbjerg Plads 3, 2000 Frederiksberg, Denmark
Phone: (+45)38153501
Email: pd@cbs.dk  Priv: pda...@gmail.com
Døden skal tape! --- Nordahl Grieg

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

2011-09-09 Thread Petr PIKAL
Hi

Isn't it something for merge is designed?

 merge(Doctors, DeptCodes, by.x=DocDepts, by.y=Depts)
  DocDeptsDocs  DeptNames
1  Christian\nChristianson  Heart
2    Bob Smith  Brain
3   Greg Jones Anesthesia
4  Al Franklin Anesthesia

It is easy to get rid of the first column.

Regards
Petr


 Re: [R] Factors? I think?
 
 It's probably easiest to think of this as a compound map (doctor - dept
 code - factor - character - integer - dept code - dept name as
 character) and to treat the code as such: if you already have R objects 
with
 the codes in them, it shouldn't be hard to do the transformation.
 
 Consider the following toy set up
 
 Docs = factor(c(Greg Jones,Bob Smith,Al Franklin,Christian
 Christianson))
 DocDepts = factor(c(,,,))
 Doctors = data.frame(Docs, DocDepts)
 
 Depts = factor(1:9 * )
 DeptNames =
 factor(c
 
(Heart,Kidney,Feet,Teeth,Brain,Digestive,Diagnostic,Surgery,Anesthesia))
 DeptCodes = data.frame(Depts,DeptNames)
 # Everything in our data frames is now some sort of factor so we can't 
match
 things up in the normal ways
 
 # Now, you have to do some unpleasantly long but pretty straightforward 
code
 to convert the factors in a way that makes the match properly:
 
 Doctors$numbers - as.numeric(as.character(Doctors[,2])) ## Will extract 
the
  as a real , rather than the internal factor code
 DeptCodes$values - as.numeric(as.character(DeptCodes[,1]))
 
 match(Doctors$numbers, DeptCodes$values) ## Will map the department 
numbers
 onto the correct rows of the DeptCodes df
 
 # Now we get the correct names using those row numbers
 DeptAssignments = as.character(DeptCodes[match(Doctors$numbers,
 DeptCodes$values),2])
 
 # Combine with doctor names to finish
 NamesandTitles = cbind(as.character(Doctors[,1]),DeptAssignments)
 
 It's not the most elegant way of doing it, but hopefully it gives some
 insight into how to work with factors. If you can send a little more
 information about how your data is currently stored we can optimize this
 into something easily repeatable but without specifics, I have to work 
in
 generalities.
 
 Hope this helps,
 
 Michael Weylandt
 
 On Thu, Sep 8, 2011 at 6:36 PM, Totally Inept kramer...@gmail.com 
wrote:
 
  First of all, let me apologize, as this is probably an absurdly basic
  question. I did search before asking, but perhaps my ineptitude didn't
  allow
  me to apply what I read to what I'm doing. Totally new to R, and 
haven't
  done any code in any language in a long time.
 
  Basically I've got categories. They're department codes for doctors 
(say,
   for radiology or  for endocrinology), which of course means 
that
  there are a good number of them, i.e. it's not practical for me to 
write
  them all out as I usually see in examples of categorical variables
  (factors).
 
  And then I've got a list of doctors that I'm actually interested in. I 
have
  the department codes associated with each, but I need to map the 
department
  name to the doctor name. So I might have Greg Jones, Bob Smith, Tom 
Wilson,
  etc... to go with 1234, , , etc.
 
  I need to turn Greg Jones, Bob Smith, ... and 1234, , ... into 
Greg
  Jones, Bob Smith, ... Cardiology, Radiology, 
 
  Obviously I could just search and replace within the csv files but I 
need
  something durable that I can run things through repeatedly.
 
  Anyhow, thanks to anyone willing to humor me with an answer.
 
  --
  View this message in context:
  http://r.789695.n4.nabble.com/Factors-I-think-tp3800413p3800413.html
  Sent from the R help mailing list archive at Nabble.com.
 
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
  http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
 
 
[[alternative HTML version deleted]]
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide 
http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] problems with function read.table

2011-09-09 Thread Petr PIKAL
Hi


 
 Hi,
 
 If you read carefully the help pages for read.table you get this:
 
 
 na.stringsa character vector of strings which are to be interpreted as
 NA../../utils/help/NA values.
 Blank fields are also considered to be missing values in logical, 
integer,
 numeric and complex fields.
 
 So, both NAs and blank fields are considered as NAs directly by 
read.table.
 
 Once you have imported your data, you can modify with any of the string
 manipulation functions (sub() or gsub()) to change your #DIV/0! to the
 string NAs. Another option is to manipulate your Excel file and 
consider
 the division by cero with a IF and get back a NA if that happens.

The only problem is that in such case all columns which has #DIV/0! are 
converted to factors and you need to consider changing it back to numeric.

read.* functions accept as na.string definition not only one value but 
also vector of values and you can get rid of all non numeric and other 
weird Excel values by defining it as a na.strings in read.table call.

 x - read.delim(clipboard)

 str(x)
'data.frame':   6 obs. of  3 variables:
 $ a: int  1 5 9 8 6 3
 $ b: int  3 5 7 0 NA 6
 $ r: Factor w/ 5 levels #DIV/0!,0.3,..: 2 4 5 1 1 3

 y-read.delim(clipboard, na.strings=c(NA, #DIV/0!))
 str(y)
'data.frame':   6 obs. of  3 variables:
 $ a: int  1 5 9 8 6 3
 $ b: int  3 5 7 0 NA 6
 $ r: num  0.333 1 1.286 NA NA ...


Regards
Petr


 
 And finally, instead of using na.omits use option na.rm=T to get done 
your
 calculations:
 
  mean(c(12,23,24,45,67,NA), na.rm=T)[1] 34.2
 
 
 
 Regards,
 Carlos Ortega
 www.qualityexcellence.es
 
 On Thu, Sep 8, 2011 at 4:23 PM, Samir Benzerfa benze...@gmx.ch wrote:
 
  Hello everyone
 
 
 
  I have a couple of questions about the usage of the R function
  read.table(.). My point of departure is that I want to import a 
matrix
  (consisting of time and daily stock returns of many stocks) in R. Most 
of
  the data is numeric, however some values are missing (blanks) and in 
other
  cases I have the character #DIV/0! (from excel). My goal is to do 
some
  regression analysis with this matrix. My questions now are the 
following
  ones:
 
 
 
  1.   How can I in general tell R to automatically replace some 
specific
  numbers or characters in tables by others? (for example to replace all
  characters #DIV/0! by the number 0 or simply NA)
 
  2.   How can I tell R to fill blanks with a number 0 or NA?
 
  3.   How can I tell R to omit the NA fields in the calculations 
but
  not the whole row or column? (I realized that the function na.omit 
omits
  the whole row)
 
 
 
  Many thanks for your help!
 
 
 
  Sincerely,
 
  Samir
 
 
 [[alternative HTML version deleted]]
 
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
  http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
 
 
[[alternative HTML version deleted]]
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide 
http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

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

2011-09-09 Thread Achim Zeileis

On Thu, 8 Sep 2011, Carlos Ortega wrote:


Hi,

Use packages rpart.plot or maptree to enhance the tree drawing.
Or another alternative, use party package that offers much more graphing
possibilities.


There is also the partykit package on R-Forge which offers party-style 
graphics for rpart objects by saying plot(as.party(rpart_object)). See
http://R-Forge.R-project.org/projects/partykit/. A CRAN release is planned 
for later this month.



Regards,
Carlos Ortega
www.qualityexcellence.es


On Thu, Sep 8, 2011 at 5:27 PM, Brian Jensvold brnjns...@gmail.com wrote:


I am trying to create a classification tree using either tree or rpart
functions but when it comes to plotting the results the formatting I get is
different than what I see in all the tutorials (like
http://www.youtube.com/watch?v=9XNhqO1bu0A or
http://www.youtube.com/watch?v=m3mLNpeke0Ifeature=related or
http://www.statmethods.net/advstats/cart.html tree for kyphosis).  I am
trying to take a large demographic population and create a tree which
systematically and accurately divides them into 2 pre-defined
classifications using multiple predictor variables.  What I would like to
see is what I have seen in the tutorials similar to to the ones provided
above, where it shows something like 12/2 in each leaf or stem or step of
the tree which is meant to be interpreded as at this step/stage there are
12 people of this category and 2 of the other.  Instead I am presently
getting something like .342524 which I guess could be the result of
divideing the two groups inorder to find some porportion but im not sure,
and would that be 34% group a or group b?  I was also wondering how you
know
which side is the side of each test where the individual passes or fails?
Is yes/pass always on the left?

I'm pretty sure this is just a mater of changing some sort of default
settings and is in no way a critique of the fine work the designers of R
have freely provided us.


Thank you

   [[alternative HTML version deleted]]

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



[[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] function censReg in panel data setting

2011-09-09 Thread Igors
Dear Arne,

Does censReg support unbalanced panel data? 



/Igors

--
View this message in context: 
http://r.789695.n4.nabble.com/function-censReg-in-panel-data-setting-tp3792227p3801030.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] The elegant way to test if a number is a whole number

2011-09-09 Thread Martin Maechler
 Marc Schwartz marc_schwa...@me.com
 on Thu, 8 Sep 2011 15:41:26 -0500 writes:

 On Sep 8, 2011, at 3:09 PM, David Winsemius wrote:

 
 On Sep 8, 2011, at 3:54 PM, Marc Schwartz wrote:
 
 
 On Sep 8, 2011, at 2:42 PM, David Winsemius wrote:
 
 
 On Sep 8, 2011, at 3:35 PM, Alexander Engelhardt wrote:
 
 Am 08.09.2011 20:48, schrieb Marc Schwartz:
 There was a post from Martin Maechler some years ago and I had to
 search a bit to find it. For these sorts of issues, I typically trust
 his judgement.

thank you, Marc!   :-) 

 The post is here:
 
  https://stat.ethz.ch/pipermail/r-help/2003-April/032471.html
 
 His solution also handles complex numbers.
 
 For those too lazy to follow
 He is basically creating the function is.whole:
 
 is.whole - function(x)
 is.numeric(x)  floor(x)==x
 
 Are you sure? I thought the test would have been all.equal(x,
 round(x,0) )

Yes, I was *very* sure that  the is.whole() above should be
exactly what it was:

The purpose was really different.
There, the discussion started from the fact that  
  as.integer(1) gives FALSE in R
and that some people where making a fuzz about the allegedly
poor design of S and R in that regard.

The purpose of the above  is.whole()  
was solely to be a 
concise  *alternative* to is.integer(),
-- which even would work with complex numbers, as Marc
mentioned...
and BTW, instead of floor(), trunc(), or ceiling() or even round()
would have been equally valid.  I had chosen floor() because  in
mathematics, it is the most frequently used of those four functions.

Your topic below ... about accidental rounding or not with
floating point arithmetic is an entirely different issue,
*and* I agree that a solution to your problem should be based on
the same calculations as all.equal.numeric() .. and hence will
be considerably more involved if it should work in all border
cases.

Martin Maechler, ETH Zurich


 My reasoning was that 1. should be considered 2 but
 floor would make it 1.
 
 David, I am confuzzled. Why would that be equal to 2?
 
 So that sqrt(3) * sqrt(3) would be a whole number. (It is true the
 the floor based wholeness criterion would make sqrt(2)*sqrt(2)
 
 Somehow it doesn't see right that only half of square roots of
 integers that have been squared should pass the wholeness test:
 
  is.whole - function(a, tol = 1e-7) {
 +is.eq - function(x,y) {
 + r - all.equal(x,y, tol=tol)
 + is.logical(r)  r
 +}
 +(is.numeric(a)  is.eq(a, floor(a))) ||
 +(is.complex(a)  {ri - c(Re(a),Im(a)); is.eq(ri, floor(ri))})
 + }
  is.whole( sqrt(2)^2 )
 [1] TRUE
  is.whole( sqrt(3)^2  )
 [1] FALSE
 

 snip content

 OK. I suspect it may down to what assumptions one is willing to make,
 including the level of tolerance used for the comparison.

 is.whole() of course works for 2 because:

 print(sqrt(2) ^ 2, 20)
 [1] 2.0004441

 is slightly larger than 2, so:

 floor(sqrt(2) ^ 2)
 [1] 2

 works, as does:

 round(sqrt(2) ^ 2, 0)
 [1] 2


 On the other hand:

 print(sqrt(3) ^ 2, 20)
 [1] 2.9995559

 is slightly smaller than 3, so:

 floor(sqrt(3) ^ 2)
 [1] 2

 versus:

 round(sqrt(3) ^ 2, 0)
 [1] 3


 Not sure if Martin (cc'd now) might have any comments 8 plus years
 later relative to this issue, as I would again defer to his judgement
 here.

 The other solution proposed, using modulo division, would logically
 fail in both cases:

 (sqrt(3) ^ 2) %% 1 == 0
 [1] FALSE

 (sqrt(2) ^ 2) %% 1 == 0
 [1] FALSE


 Regards,

 Marc

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


[R] prcomp: results with reversed sign in output?

2011-09-09 Thread René Mayer

Dear All,

when I'm running a PCA with

prcomp(USArrests, scale = TRUE)

I get the right principal components, but with the wrong sign infront

Rotation:
PC1 PC2 PC3 PC4
Murder 0.5358995 -0.4181809 0.3412327 0.64922780
Assault 0.5831836 -0.1879856 0.2681484 -0.74340748
UrbanPop 0.2781909 0.8728062 0.3780158 0.13387773
Rape 0.5434321 0.1673186 -0.819 0.08902432

instead of

PC1 PC2 PC3 PC4
Murder -0.5358995 0.4181809 -0.3412327 0.64922780
Assault -0.5831836 0.1879856 -0.2681484 -0.74340748
UrbanPop -0.2781909 -0.8728062 -0.3780158 0.13387773
Rape -0.5434321 -0.1673186 0.819 0.08902432

what is happening here?
any ideas?

thanks,
René

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] prcomp: results with reversed sign in output?

2011-09-09 Thread Paul Hiemstra
 Hi,

If all the signs are switched the PC's are still the same. The principal
vectors are along the same axis, only in a different direction. So there
is no problem :).

hope this helps,
Paul

On 09/09/2011 09:01 AM, René Mayer wrote:
 Dear All,

 when I'm running a PCA with

 prcomp(USArrests, scale = TRUE)

 I get the right principal components, but with the wrong sign infront

 Rotation:
 PC1 PC2 PC3 PC4
 Murder 0.5358995 -0.4181809 0.3412327 0.64922780
 Assault 0.5831836 -0.1879856 0.2681484 -0.74340748
 UrbanPop 0.2781909 0.8728062 0.3780158 0.13387773
 Rape 0.5434321 0.1673186 -0.819 0.08902432

 instead of

 PC1 PC2 PC3 PC4
 Murder -0.5358995 0.4181809 -0.3412327 0.64922780
 Assault -0.5831836 0.1879856 -0.2681484 -0.74340748
 UrbanPop -0.2781909 -0.8728062 -0.3780158 0.13387773
 Rape -0.5434321 -0.1673186 0.819 0.08902432

 what is happening here?
 any ideas?

 thanks,
 René

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


-- 
Paul Hiemstra, Ph.D.
Global Climate Division
Royal Netherlands Meteorological Institute (KNMI)
Wilhelminalaan 10 | 3732 GK | De Bilt | Kamer B 3.39
P.O. Box 201 | 3730 AE | De Bilt
tel: +31 30 2206 494

http://intamap.geo.uu.nl/~paul
http://nl.linkedin.com/pub/paul-hiemstra/20/30b/770

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

2011-09-09 Thread Paul Hiemstra
 Hi,

You could also consider following a statistics course at a university
near you, if you can convince your boss to give you the time. Maybe an
Introductory Statistics course is all you need.

cheers,
Paul

On 09/08/2011 04:49 PM, kensuguro wrote:
 I understand this isn't a r specific question.  I'm switching departments to
 work with the analytics team at my company as a service side manager to
 better incorporate the analytics process into product design / production. 
 We're an online gaming company.

 As I'm going through tools like R, rapidminer, tableau, I was also thinking
 that I should get some formal training in statistics, since my formal
 background is in media arts - it's quite a jump..  I was wondering if there
 was a single go to place for online courses for statistics.  I've checked
 out statistics.com and have found a bunch of online courses from various
 universities, but it's very difficult to find info about how the courses are
 different from each other, which ones aren't worth it, etc..

 I'd greatly appreciate any advice.  I apologize in advance if asking here
 was inappropriate.

 --
 View this message in context: 
 http://r.789695.n4.nabble.com/General-help-online-statistics-courses-tp3799327p3799327.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.


-- 
Paul Hiemstra, Ph.D.
Global Climate Division
Royal Netherlands Meteorological Institute (KNMI)
Wilhelminalaan 10 | 3732 GK | De Bilt | Kamer B 3.39
P.O. Box 201 | 3730 AE | De Bilt
tel: +31 30 2206 494

http://intamap.geo.uu.nl/~paul
http://nl.linkedin.com/pub/paul-hiemstra/20/30b/770

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

2011-09-09 Thread Paul Hiemstra
 And do take a look at Open Courseware:

http://ocw.mit.edu/index.htm

cheers,
Paul

On 09/09/2011 09:30 AM, Paul Hiemstra wrote:
  Hi,

 You could also consider following a statistics course at a university
 near you, if you can convince your boss to give you the time. Maybe an
 Introductory Statistics course is all you need.

 cheers,
 Paul

 On 09/08/2011 04:49 PM, kensuguro wrote:
 I understand this isn't a r specific question.  I'm switching departments to
 work with the analytics team at my company as a service side manager to
 better incorporate the analytics process into product design / production. 
 We're an online gaming company.

 As I'm going through tools like R, rapidminer, tableau, I was also thinking
 that I should get some formal training in statistics, since my formal
 background is in media arts - it's quite a jump..  I was wondering if there
 was a single go to place for online courses for statistics.  I've checked
 out statistics.com and have found a bunch of online courses from various
 universities, but it's very difficult to find info about how the courses are
 different from each other, which ones aren't worth it, etc..

 I'd greatly appreciate any advice.  I apologize in advance if asking here
 was inappropriate.

 --
 View this message in context: 
 http://r.789695.n4.nabble.com/General-help-online-statistics-courses-tp3799327p3799327.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.



-- 
Paul Hiemstra, Ph.D.
Global Climate Division
Royal Netherlands Meteorological Institute (KNMI)
Wilhelminalaan 10 | 3732 GK | De Bilt | Kamer B 3.39
P.O. Box 201 | 3730 AE | De Bilt
tel: +31 30 2206 494

http://intamap.geo.uu.nl/~paul
http://nl.linkedin.com/pub/paul-hiemstra/20/30b/770

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] prcomp: results with reversed sign in output?

2011-09-09 Thread René Mayer

thanks for pointing out Paul,
but the thing which is annoying me in the first place IS this  
direction reversal.

this makes no sense for me
why could this be?

Zitat von Paul Hiemstra paul.hiems...@knmi.nl:


 Hi,

If all the signs are switched the PC's are still the same. The principal
vectors are along the same axis, only in a different direction. So there
is no problem :).

hope this helps,
Paul

On 09/09/2011 09:01 AM, René Mayer wrote:

Dear All,

when I'm running a PCA with

prcomp(USArrests, scale = TRUE)

I get the right principal components, but with the wrong sign infront

Rotation:
PC1 PC2 PC3 PC4
Murder 0.5358995 -0.4181809 0.3412327 0.64922780
Assault 0.5831836 -0.1879856 0.2681484 -0.74340748
UrbanPop 0.2781909 0.8728062 0.3780158 0.13387773
Rape 0.5434321 0.1673186 -0.819 0.08902432

instead of

PC1 PC2 PC3 PC4
Murder -0.5358995 0.4181809 -0.3412327 0.64922780
Assault -0.5831836 0.1879856 -0.2681484 -0.74340748
UrbanPop -0.2781909 -0.8728062 -0.3780158 0.13387773
Rape -0.5434321 -0.1673186 0.819 0.08902432

what is happening here?
any ideas?

thanks,
René

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



--
Paul Hiemstra, Ph.D.
Global Climate Division
Royal Netherlands Meteorological Institute (KNMI)
Wilhelminalaan 10 | 3732 GK | De Bilt | Kamer B 3.39
P.O. Box 201 | 3730 AE | De Bilt
tel: +31 30 2206 494

http://intamap.geo.uu.nl/~paul
http://nl.linkedin.com/pub/paul-hiemstra/20/30b/770




__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] prcomp: results with reversed sign in output?

2011-09-09 Thread Duncan Murdoch

On 11-09-09 5:42 AM, René Mayer wrote:

thanks for pointing out Paul,
but the thing which is annoying me in the first place IS this
direction reversal.
this makes no sense for me
why could this be?


I think you need to read more about principal components.  The signs 
within a PC vector are meaningful only relative to each other (e.g. most 
of the variation in the variables is in the same direction:  as 
standardized murder rates increase, the others increase by approximately 
the same amount), but switching all of the signs within a vector gives 
the same meaning (e.g. as murder decreases, all of the others decrease).


Duncan Murdoch



Zitat von Paul Hiemstrapaul.hiems...@knmi.nl:


  Hi,

If all the signs are switched the PC's are still the same. The principal
vectors are along the same axis, only in a different direction. So there
is no problem :).

hope this helps,
Paul

On 09/09/2011 09:01 AM, René Mayer wrote:

Dear All,

when I'm running a PCA with

prcomp(USArrests, scale = TRUE)

I get the right principal components, but with the wrong sign infront

Rotation:
PC1 PC2 PC3 PC4
Murder 0.5358995 -0.4181809 0.3412327 0.64922780
Assault 0.5831836 -0.1879856 0.2681484 -0.74340748
UrbanPop 0.2781909 0.8728062 0.3780158 0.13387773
Rape 0.5434321 0.1673186 -0.819 0.08902432

instead of

PC1 PC2 PC3 PC4
Murder -0.5358995 0.4181809 -0.3412327 0.64922780
Assault -0.5831836 0.1879856 -0.2681484 -0.74340748
UrbanPop -0.2781909 -0.8728062 -0.3780158 0.13387773
Rape -0.5434321 -0.1673186 0.819 0.08902432

what is happening here?
any ideas?

thanks,
René

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



--
Paul Hiemstra, Ph.D.
Global Climate Division
Royal Netherlands Meteorological Institute (KNMI)
Wilhelminalaan 10 | 3732 GK | De Bilt | Kamer B 3.39
P.O. Box 201 | 3730 AE | De Bilt
tel: +31 30 2206 494

http://intamap.geo.uu.nl/~paul
http://nl.linkedin.com/pub/paul-hiemstra/20/30b/770




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

2011-09-09 Thread Göran Broström
I have a data frame 'tmp' and a vector 'name' containing 'd2'.
I want to save 'tmp' under the name hidden in 'name', and the file must have
the same name, plus the extension '.rda'.
So I try

 tmp
 x y
1 1 3
2 2 4
 name
[1] d2
 assign(name, tmp)
 summary(get(name))
  x  y
Min.   :1.00   Min.   :3.00
1st Qu.:1.25   1st Qu.:3.25
Median :1.50   Median :3.50
Mean   :1.50   Mean   :3.50
3rd Qu.:1.75   3rd Qu.:3.75
Max.   :2.00   Max.   :4.00
 save(get(name), file = paste(name, rda, sep = .))
Error in save(get(name), file = paste(name, rda, sep = .)) :
 object ‘get(name)’ not found
++
I can't figure out where I 'get' it wrong.

In real life I have a bunch of text files named like 'd1.txt', ...,
'd100.txt'.
I want to convert all of them to R data files, with one data frame in each
named d1, ..., d100.

Any suggestion is much appreciated!
-- 
Göran Broström

[[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] Factors? I think?

2011-09-09 Thread peter dalgaard

On Sep 9, 2011, at 09:13 , Petr PIKAL wrote:

 Hi
 
 Isn't it something for merge is designed?

Sort of. (You'd need to think carefully about what happens with non-matched 
codes.)

Wouldn't this do the trick as well?

in - as.character(DeptCodes$DeptCodes)
out - as.character(DeptCodes$DeptNames)
Doctors - within(Doctors, DeptNames - factor(DocDepts, levels=in, labels=out))

 
 merge(Doctors, DeptCodes, by.x=DocDepts, by.y=Depts)
  DocDeptsDocs  DeptNames
 1  Christian\nChristianson  Heart
 2    Bob Smith  Brain
 3   Greg Jones Anesthesia
 4  Al Franklin Anesthesia
 
 It is easy to get rid of the first column.
 
 Regards
 Petr
 
 
 Re: [R] Factors? I think?
 
 It's probably easiest to think of this as a compound map (doctor - dept
 code - factor - character - integer - dept code - dept name as
 character) and to treat the code as such: if you already have R objects 
 with
 the codes in them, it shouldn't be hard to do the transformation.
 
 Consider the following toy set up
 
 Docs = factor(c(Greg Jones,Bob Smith,Al Franklin,Christian
 Christianson))
 DocDepts = factor(c(,,,))
 Doctors = data.frame(Docs, DocDepts)
 
 Depts = factor(1:9 * )
 DeptNames =
 factor(c
 
 (Heart,Kidney,Feet,Teeth,Brain,Digestive,Diagnostic,Surgery,Anesthesia))
 DeptCodes = data.frame(Depts,DeptNames)
 # Everything in our data frames is now some sort of factor so we can't 
 match
 things up in the normal ways
 
 # Now, you have to do some unpleasantly long but pretty straightforward 
 code
 to convert the factors in a way that makes the match properly:
 
 Doctors$numbers - as.numeric(as.character(Doctors[,2])) ## Will extract 
 the
  as a real , rather than the internal factor code
 DeptCodes$values - as.numeric(as.character(DeptCodes[,1]))
 
 match(Doctors$numbers, DeptCodes$values) ## Will map the department 
 numbers
 onto the correct rows of the DeptCodes df
 
 # Now we get the correct names using those row numbers
 DeptAssignments = as.character(DeptCodes[match(Doctors$numbers,
 DeptCodes$values),2])
 
 # Combine with doctor names to finish
 NamesandTitles = cbind(as.character(Doctors[,1]),DeptAssignments)
 
 It's not the most elegant way of doing it, but hopefully it gives some
 insight into how to work with factors. If you can send a little more
 information about how your data is currently stored we can optimize this
 into something easily repeatable but without specifics, I have to work 
 in
 generalities.
 
 Hope this helps,
 
 Michael Weylandt
 
 On Thu, Sep 8, 2011 at 6:36 PM, Totally Inept kramer...@gmail.com 
 wrote:
 
 First of all, let me apologize, as this is probably an absurdly basic
 question. I did search before asking, but perhaps my ineptitude didn't
 allow
 me to apply what I read to what I'm doing. Totally new to R, and 
 haven't
 done any code in any language in a long time.
 
 Basically I've got categories. They're department codes for doctors 
 (say,
  for radiology or  for endocrinology), which of course means 
 that
 there are a good number of them, i.e. it's not practical for me to 
 write
 them all out as I usually see in examples of categorical variables
 (factors).
 
 And then I've got a list of doctors that I'm actually interested in. I 
 have
 the department codes associated with each, but I need to map the 
 department
 name to the doctor name. So I might have Greg Jones, Bob Smith, Tom 
 Wilson,
 etc... to go with 1234, , , etc.
 
 I need to turn Greg Jones, Bob Smith, ... and 1234, , ... into 
 Greg
 Jones, Bob Smith, ... Cardiology, Radiology, 
 
 Obviously I could just search and replace within the csv files but I 
 need
 something durable that I can run things through repeatedly.
 
 Anyhow, thanks to anyone willing to humor me with an answer.
 
 --
 View this message in context:
 http://r.789695.n4.nabble.com/Factors-I-think-tp3800413p3800413.html
 Sent from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 
   [[alternative HTML version deleted]]
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide 
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

-- 
Peter Dalgaard, Professor
Center for Statistics, Copenhagen Business 

Re: [R] get and save

2011-09-09 Thread Ivan Calandra

Hi!

Why don't you just save tmp?
save(tmp, file = paste(name, rda, sep = .))

I don't think it makes much difference, at least not with your example. 
Or does it?


Ivan


Le 9/9/2011 11:55, Göran Broström a écrit :

I have a data frame 'tmp' and a vector 'name' containing 'd2'.
I want to save 'tmp' under the name hidden in 'name', and the file must have
the same name, plus the extension '.rda'.
So I try


tmp

  x y
1 1 3
2 2 4

name

[1] d2

assign(name, tmp)
summary(get(name))

   x  y
Min.   :1.00   Min.   :3.00
1st Qu.:1.25   1st Qu.:3.25
Median :1.50   Median :3.50
Mean   :1.50   Mean   :3.50
3rd Qu.:1.75   3rd Qu.:3.75
Max.   :2.00   Max.   :4.00

save(get(name), file = paste(name, rda, sep = .))

Error in save(get(name), file = paste(name, rda, sep = .)) :
  object ‘get(name)’ not found
++
I can't figure out where I 'get' it wrong.

In real life I have a bunch of text files named like 'd1.txt', ...,
'd100.txt'.
I want to convert all of them to R data files, with one data frame in each
named d1, ..., d100.

Any suggestion is much appreciated!


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


--
Ivan CALANDRA
PhD Student
University of Hamburg
Biozentrum Grindel und Zoologisches Museum
Dept. Mammalogy
Martin-Luther-King-Platz 3
D-20146 Hamburg, GERMANY
+49(0)40 42838 6231
ivan.calan...@uni-hamburg.de

**
http://www.for771.uni-bonn.de
http://webapp5.rrz.uni-hamburg.de/mammals/eng/1525_8_1.php

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] prcomp: results with reversed sign in output?

2011-09-09 Thread Ted Harding
The point is that a principal component vector is a solution,
say V, of a matrix equation A%*%V = L*V where A is the matrix
and L is a scalar..

Since this equation can be written A%*%(-V) = L*(-V), the
result is indeterminate with respect to its sign. If V is a
solution, so is (-V), and vice versa. It is not a case of
direction reversal, since neither L nor (-L) has a primary
role -- they are equivalent, and you can adopt either one.
Just make it clear which one you adopt -- or someone else
might think that they disagree!

You originally wrote I get the right principal components,
but with the wrong sign in front. You did not get the wrong
sign -- both are correct! It may be that you are comparing
your result from R with the result from some other software
(or from a textbook, or whatever) which produced the equivalent
result but with the opposite sign. Again, both are correct.

If, for some reason, you do not like the sign of the result
you get, then change its sign.

Hoping this helps,
Ted.

On 09-Sep-11 09:42:49, René Mayer wrote:
 thanks for pointing out Paul,
 but the thing which is annoying me in the first place IS this  
 direction reversal.
 this makes no sense for me
 why could this be?
 
 Zitat von Paul Hiemstra paul.hiems...@knmi.nl:
 
  Hi,

 If all the signs are switched the PC's are still the same. The
 principal
 vectors are along the same axis, only in a different direction. So
 there
 is no problem :).

 hope this helps,
 Paul

 On 09/09/2011 09:01 AM, René Mayer wrote:
 Dear All,

 when I'm running a PCA with

 prcomp(USArrests, scale = TRUE)

 I get the right principal components, but with the wrong sign infront

 Rotation:
 PC1 PC2 PC3 PC4
 Murder 0.5358995 -0.4181809 0.3412327 0.64922780
 Assault 0.5831836 -0.1879856 0.2681484 -0.74340748
 UrbanPop 0.2781909 0.8728062 0.3780158 0.13387773
 Rape 0.5434321 0.1673186 -0.819 0.08902432

 instead of

 PC1 PC2 PC3 PC4
 Murder -0.5358995 0.4181809 -0.3412327 0.64922780
 Assault -0.5831836 0.1879856 -0.2681484 -0.74340748
 UrbanPop -0.2781909 -0.8728062 -0.3780158 0.13387773
 Rape -0.5434321 -0.1673186 0.819 0.08902432

 what is happening here?
 any ideas?

 thanks,
 René

 --
 Paul Hiemstra, Ph.D.
 Global Climate Division
 Royal Netherlands Meteorological Institute (KNMI)
 Wilhelminalaan 10 | 3732 GK | De Bilt | Kamer B 3.39
 P.O. Box 201 | 3730 AE | De Bilt
 tel: +31 30 2206 494

 http://intamap.geo.uu.nl/~paul
 http://nl.linkedin.com/pub/paul-hiemstra/20/30b/770


E-Mail: (Ted Harding) ted.hard...@wlandres.net
Fax-to-email: +44 (0)870 094 0861
Date: 09-Sep-11   Time: 11:04:16
-- 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.


Re: [R] get and save

2011-09-09 Thread Duncan Murdoch

On 11-09-09 5:55 AM, Göran Broström wrote:

I have a data frame 'tmp' and a vector 'name' containing 'd2'.
I want to save 'tmp' under the name hidden in 'name', and the file must have
the same name, plus the extension '.rda'.
So I try


tmp

  x y
1 1 3
2 2 4

name

[1] d2

assign(name, tmp)
summary(get(name))

   x  y
Min.   :1.00   Min.   :3.00
1st Qu.:1.25   1st Qu.:3.25
Median :1.50   Median :3.50
Mean   :1.50   Mean   :3.50
3rd Qu.:1.75   3rd Qu.:3.75
Max.   :2.00   Max.   :4.00

save(get(name), file = paste(name, rda, sep = .))

Error in save(get(name), file = paste(name, rda, sep = .)) :
  object ‘get(name)’ not found
++
I can't figure out where I 'get' it wrong.


save() does funny things with its arguments:  the first one is not 
evaluated in the standard way.  It wants a name or a character string. 
So you would want to pass name instead of get(name), but that would end 
up being interpreted the wrong way and you'd just save the name 
variable, not the thing it refers to.


You will probably have to use do.call to do this, i.e.

do.call(save, list(name, file=paste(name, rda, sep = .)))

This will evaluate the name before passing it to save().

Duncan Murdoch




In real life I have a bunch of text files named like 'd1.txt', ...,
'd100.txt'.
I want to convert all of them to R data files, with one data frame in each
named d1, ..., d100.

Any suggestion is much appreciated!



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

2011-09-09 Thread Duncan Murdoch

On 11-09-09 6:01 AM, Ivan Calandra wrote:

Hi!

Why don't you just save tmp?
save(tmp, file = paste(name, rda, sep = .))

I don't think it makes much difference, at least not with your example.
Or does it?


That saves it with the name tmp, not d2 as Göran wanted.

Duncan Murdoch



Ivan


Le 9/9/2011 11:55, Göran Broström a écrit :

I have a data frame 'tmp' and a vector 'name' containing 'd2'.
I want to save 'tmp' under the name hidden in 'name', and the file must have
the same name, plus the extension '.rda'.
So I try


tmp

   x y
1 1 3
2 2 4

name

[1] d2

assign(name, tmp)
summary(get(name))

x  y
Min.   :1.00   Min.   :3.00
1st Qu.:1.25   1st Qu.:3.25
Median :1.50   Median :3.50
Mean   :1.50   Mean   :3.50
3rd Qu.:1.75   3rd Qu.:3.75
Max.   :2.00   Max.   :4.00

save(get(name), file = paste(name, rda, sep = .))

Error in save(get(name), file = paste(name, rda, sep = .)) :
   object ‘get(name)’ not found
++
I can't figure out where I 'get' it wrong.

In real life I have a bunch of text files named like 'd1.txt', ...,
'd100.txt'.
I want to convert all of them to R data files, with one data frame in each
named d1, ..., d100.

Any suggestion is much appreciated!


__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] rgl: axis/viewport/box problems in persp3d()

2011-09-09 Thread Marius Hofert
Dear expeRts,

I am a new user of rgl, below is my first trial to plot a simple function in 
3d. 
I managed to put the axes in the right locations, but:
(1) The xlab, ylab, and zlab arguments are ignored; how can I put in axes 
labels?
(2) Since I removed the axes in persp3d() the viewport is too small; is it 
possible
to keep the size of the viewport?
(3) The box is not correctly drawn, there are two holes, one in (0,0,1) and 
one
in (1,1,0); how can I fix that?

Cheers,

Marius


require(rgl)
s - seq(0, 1, length.out=21)
M - function(u) apply(u, 1, min)
u - s
v - s
z - outer(u, v, function(u,v) M(cbind(u,v)))
persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=, 
ylab=, zlab=) 
axes3d(edges=c('x--','y--','z+-'), xlab=x, ylab=y, zlab=z)
par3d(windowRect=c(0,0,480,480))

R1 - rotationMatrix(-55*pi/180, 1,0,0) 
R3 - rotationMatrix(50*pi/180, 0,0,1)
R - R1 %*% R3
rgl.viewpoint(interactive=TRUE, userMatrix=R) # rotate
rgl.postscript(myplot.pdf, fmt=pdf)
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] NMDS plot and Adonis (PerMANOVA) of community composition with presence absence and relative intensity

2011-09-09 Thread Markus Lindh
Hi!

Thanks for providing great help in R-related statistics. Now, however I'm
stuck. I'm not a statistics person but I was recommended to use R to perform
a nmds plot and PerMANOVA of my dataset.

Sample(treatment) in the columns and species (OTU) in the rows. I have 4
treatments (Ambient Temperature, Ambient temperature+Low pH, High
temperature, High temperature+low pH), and I have 16 different species
(OTUs).

Firstly I want to make a nmds plot and see how different the treatments are.
Meaning how the species covaries with the environment, I can produce a plot
but don't understand how to interpret it and I get an error message. Then I
want to make a PerMANOVA analysis comparing the different treatments with
each other.

I tried to look at similar problems but I get error messages.

Any ideas?

#NMDS
step1-read.delim2(day20.txt, row.names=1)
library(clusterSim)
step2-data.Normalization(step1,type=n10)
step3-asin(sqrt(step2))*57.3
step4-t(step3)
library(vegan)
step5-data.matrix(vegdist(step4,method=bray))
step6-metaMDS(step5, autotransform=FALSE)
plot(step5)
Warning message:
In ordiplot(x, choices = choices, type = type, display = display,  :
  Species scores not available


#PerMANOVA
step1-read.delim2(day20.txt, row.names=1)
library(clusterSim)
step2-data.Normalization(step1,type=n10)
step3-asin(sqrt(step2))*57.3
step4-t(step3)
library(vegan)
step5-data.matrix(vegdist(step4,method=bray))
attach(step5)
step6-adonis(step5~Sample1+Sample2+Sample3,data=step5)

[[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] General help - online statistics courses?

2011-09-09 Thread Michael Griffiths
Hi Kensuguro,

I would try the Royal Statistical Society, based in London. They have a
range of courses from relatively basic through to graduate level, all of
which are based on distance learning. However, what they also have now are
modularised programmes, which might be ideal for you, as you are able to
select those areas in which you are particularly interested.

Good Luck!

Mike Griffiths

On Thu, Sep 8, 2011 at 5:49 PM, kensuguro magronb...@gmail.com wrote:

 I understand this isn't a r specific question.  I'm switching departments
 to
 work with the analytics team at my company as a service side manager to
 better incorporate the analytics process into product design / production.
 We're an online gaming company.

 As I'm going through tools like R, rapidminer, tableau, I was also thinking
 that I should get some formal training in statistics, since my formal
 background is in media arts - it's quite a jump..  I was wondering if there
 was a single go to place for online courses for statistics.  I've checked
 out statistics.com and have found a bunch of online courses from various
 universities, but it's very difficult to find info about how the courses
 are
 different from each other, which ones aren't worth it, etc..

 I'd greatly appreciate any advice.  I apologize in advance if asking here
 was inappropriate.

 --
 View this message in context:
 http://r.789695.n4.nabble.com/General-help-online-statistics-courses-tp3799327p3799327.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.




-- 

*Michael Griffiths, Ph.D
*Statistician

*Upstream Systems
*23 Old Bond Street, London, W1S 4PZ
Mob +44 750 011 2873
DL   +44 (0) 20 7869 5147
Fax  +44 207 290 1321

*griffi...@upstreamsystems.com einst...@upstreamsystems.com*
www.upstreamsystems.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] rgl: axis/viewport/box problems in persp3d()

2011-09-09 Thread Duncan Murdoch

On 11-09-09 6:18 AM, Marius Hofert wrote:

Dear expeRts,

I am a new user of rgl, below is my first trial to plot a simple function in 3d.
I managed to put the axes in the right locations, but:
(1) The xlab, ylab, and zlab arguments are ignored; how can I put in axes 
labels?


Those are documented on the axes3d page, but are arguments to title3d, 
not axes3d.  So add title3d(xlab=x, etc.



(2) Since I removed the axes in persp3d() the viewport is too small; is it 
possible
 to keep the size of the viewport?


You can manually adjust it to your taste, then write down the value of 
par3d(zoom).  Later you can reproduce the resizing by calling 
par3d(zoom= saved value ).




(3) The box is not correctly drawn, there are two holes, one in (0,0,1) and 
one
 in (1,1,0); how can I fix that?


That happens because OpenGL has a limit on the range of depths that can 
be displayed, and the corners of the box have been adjusted to be too 
close or far.  This is arguably a bug in rgl, but it's sometimes a feature.


What I'd suggest is that you don't use rgl.viewpoint, you just manually 
adjust the display as you like, without making it quite as extreme, then 
record the values of par3d(c(userMatrix, zoom, FOV)); those 
control the viewpoint.


Duncan Murdoch



Cheers,

Marius


require(rgl)
s- seq(0, 1, length.out=21)
M- function(u) apply(u, 1, min)
u- s
v- s
z- outer(u, v, function(u,v) M(cbind(u,v)))
persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=,
 ylab=, zlab=)
axes3d(edges=c('x--','y--','z+-'), xlab=x, ylab=y, zlab=z)
par3d(windowRect=c(0,0,480,480))

R1- rotationMatrix(-55*pi/180, 1,0,0)
R3- rotationMatrix(50*pi/180, 0,0,1)
R- R1 %*% R3
rgl.viewpoint(interactive=TRUE, userMatrix=R) # rotate
rgl.postscript(myplot.pdf, fmt=pdf)
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] prcomp: results with reversed sign in output?

2011-09-09 Thread René Mayer

thanks for explaining Duncan and Ted,
Indeed, I did compare my results from a textbook
and noticed that I consitenly get flipped signs and biplots.
regards
René



Zitat von ted.hard...@wlandres.net:


The point is that a principal component vector is a solution,
say V, of a matrix equation A%*%V = L*V where A is the matrix
and L is a scalar..

Since this equation can be written A%*%(-V) = L*(-V), the
result is indeterminate with respect to its sign. If V is a
solution, so is (-V), and vice versa. It is not a case of
direction reversal, since neither L nor (-L) has a primary
role -- they are equivalent, and you can adopt either one.
Just make it clear which one you adopt -- or someone else
might think that they disagree!

You originally wrote I get the right principal components,
but with the wrong sign in front. You did not get the wrong
sign -- both are correct! It may be that you are comparing
your result from R with the result from some other software
(or from a textbook, or whatever) which produced the equivalent
result but with the opposite sign. Again, both are correct.

If, for some reason, you do not like the sign of the result
you get, then change its sign.

Hoping this helps,
Ted.

On 09-Sep-11 09:42:49, René Mayer wrote:

thanks for pointing out Paul,
but the thing which is annoying me in the first place IS this
direction reversal.
this makes no sense for me
why could this be?

Zitat von Paul Hiemstra paul.hiems...@knmi.nl:


 Hi,

If all the signs are switched the PC's are still the same. The
principal
vectors are along the same axis, only in a different direction. So
there
is no problem :).

hope this helps,
Paul

On 09/09/2011 09:01 AM, René Mayer wrote:

Dear All,

when I'm running a PCA with

prcomp(USArrests, scale = TRUE)

I get the right principal components, but with the wrong sign infront

Rotation:
PC1 PC2 PC3 PC4
Murder 0.5358995 -0.4181809 0.3412327 0.64922780
Assault 0.5831836 -0.1879856 0.2681484 -0.74340748
UrbanPop 0.2781909 0.8728062 0.3780158 0.13387773
Rape 0.5434321 0.1673186 -0.819 0.08902432

instead of

PC1 PC2 PC3 PC4
Murder -0.5358995 0.4181809 -0.3412327 0.64922780
Assault -0.5831836 0.1879856 -0.2681484 -0.74340748
UrbanPop -0.2781909 -0.8728062 -0.3780158 0.13387773
Rape -0.5434321 -0.1673186 0.819 0.08902432

what is happening here?
any ideas?

thanks,
René


--
Paul Hiemstra, Ph.D.
Global Climate Division
Royal Netherlands Meteorological Institute (KNMI)
Wilhelminalaan 10 | 3732 GK | De Bilt | Kamer B 3.39
P.O. Box 201 | 3730 AE | De Bilt
tel: +31 30 2206 494

http://intamap.geo.uu.nl/~paul
http://nl.linkedin.com/pub/paul-hiemstra/20/30b/770



E-Mail: (Ted Harding) ted.hard...@wlandres.net
Fax-to-email: +44 (0)870 094 0861
Date: 09-Sep-11   Time: 11:04:16
-- 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.


Re: [R] get and save

2011-09-09 Thread Göran Broström
On Fri, Sep 9, 2011 at 12:09 PM, Duncan Murdoch murdoch.dun...@gmail.comwrote:

 On 11-09-09 6:01 AM, Ivan Calandra wrote:

 Hi!

 Why don't you just save tmp?
 save(tmp, file = paste(name, rda, sep = .))

 I don't think it makes much difference, at least not with your example.
 Or does it?


 That saves it with the name tmp, not d2 as Göran wanted.

 Exactly; that was one of my 'failures'. Thank you both for the help!

Göran


 Duncan Murdoch



 Ivan


 Le 9/9/2011 11:55, Göran Broström a écrit :

 I have a data frame 'tmp' and a vector 'name' containing 'd2'.
 I want to save 'tmp' under the name hidden in 'name', and the file must
 have
 the same name, plus the extension '.rda'.
 So I try

  tmp

   x y
 1 1 3
 2 2 4

 name

 [1] d2

 assign(name, tmp)
 summary(get(name))

x  y
 Min.   :1.00   Min.   :3.00
 1st Qu.:1.25   1st Qu.:3.25
 Median :1.50   Median :3.50
 Mean   :1.50   Mean   :3.50
 3rd Qu.:1.75   3rd Qu.:3.75
 Max.   :2.00   Max.   :4.00

 save(get(name), file = paste(name, rda, sep = .))

 Error in save(get(name), file = paste(name, rda, sep = .)) :
   object ‘get(name)’ not found
 ++
 I can't figure out where I 'get' it wrong.

 In real life I have a bunch of text files named like 'd1.txt', ...,
 'd100.txt'.
 I want to convert all of them to R data files, with one data frame in
 each
 named d1, ..., d100.

 Any suggestion is much appreciated!


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




-- 
Göran Broström

[[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] Spider (Radar) Plot

2011-09-09 Thread Jim Lemon

On 09/09/2011 06:31 AM, XINLI LI wrote:

Dear R Group:

Based on the following data, how to do a great Spider (Radar) Plot? Any
advice is greatly appreciated.

  HospID Rate Age Charlson NIHSS 1 0.2 49 3.5 0 2 0.1 48 1.8 12 3 0.4 56
2.1 5 4 0.3 77 0 7 5 0.2 67 6.5 3 6 0.1 62 4.8 4.6 7 0.1 64 12 5.2 8 0.3 61
3 2.8 9 0.15 69 4.5 1.9 10 0.22 80 0 6.7 11 0.34 61 6 4.2 12 0.18 63 3 6.1
13 0.09 64 8 15 12 0 56 10 11 15 0.1 70 11 7 16 0.21 43 9 9


Hi XINLI,

How about this:

hosp-read.table(hosp.dat,header=TRUE)
library(plotrix)
radial.plot(t(hosp[,4:5]),rp.type=p,labels=hosp$HospID,
 show.radial.grid=FALSE)
par(xpd=TRUE)
legend(10,17,c(NIHSS,Charlson),lty=1,col=c(2,1))
par(xpd=FALSE)

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] NMDS plot and Adonis (PerMANOVA) of community composition with presence absence and relative intensity

2011-09-09 Thread Jari Oksanen
Markus Lindh markusvlindh at gmail.com writes:

 
 Hi!
 
 Thanks for providing great help in R-related statistics. Now, however I'm
 stuck. I'm not a statistics person but I was recommended to use R to perform
 a nmds plot and PerMANOVA of my dataset.
 
 Sample(treatment) in the columns and species (OTU) in the rows. I have 4
 treatments (Ambient Temperature, Ambient temperature+Low pH, High
 temperature, High temperature+low pH), and I have 16 different species
 (OTUs).
 
 Firstly I want to make a nmds plot and see how different the treatments are.
 Meaning how the species covaries with the environment, I can produce a plot
 but don't understand how to interpret it and I get an error message. Then I
 want to make a PerMANOVA analysis comparing the different treatments with
 each other.
 
 I tried to look at similar problems but I get error messages.
 
 Any ideas?

What kind of error messages do you get? I don't see any in your message.
You get one Warning message, but that is not an error (difference 
between warnings and errors are, among other things, that warnings are
called warnings, and errors are called errors). Of course, I cannot
reproduce your analys because I have no clue what is in day20.txt. 
A reproducible example could be useful. Here some comments, though.

1) do not use data.matrix() for the the result of dist() or vegdist().
Several functions expect dissimilarities like returned by dist() and
vegdist() and they may get confused if they get a square data matrix
instead. In this case metaMDS is clever enough to reallize that your
square matrix may be dissimilarities and works correctly, but adonis()
was not quite as clever: it was fooled to think that input was a 
data.matrix and calculated dissimilarities on that. So it analysed
dissimilarities of dissimilarities. This may give correct looking
results but they are all wrong. Either give the result of vegdist() 
unchanged or give the original data matrix -- both adonis and 
metaMDS can handle either. Obviously we must change adonis() so
that it realized when the input is dissimilarity data in disguise
of data.

2) I did not see any error messages, but I saw one warning.
Because you gave dissimilrities between OTUs (?) to metaMDS instead 
of the original data, metaMDS() did not know anything about
column variables (species), and therefore plot() could
not plot them even if it tried. It warned you about this.
If you want to have even those in your graph, you should use
your data as input (either step3 or step4 -- I have no idea which
one actually makes sense).


HTH, Jari Oksanen
 
 #NMDS
 step1-read.delim2(day20.txt, row.names=1)
 library(clusterSim)
 step2-data.Normalization(step1,type=n10)
 step3-asin(sqrt(step2))*57.3
 step4-t(step3)
 library(vegan)
 step5-data.matrix(vegdist(step4,method=bray))
 step6-metaMDS(step5, autotransform=FALSE)
 plot(step5)
 Warning message:
 In ordiplot(x, choices = choices, type = type, display = display,  :
   Species scores not available
 
 #PerMANOVA
 step1-read.delim2(day20.txt, row.names=1)
 library(clusterSim)
 step2-data.Normalization(step1,type=n10)
 step3-asin(sqrt(step2))*57.3
 step4-t(step3)
 library(vegan)
 step5-data.matrix(vegdist(step4,method=bray))
 attach(step5)
 step6-adonis(step5~Sample1+Sample2+Sample3,data=step5)
 
   [[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] pie chart

2011-09-09 Thread Jim Lemon

On 09/09/2011 03:35 AM, Mohan L wrote:

Hi All,

I have txt file like :

$ cat data.txt
US  10
UK  12
Ind 4
Germany 14
France  8


rawdata- read.table(file='data.txt',sep='\t' , header=FALSE)



rawdata

V1 V2
1  US 10
2  UK 12
3 Ind  4
4 Germany 14
5  France  8

I want to draw pie chart for the above data.

How to  split rawdata into :
con- c(US,UK,Ind,Germany,France);
total- (10,12,4,14,8)


Hi Mohan,
Anything wrong with:

pie(rawdata$V2,labels=rawdata$V1)

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] envfit vector labels with ordiplot3d

2011-09-09 Thread Jari Oksanen
Briony brionynorton at gmail.com writes:

 
 Hi R experts,
 
 I'm looking for some help with plotting vectors from envfit in vegan, onto a
 3d plot using ordiplot3d. So far I have
 
 data.mds - metaMDS(data, k=3,trace = FALSE)
 vect_data-envfit(data.mds,vegdata[,3:21],choices=1:3,permu=)
 ordiplot3d(data.mds,envfit=vect_data)
 ordixyplot(data.mds,pch=pts,envfit=vect_data)
 
 (my data's not really called data, I thought it might be easier to
 communicate this way)
 
 These display the vectors as arrows, but what I would really like is for the
 arrows to be labelled, like what comes up automatically in ordirgl or with a
 2D ordiplot.
 
 I've gone through the help and tried everything I can work out, but I must
 be missing something important, because nothing's worked so far. I would be
 happy to use ordixyplot and show a series of 2D plots, but I can't get
 labels on those arrows either. 
 
 Any pointers in the right direction would be gratefully received.
 Briony

Briony,

There really is no way to do this automatically, but if someone fixes the
functions, we are happy to incorporate those changes in vegan. 

You may be able to achieve something like that with ordiplot3d, but I am
not sure it looks completely satisfactory. Function ordiplot3d returns 
invisibly the plotting object which contains, among other items. the
coordinates of arrow heads in the  flattened graph. So this could work:

pl - ordiplot3d(data.mds,envfit=vect_data)
ordilabel(pl$arrows)

Function ordiplot3d uses scatterplot3d, and it returns also all 
scatterplot3d items, like functions xyz.converter and points3d that
can be used for tuning labels.

With ordixyplot I can see no other choice than that you edit the
function and preferably contribute your edited function to vegan 
(and will be credited with the function help).

Cheers, Jari Oksanen

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] rgl: axis/viewport/box problems in persp3d()

2011-09-09 Thread Marius Hofert
Dear Duncan,

thanks for your quick response. 
Below is my second trial. I had to use mtext3d to place the label for the 
z-axis 
at the new axis where the ticks are drawn (if there's a simpler solution, please
let me know). Was the usage of rgl.viewpoint meant this way? It is nice to 
adjust the rotation of the plotted object by hand but then I want to make sure
the subsequent plot(s) have precisely the same rotation. 

Okay, great. 

One more thing I am wondering is: I tried to pass through arguments like
marklen or expand to rgl.bbox/bbox3d. Is anything like this possible? I would 
like
to change the length of the axis ticks.

Cheers  many thanks,

Marius

PS: I read somewhere that plotmath-expressions are not available in rgl. Is 
there
an update on this? I know it may be very difficult to implement this, I'm just 
wondering if there is an update/workaround on this (?)

require(rgl)
s - seq(0, 1, length.out=21)
M - function(u) apply(u, 1, min)
u - s
v - s
z - outer(u, v, function(u,v) M(cbind(u,v)))
persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=, 
ylab=, zlab=) 
axes3d(edges=c('x--','y--','z+-')) # label the right axes
title3d(xlab=x, ylab=y, zlab=) # put in axes labels [z is wrong]
mtext3d(z, edge='z+-', line=2) # put in z-axis label by hand
par3d(windowRect=c(0,0,480,480), zoom=1.2) # use zoom to get everything on the 
viewport; then adjust rotation by hand
pl - par3d(c(userMatrix, zoom, FOV)) # record for use in other plots
rgl.postscript(myplot.pdf, fmt=pdf) # print to file
rgl.viewpoint(zoom=pl$zoom, fov=pl$FOV, userMatrix=pl$userMatrix, 
interactive=FALSE) # set the viewpoint for the next plot to make sure it looks 
the same

On 2011-09-09, at 12:41 , Duncan Murdoch wrote:

 On 11-09-09 6:18 AM, Marius Hofert wrote:
 Dear expeRts,
 
 I am a new user of rgl, below is my first trial to plot a simple function in 
 3d.
 I managed to put the axes in the right locations, but:
 (1) The xlab, ylab, and zlab arguments are ignored; how can I put in axes 
 labels?
 
 Those are documented on the axes3d page, but are arguments to title3d, not 
 axes3d.  So add title3d(xlab=x, etc.
 
 (2) Since I removed the axes in persp3d() the viewport is too small; is it 
 possible
 to keep the size of the viewport?
 
 You can manually adjust it to your taste, then write down the value of 
 par3d(zoom).  Later you can reproduce the resizing by calling par3d(zoom= 
 saved value ).
 
 
 (3) The box is not correctly drawn, there are two holes, one in (0,0,1) 
 and one
 in (1,1,0); how can I fix that?
 
 That happens because OpenGL has a limit on the range of depths that can be 
 displayed, and the corners of the box have been adjusted to be too close or 
 far.  This is arguably a bug in rgl, but it's sometimes a feature.
 
 What I'd suggest is that you don't use rgl.viewpoint, you just manually 
 adjust the display as you like, without making it quite as extreme, then 
 record the values of par3d(c(userMatrix, zoom, FOV)); those control the 
 viewpoint.
 
 Duncan Murdoch
 
 
 Cheers,
 
 Marius
 
 
 require(rgl)
 s- seq(0, 1, length.out=21)
 M- function(u) apply(u, 1, min)
 u- s
 v- s
 z- outer(u, v, function(u,v) M(cbind(u,v)))
 persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=,
 ylab=, zlab=)
 axes3d(edges=c('x--','y--','z+-'), xlab=x, ylab=y, zlab=z)
 par3d(windowRect=c(0,0,480,480))
 
 R1- rotationMatrix(-55*pi/180, 1,0,0)
 R3- rotationMatrix(50*pi/180, 0,0,1)
 R- R1 %*% R3
 rgl.viewpoint(interactive=TRUE, userMatrix=R) # rotate
 rgl.postscript(myplot.pdf, fmt=pdf)
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/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] reshape data from long to wide format

2011-09-09 Thread maxbre
This is my reproducible example:

example-structure(list(SENSOR = structure(1:6, .Label = c(A, B, C, 
D, E, F), class = factor), VALUE = c(270, 292.5, 0, 45, 
247.5, 315), DATE = structure(1:6, .Label = c( 01/01/2010 1, 
 01/01/2010 2,  01/01/2010 3,  01/01/2010 4,  01/01/2010 5, 
 01/01/2010 6), class = factor)), .Names = c(SENSOR, VALUE, 
DATE), class = data.frame, row.names = c(1, 2, 3, 4, 
5, 6))

I need to resahpe example  in a wide format so that “SENSOR” appear as
columns and “DATE” as rows with corresponding “VALUE” as value;

I thought it was very simple so that I’ve been trying this:

dcast(example,DATE~SENSOR)

But I've got this message:

Using DATE as value column.  Use the value argument to cast to override this
choice
Errore in `[.data.frame`(data, , variables, drop = FALSE) : 
  undefined columns selected

and even if by using the value argument sorted out any good result…

sorry for the very trivial question, I've been looking at documentation and
forums anywhere but I was not successful at all and now I’m somehow in the
right middle of nowhere…

any help or hint for this?

thank you

max

--
View this message in context: 
http://r.789695.n4.nabble.com/reshape-data-from-long-to-wide-format-tp3801381p3801381.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] Subset function

2011-09-09 Thread stat.kk
Hi,

can anyone help me how to use 'subset' function on my data frame?
I have created data frame 'data' with a few variables and with row names.
Now I would like to subset rows with concrete row names. 
Using data[] I know how to do it. But I dont know how to formulate the
subset condition:
subset(data, subset = ?, select = c(var1, var2))

Thank you very much,
stat.kk

--
View this message in context: 
http://r.789695.n4.nabble.com/Subset-function-tp3801397p3801397.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] Survival Analysis for soccer scoring process

2011-09-09 Thread ryusuke
6.4.1 Estimation of fixed effects
Heterogeneous team ability is a possible explanation for the result in
Section 6.3. That
result simply indicates that the more goals a team scores, the higher the
probability
that it will score more. However, teams that can score more goals also
indicate teams
with greater ability, or just greater scoring ability, than their opponents.
As mentioned
in Section 6.1.1, it is very natural that the higher ability team will score
more, even
though that stronger team has decided to defend more. Therefore, we need to
control
for this in the estimation.

This problem is similar, but not identical, to the problem identified by
Heckman and
Singer (1984), that heterogeneous subjects will tend to yield decreasing
hazards
because a longer time of no failure indicates that the particular sample is
intrinsically
not likely to fail. Here, heterogeneous teams tend to yield a higher hazard
of scoring
for the leading team because such teams are intrinsically more likely to
score.

*For the estimation controlled for fixed effects, a series of dummy
variables are added.
Two dummy variables are added for each team, one for home ground and another
for
away ground (with the exception of Manchester City to prevent colinearity
and act as
baseline hazard). * As there are in total 25 teams6 involving in the English
Premiership,
there are 48 dummy variables added. Similarly, there are 22 and 25 teams in
the
German Bundesliga and Spanish Primera Liga and therefore I add 42 and 48
dummy
variables respectively.

Fitting the model by these dummy variables only can yield the ability
difference
between teams. Manchester City’s hazard ratio is equal to 1 in all columns
because it is the
benchmark. Arsenal is obviously a strong team as it has very high hazard
ratios for
scoring (1.40 of Manchester City’s hazard at home) and very low hazard
ratios for
being scored against (0.81 of Manchester City’s hazard at home)


I am currently modelling soccer scoring process, may I know how do I set a
variable as baseline hazard benchmark (as Bold characters inside above
articles)? Hope some body willing to share your precious ideas. Thanks.

--
View this message in context: 
http://r.789695.n4.nabble.com/Survival-Analysis-for-soccer-scoring-process-tp3801400p3801400.html
Sent from the R help mailing list archive at Nabble.com.
[[alternative HTML version deleted]]

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


[R] How to translate the 2D-density matrix (the output of bkde2D function) into matrix of datapoints' amounts?

2011-09-09 Thread Юрий Посохов
It is known that function bkde2D (package KernSmooth) returns a
matrix of density estimates over the mesh induced by x1 and x2. In
Details it is written that ... heights of the kernel, scaled by the
bandwidths, at each datapoint are summed. This sum, after a
normalization, is the corresponding fhat value in the output.
There are several questions:
1) How to calculate amount (sum) of datapoints from this fhat value?
Is it a non-normalized matrix of fhat values?
2) Values in fhat matrix are greater than 1 sometimes, why? Values of
density are normalized and hence must be less than 1, isn't it?

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


Re: [R] problems with function read.table

2011-09-09 Thread Petr PIKAL
Hi

 Hello,
 
 in the meanwhile i found the problem for the #. The problem was, that # 
and
 “none” are comments by default, so I turned comments of with 
comment.char=””
 as stated in the help.

Hm. You are not telling the whole story.
I made a sample excel sheet with blank values and #DIV/0! result (simple 
two columns with missing values and division by these missing values in 
first column)

 x-read.delim(clipboard)
 x
r a  b
11.25 5  4
2 1.16667 7  6
3 0.142857143 1  7
4 #DIV/0! 5 NA
5 4.5 9  2
6   2 6  3
7 #DIV/0! 4 NA

simple read.delim gave me correctly NA in place of blank values and 
#DIV/0! in first column, which therefore became factor.

 y-read.delim(clipboard, na.strings=c(NA, #DIV/0!))
 y
  r a  b
1 1.250 5  4
2 1.167 7  6
3 0.1428571 1  7
4NA 5 NA
5 4.500 9  2
6 2.000 6  3
7NA 4 NA

When I defined vector of na.strings I correctly got NA values for both 
blank and #DIV/0! in Excel sheet. So I wonder where you get problems with 
blank values 

 
 For the problem with blanks I still didn’t find the solution, even it 
seems
 to be a frequent problem. I got the solution with fill=T, but this fills 
my
 row at the end with a NA. I want to have the NA exactly in the field 
where
 there is a blank in the .txt file.

It seems to me that you do not use correct delimiter. For Excel it is tab 
- \t. For csv it can be either ; or , depending on your locale. You 
need to show us some small part of a file, preferably together with your 
read.* command and the result you got.

Regards
Petr

 
 Date  rtn   vwretdewretd
 sprtrn
 19700102 0.000686  0.00547 0.033450
 0.010211
 19700105 0.009596 0.018947
 0.004946
 19700106 #DIV0!-0.007233
 -0.006848
 19700107 0.000678  -0.001272
 0.003559  -0.002047
 19700108 0.002034  0.0005640.11
 0.000540
 19700109  -0.002797   0
 -0.003021
 19700113 0.017335  0.000737 -0.001090
 0
 
 Can you provide a solution?
 
 Thanks,
 Samir
 
 p.s. I know that my questions seem obvious to you, I’m sorry for that, 
but I
 just started to work with R ;)
 
 -Ursprüngliche Nachricht-
 Von: Petr PIKAL [mailto:petr.pi...@precheza.cz] 
 Gesendet: Freitag, 9. September 2011 09:23
 An: Carlos Ortega
 Cc: Samir Benzerfa; r-help@r-project.org
 Betreff: Re: [R] problems with function read.table
 
 Hi
 
 
  
  Hi,
  
  If you read carefully the help pages for read.table you get this:
  
  
  na.stringsa character vector of strings which are to be interpreted as
  NA../../utils/help/NA values.
  Blank fields are also considered to be missing values in logical, 
 integer,
  numeric and complex fields.
  
  So, both NAs and blank fields are considered as NAs directly by 
 read.table.
  
  Once you have imported your data, you can modify with any of the 
string
  manipulation functions (sub() or gsub()) to change your #DIV/0! to 
the
  string NAs. Another option is to manipulate your Excel file and 
 consider
  the division by cero with a IF and get back a NA if that happens.
 
 The only problem is that in such case all columns which has #DIV/0! 
are 
 converted to factors and you need to consider changing it back to 
numeric.
 
 read.* functions accept as na.string definition not only one value but 
 also vector of values and you can get rid of all non numeric and other 
 weird Excel values by defining it as a na.strings in read.table call.
 
  x - read.delim(clipboard)
 
  str(x)
 'data.frame':   6 obs. of  3 variables:
  $ a: int  1 5 9 8 6 3
  $ b: int  3 5 7 0 NA 6
  $ r: Factor w/ 5 levels #DIV/0!,0.3,..: 2 4 5 1 1 3
 
  y-read.delim(clipboard, na.strings=c(NA, #DIV/0!))
  str(y)
 'data.frame':   6 obs. of  3 variables:
  $ a: int  1 5 9 8 6 3
  $ b: int  3 5 7 0 NA 6
  $ r: num  0.333 1 1.286 NA NA ...
 
 
 Regards
 Petr
 
 
  
  And finally, instead of using na.omits use option na.rm=T to get done 
 your
  calculations:
  
   mean(c(12,23,24,45,67,NA), na.rm=T)[1] 34.2
  
  
  
  Regards,
  Carlos Ortega
  www.qualityexcellence.es
  
  On Thu, Sep 8, 2011 at 4:23 PM, Samir Benzerfa benze...@gmx.ch 
wrote:
  
   Hello everyone
  
  
  
   I have a couple of questions about the usage of the R function
   read.table(.). My point of departure is that I want to import a 
 matrix
   (consisting of time and daily stock returns of many stocks) in R. 
Most 
 of
   the data is numeric, however some values are missing (blanks) and in 

 other
   cases I have the character #DIV/0! (from excel). My goal is to do 
 some
   regression analysis with this matrix. My questions now are the 
 following
   ones:
  
  
  
   1.   How can I in general tell R to automatically replace some 
 specific
   numbers or characters in tables by others? (for example 

Re: [R] rgl: axis/viewport/box problems in persp3d()

2011-09-09 Thread Duncan Murdoch

On 09/09/2011 8:02 AM, Marius Hofert wrote:

Dear Duncan,

thanks for your quick response.
Below is my second trial. I had to use mtext3d to place the label for the z-axis
at the new axis where the ticks are drawn (if there's a simpler solution, please
let me know). Was the usage of rgl.viewpoint meant this way?


I would have used par3d(zoom=pl$zoom, ... ) instead of rgl.viewpoint, 
but it probably gives the same result.  You'd have to check the sources.



  It is nice to
adjust the rotation of the plotted object by hand but then I want to make sure
the subsequent plot(s) have precisely the same rotation.

Okay, great.

One more thing I am wondering is: I tried to pass through arguments like
marklen or expand to rgl.bbox/bbox3d. Is anything like this possible? I would 
like
to change the length of the axis ticks.


You can use marklen in bbox3d, but it won't affect the axes that were 
drawn with axes3d.  This is something that has been on my todo-list for 
a long time, uniting the two different ways of drawing axes (static ones 
from axes3d, dynamic ones from bbox3d).  To modify the static tick 
length, you'll need to modify the source to axis3d.

Cheers  many thanks,

Marius

PS: I read somewhere that plotmath-expressions are not available in rgl. Is 
there
an update on this? I know it may be very difficult to implement this, I'm just
wondering if there is an update/workaround on this (?)


No, unfortunately not.   That's also on my todo-list, but fairly far 
down.  (I wouldn't implement plotmath in rgl; I would implement a 
standard graphics driver in rgl, so all graphics functions could draw to 
a plane in 3d space.  But that's really a lot of work.)


Duncan Murdoch


require(rgl)
s- seq(0, 1, length.out=21)
M- function(u) apply(u, 1, min)
u- s
v- s
z- outer(u, v, function(u,v) M(cbind(u,v)))
persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=,
 ylab=, zlab=)
axes3d(edges=c('x--','y--','z+-')) # label the right axes
title3d(xlab=x, ylab=y, zlab=) # put in axes labels [z is wrong]
mtext3d(z, edge='z+-', line=2) # put in z-axis label by hand
par3d(windowRect=c(0,0,480,480), zoom=1.2) # use zoom to get everything on the 
viewport; then adjust rotation by hand
pl- par3d(c(userMatrix, zoom, FOV)) # record for use in other plots
rgl.postscript(myplot.pdf, fmt=pdf) # print to file
rgl.viewpoint(zoom=pl$zoom, fov=pl$FOV, userMatrix=pl$userMatrix, 
interactive=FALSE) # set the viewpoint for the next plot to make sure it looks 
the same

On 2011-09-09, at 12:41 , Duncan Murdoch wrote:

  On 11-09-09 6:18 AM, Marius Hofert wrote:
  Dear expeRts,

  I am a new user of rgl, below is my first trial to plot a simple function 
in 3d.
  I managed to put the axes in the right locations, but:
  (1) The xlab, ylab, and zlab arguments are ignored; how can I put in axes 
labels?

  Those are documented on the axes3d page, but are arguments to title3d, not axes3d.  So 
add title3d(xlab=x, etc.

  (2) Since I removed the axes in persp3d() the viewport is too small; is it 
possible
  to keep the size of the viewport?

  You can manually adjust it to your taste, then write down the value of 
par3d(zoom).  Later you can reproduce the resizing by calling par3d(zoom=saved 
value  ).


  (3) The box is not correctly drawn, there are two holes, one in (0,0,1) 
and one
  in (1,1,0); how can I fix that?

  That happens because OpenGL has a limit on the range of depths that can be 
displayed, and the corners of the box have been adjusted to be too close or far.  
This is arguably a bug in rgl, but it's sometimes a feature.

  What I'd suggest is that you don't use rgl.viewpoint, you just manually adjust the display as you like, 
without making it quite as extreme, then record the values of par3d(c(userMatrix, zoom, 
FOV)); those control the viewpoint.

  Duncan Murdoch


  Cheers,

  Marius


  require(rgl)
  s- seq(0, 1, length.out=21)
  M- function(u) apply(u, 1, min)
  u- s
  v- s
  z- outer(u, v, function(u,v) M(cbind(u,v)))
  persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=,
  ylab=, zlab=)
  axes3d(edges=c('x--','y--','z+-'), xlab=x, ylab=y, zlab=z)
  par3d(windowRect=c(0,0,480,480))

  R1- rotationMatrix(-55*pi/180, 1,0,0)
  R3- rotationMatrix(50*pi/180, 0,0,1)
  R- R1 %*% R3
  rgl.viewpoint(interactive=TRUE, userMatrix=R) # rotate
  rgl.postscript(myplot.pdf, fmt=pdf)
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/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] How to translate the 2D-density matrix (the output of bkde2D function) into matrix of datapoints' amounts?

2011-09-09 Thread omega1x
It is known that function bkde2D (package KernSmooth) returns a matrix of
density estimates over the mesh induced by x1 and x2. In Details it is
written that ... heights of the kernel, scaled by the bandwidths, at each
datapoint are summed. This sum, after a normalization, is the corresponding
fhat value in the output. 
There are several questions:
1) How to calculate amount (sum) of datapoints from this fhat value? Is it a
non-normalized matrix of fhat values?
2) Values in fhat matrix are greater than 1 sometimes, why? Values of
density are normalized and hence must be less than 1, isn't it?

--
View this message in context: 
http://r.789695.n4.nabble.com/How-to-translate-the-2D-density-matrix-the-output-of-bkde2D-function-into-matrix-of-datapoints-amoun-tp3801482p3801482.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] function censReg in panel data setting

2011-09-09 Thread Arne Henningsen
On 8 September 2011 09:56, Igors igors.lahanc...@gmail.com wrote:
 Does censReg expect from panel data to be balanced?

No. censReg() can estimate models with unbalanced panel data. The
estimation in the code that you sent me indeed does not work but if I
remove the user-specified starting values (argument start), it works
fine.

R summary(UpC)

Call:
censReg(formula = Power ~ Windspeed, left = -Inf, right = 2000,
data = PData_In, nGHQ = 4, method = BHHH, iterlim = 200)

Observations:
 Total  Left-censored Uncensored Right-censored
   874  0847 27

Coefficients:
  Estimate Std. error t value Pr( t)
(Intercept) -462.26676   19.89517  -23.23  2e-16 ***
Windspeed188.387961.72492  109.22  2e-16 ***
logSigmaMu 4.540530.03352  135.45  2e-16 ***
logSigmaNu 5.087330.02706  188.00  2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

BHHH maximisation, 131 iterations
Return code 2: successive function values within tolerance limit
Log-likelihood: -5538.124 on 4 Df


/Arne


-- 
Arne Henningsen
http://www.arne-henningsen.name

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

2011-09-09 Thread Rainer Schuermann
Does that help:

 x
  xin xout
1   1   14
2   85
3  16  884
4   1   14
5   85  
  
6  16  884

 subset( x, x$xin  7, select = xout ) 
  xout  
  
25  
  
3  884  
  
55  
  
6  884 

Rgds,
Rainer


On Friday 09 September 2011 04:38:44 stat.kk wrote:
 Hi,
 
 can anyone help me how to use 'subset' function on my data frame?
 I have created data frame 'data' with a few variables and with row names.
 Now I would like to subset rows with concrete row names.
 Using data[] I know how to do it. But I dont know how to formulate the
 subset condition:
 subset(data, subset = ?, select = c(var1, var2))
 
 Thank you very much,
 stat.kk
 
 --
 View this message in context:
 http://r.789695.n4.nabble.com/Subset-function-tp3801397p3801397.html Sent
 from the R help mailing list archive at Nabble.com.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] Subset function

2011-09-09 Thread Jean V Adams
Rainer provided an example of subsetting by the value of a variable in the 
data frame.  Below is an example of subsetting by the value of the row 
name of the data frame.

df - data.frame(var1=1:10, var2=letters[1:10], var3=sample(10), 
 row.names=month.abb[1:10])
subset(df, subset = row.names(df) %in% c(Feb, Jul, Sep), 
 select=c(var1, var2))

var1 var2
Feb2b
Jul7g
Sep9i

Jean


Rainer Schuermann wrote on 09/09/2011 07:44:11 AM:
 
 Does that help:
 
  x
   xin xout
 1   1   14
 2   85
 3  16  884
 4   1   14
 5   85  
 6  16  884
 
  subset( x, x$xin  7, select = xout )  
   xout  
 25  
 3  884  
 55  
 6  884 
 
 Rgds,
 Rainer
 
 
 On Friday 09 September 2011 04:38:44 stat.kk wrote:
  Hi,
  
  can anyone help me how to use 'subset' function on my data frame?
  I have created data frame 'data' with a few variables and with row 
names.
  Now I would like to subset rows with concrete row names.
  Using data[] I know how to do it. But I dont know how to formulate the
  subset condition:
  subset(data, subset = ?, select = c(var1, var2))
  
  Thank you very much,
  stat.kk
  

[[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] reshape data from long to wide format

2011-09-09 Thread Jean V Adams
maxbre wrote on 09/09/2011 06:28:15 AM:
 
 This is my reproducible example:
 
 example-structure(list(SENSOR = structure(1:6, .Label = c(A, B, 
C, 
 D, E, F), class = factor), VALUE = c(270, 292.5, 0, 45, 
 247.5, 315), DATE = structure(1:6, .Label = c( 01/01/2010 1, 
  01/01/2010 2,  01/01/2010 3,  01/01/2010 4,  01/01/2010 5, 
  01/01/2010 6), class = factor)), .Names = c(SENSOR, VALUE, 
 DATE), class = data.frame, row.names = c(1, 2, 3, 4, 
 5, 6))
 
 I need to resahpe example  in a wide format so that ?SENSOR? appear as
 columns and ?DATE? as rows with corresponding ?VALUE? as value;


It's not clear to me what the wide format of your example would look like. 
 You start with this:

  SENSOR VALUE  DATE
1  A 270.0  01/01/2010 1
2  B 292.5  01/01/2010 2
3  C   0.0  01/01/2010 3
4  D  45.0  01/01/2010 4
5  E 247.5  01/01/2010 5
6  F 315.0  01/01/2010 6

And you want to change it to ... something like this?

 DATE   A B  C  D E   F
 01/01/2010 1 270NA NA NANA  NA
 01/01/2010 2  NA 292.5 NA NANA  NA
 01/01/2010 3  NANA  0 NANA  NA
 01/01/2010 4  NANA NA 45NA  NA
 01/01/2010 5  NANA NA NA 247.5  NA
 01/01/2010 6  NANA NA NANA 315

Jean


 
 I thought it was very simple so that I?ve been trying this:
 
 dcast(example,DATE~SENSOR)
 
 But I've got this message:
 
 Using DATE as value column.  Use the value argument to cast to override 
this
 choice
 Errore in `[.data.frame`(data, , variables, drop = FALSE) : 
   undefined columns selected
 
 and even if by using the value argument sorted out any good result?
 
 sorry for the very trivial question, I've been looking at documentation 
and
 forums anywhere but I was not successful at all and now I?m somehow in 
the
 right middle of nowhere?
 
 any help or hint for this?
 
 thank you
 
 max
 
[[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] rgl: axis/viewport/box problems in persp3d()

2011-09-09 Thread Marius Hofert
Dear Duncan,

thanks a lot. 
Is it possible to rotate the label drawn by mtext3d, say, by 90 degrees? [a 
rot=90 did not help]

Cheers,

Marius

On 2011-09-09, at 14:32 , Duncan Murdoch wrote:

 On 09/09/2011 8:02 AM, Marius Hofert wrote:
 Dear Duncan,
 
 thanks for your quick response.
 Below is my second trial. I had to use mtext3d to place the label for the 
 z-axis
 at the new axis where the ticks are drawn (if there's a simpler solution, 
 please
 let me know). Was the usage of rgl.viewpoint meant this way?
 
 I would have used par3d(zoom=pl$zoom, ... ) instead of rgl.viewpoint, but it 
 probably gives the same result.  You'd have to check the sources.
 
  It is nice to
 adjust the rotation of the plotted object by hand but then I want to make 
 sure
 the subsequent plot(s) have precisely the same rotation.
 
 Okay, great.
 
 One more thing I am wondering is: I tried to pass through arguments like
 marklen or expand to rgl.bbox/bbox3d. Is anything like this possible? I 
 would like
 to change the length of the axis ticks.
 
 You can use marklen in bbox3d, but it won't affect the axes that were drawn 
 with axes3d.  This is something that has been on my todo-list for a long 
 time, uniting the two different ways of drawing axes (static ones from 
 axes3d, dynamic ones from bbox3d).  To modify the static tick length, you'll 
 need to modify the source to axis3d.
 Cheers  many thanks,
 
 Marius
 
 PS: I read somewhere that plotmath-expressions are not available in rgl. Is 
 there
 an update on this? I know it may be very difficult to implement this, I'm 
 just
 wondering if there is an update/workaround on this (?)
 
 No, unfortunately not.   That's also on my todo-list, but fairly far down.  
 (I wouldn't implement plotmath in rgl; I would implement a standard graphics 
 driver in rgl, so all graphics functions could draw to a plane in 3d space.  
 But that's really a lot of work.)
 
 Duncan Murdoch
 
 require(rgl)
 s- seq(0, 1, length.out=21)
 M- function(u) apply(u, 1, min)
 u- s
 v- s
 z- outer(u, v, function(u,v) M(cbind(u,v)))
 persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=,
 ylab=, zlab=)
 axes3d(edges=c('x--','y--','z+-')) # label the right axes
 title3d(xlab=x, ylab=y, zlab=) # put in axes labels [z is wrong]
 mtext3d(z, edge='z+-', line=2) # put in z-axis label by hand
 par3d(windowRect=c(0,0,480,480), zoom=1.2) # use zoom to get everything on 
 the viewport; then adjust rotation by hand
 pl- par3d(c(userMatrix, zoom, FOV)) # record for use in other plots
 rgl.postscript(myplot.pdf, fmt=pdf) # print to file
 rgl.viewpoint(zoom=pl$zoom, fov=pl$FOV, userMatrix=pl$userMatrix, 
 interactive=FALSE) # set the viewpoint for the next plot to make sure it 
 looks the same
 
 On 2011-09-09, at 12:41 , Duncan Murdoch wrote:
 
   On 11-09-09 6:18 AM, Marius Hofert wrote:
   Dear expeRts,
 
   I am a new user of rgl, below is my first trial to plot a simple 
  function in 3d.
   I managed to put the axes in the right locations, but:
   (1) The xlab, ylab, and zlab arguments are ignored; how can I put in 
  axes labels?
 
   Those are documented on the axes3d page, but are arguments to title3d, 
  not axes3d.  So add title3d(xlab=x, etc.
 
   (2) Since I removed the axes in persp3d() the viewport is too small; is 
  it possible
   to keep the size of the viewport?
 
   You can manually adjust it to your taste, then write down the value of 
  par3d(zoom).  Later you can reproduce the resizing by calling 
  par3d(zoom=saved value  ).
 
 
   (3) The box is not correctly drawn, there are two holes, one in 
  (0,0,1) and one
   in (1,1,0); how can I fix that?
 
   That happens because OpenGL has a limit on the range of depths that can 
  be displayed, and the corners of the box have been adjusted to be too 
  close or far.  This is arguably a bug in rgl, but it's sometimes a feature.
 
   What I'd suggest is that you don't use rgl.viewpoint, you just manually 
  adjust the display as you like, without making it quite as extreme, then 
  record the values of par3d(c(userMatrix, zoom, FOV)); those control 
  the viewpoint.
 
   Duncan Murdoch
 
 
   Cheers,
 
   Marius
 
 
   require(rgl)
   s- seq(0, 1, length.out=21)
   M- function(u) apply(u, 1, min)
   u- s
   v- s
   z- outer(u, v, function(u,v) M(cbind(u,v)))
   persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, 
  xlab=,
   ylab=, zlab=)
   axes3d(edges=c('x--','y--','z+-'), xlab=x, ylab=y, zlab=z)
   par3d(windowRect=c(0,0,480,480))
 
   R1- rotationMatrix(-55*pi/180, 1,0,0)
   R3- rotationMatrix(50*pi/180, 0,0,1)
   R- R1 %*% R3
   rgl.viewpoint(interactive=TRUE, userMatrix=R) # rotate
   rgl.postscript(myplot.pdf, fmt=pdf)
   __
   R-help@r-project.org mailing list
   https://stat.ethz.ch/mailman/listinfo/r-help
   PLEASE do read the posting guide 
  http://www.R-project.org/posting-guide.html
   and provide commented, minimal, self-contained, 

Re: [R] Subset function

2011-09-09 Thread stat.kk
No :(
I dont have a condition on variable but on row name. My data frame looks
like (from your example):
 x
 xin xout
Peter1   14
Tom 8 5
Jane   16  884
Paul 114
Cathy   8 5 
   
Fred   16  884 

I would like to subset xout for Jane and Paul for example.

--
View this message in context: 
http://r.789695.n4.nabble.com/Subset-function-tp3801397p3801605.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] conditional Latin hypercube sampling

2011-09-09 Thread Eric Tiger
Hello,

I got one question on the Latin hypercube sampling.

suppose
 there are three variables a, b, c, all of them follow the normal 
distribution. the mean value and standard deviation for each are  a(32, 
2), b(35,5), c(37,3). I would like to use Latin hypercube sampling to 
random generate 1000 samples. but it needs to satisfy the condition that
 abc.

How can I implement this conditional sampling?

Thanks a lot.
Frank
[[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] Subset function

2011-09-09 Thread stat.kk
Thank you very much :)

stat.kk

--
View this message in context: 
http://r.789695.n4.nabble.com/Subset-function-tp3801397p3801626.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] analyse fMRI statistical maps in a mixed procedure (within-subject and between-subject error term)

2011-09-09 Thread Diederick Stoffers
Hi all,

I have a set of 3D statistical parametric maps derived from an fMRI experiment 
in which we have two levels of dependency.

- All subjects were scanned twice.
- All subject are dizygotic (non-identical) twins.

For every single scan (and thus statical map) we have a behavioural measure, 
which we would like to regress against the set of maps, while correcting for 
the dependency between scans and ideally also between twins. 

In neuroimaging software like SPM and FSL it is not possible to have a 
within-subject error term as well as a between-subject error term. The advice I 
got from multiple people on the FSL and SPM mailinglists was to just average 
the spatial maps over sessions and possibly also twins. However, that would 
mean also averaging over the behavioural measure which can actually be quite 
different over sessions as well as twins, in this way greatly reducing power. 

I was told that it might be possible to analyse these statistical parametric 
maps within R with a mixed procedure. Can anyone confirm that this is the case 
and perhaps point me in the right direction, e.g. with a publication or other 
helpful document, webpage or package?

Cheers,

Diederick
[[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] How to translate the 2D-density matrix (the output of bkde2D function) into matrix of datapoints' amounts?

2011-09-09 Thread David Winsemius


On Sep 9, 2011, at 8:12 AM, omega1x wrote:

It is known that function bkde2D (package KernSmooth) returns a  
matrix of

density estimates over the mesh induced by x1 and x2. In Details it is
written that ... heights of the kernel, scaled by the bandwidths,  
at each
datapoint are summed. This sum, after a normalization, is the  
corresponding

fhat value in the output.
There are several questions:
1) How to calculate amount (sum) of datapoints from this fhat value?  
Is it a

non-normalized matrix of fhat values?


Wouldn't you just use 'cut' or 'findInterval' using the mesh points  
and then 'table' or 'tabulate'?



2) Values in fhat matrix are greater than 1 sometimes, why? Values of
density are normalized and hence must be less than 1, isn't it?


No. One would expect that the normalization is done to get the  
integral equal to 1, not the individual points to have fhat in the  
range of [0,1]. (Failure to understand the difference is a common  
cause of questions on Rhelp, so if you need further explanation and  
many examples then do a search.)




--
View this message in context: 
http://r.789695.n4.nabble.com/How-to-translate-the-2D-density-matrix-the-output-of-bkde2D-function-into-matrix-of-datapoints-amoun-tp3801482p3801482.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.


[R] Read a list of files into named R data.frames

2011-09-09 Thread Michael Friendly
I have a collection of .csv files in a directory, and want to read them 
into R data.frames whose names

are the same as the file names, without the .csv extension

e.g., from
 (files - list.files(pattern=*.csv))
 [1] Allstar.csv AllstarFull.csv
 [3] Appearances.csv AwardsManagers.csv
 [5] AwardsPlayers.csv   AwardsShareManagers.csv
 [7] AwardsSharePlayers.csv  Batting.csv
 [9] BattingPost.csv Fielding.csv
[11] FieldingOF.csv  FieldingPost.csv
[13] HallOfFame.csv  HOFold.csv
[15] Managers.csvManagersHalf.csv
[17] Master.csv  Pitching.csv
[19] PitchingPost.csvSalaries.csv
[21] Schools.csv SchoolsPlayers.csv
[23] SeriesPost.csv  Teams.csv
[25] TeamsFranchises.csv TeamsHalf.csv

 Allstar - read.csv(Allstar.csv, header=TRUE)
  ...
 TeamsHalf - read.csv(TeamsHalf.csv, header=TRUE)

Below is what I tried, which reads all the files, but doesn't create the 
R objects in the global environment.

What is missing here?

for (i in 1:length(files)) {
inp - read.csv(file=files[i], header=TRUE)
name - sub(.csv, , files[i])
cat(Read , files[i], \trows: , nrow(inp),  cols: , ncol(inp), 
\n)

eval(paste(name, - inp))
}

Read  Allstar.csv   rows:  4475  cols:  3
Read  AllstarFull.csv   rows:  4676  cols:  8
Read  Appearances.csv   rows:  94157  cols:  20
Read  AwardsManagers.csvrows:  57  cols:  6
Read  AwardsPlayers.csv rows:  2679  cols:  6
Read  AwardsShareManagers.csv   rows:  344  cols:  7
Read  AwardsSharePlayers.csvrows:  6354  cols:  7
Read  Batting.csv   rows:  93955  cols:  24
Read  BattingPost.csv   rows:  9840  cols:  22
Read  Fielding.csv  rows:  160710  cols:  18
Read  FieldingOF.csvrows:  12028  cols:  6
Read  FieldingPost.csv  rows:  10458  cols:  17
Read  HallOfFame.csvrows:  3913  cols:  8
Read  HOFold.csvrows:  289  cols:  7
Read  Managers.csv  rows:  3238  cols:  10
Read  ManagersHalf.csv  rows:  93  cols:  10
Read  Master.csvrows:  17674  cols:  33
Read  Pitching.csv  rows:  40432  cols:  30
Read  PitchingPost.csv  rows:  4284  cols:  30
Read  Salaries.csv  rows:  21464  cols:  5
Read  Schools.csv   rows:  749  cols:  5
Read  SchoolsPlayers.csvrows:  6147  cols:  4
Read  SeriesPost.csvrows:  256  cols:  9
Read  Teams.csv rows:  2655  cols:  48
Read  TeamsFranchises.csv   rows:  120  cols:  4
Read  TeamsHalf.csv rows:  52  cols:  10
Read  Xref_Stats.csvrows:  2753  cols:  3
 ls()
[1] files i inp   name


--
Michael Friendly Email: friendly AT yorku DOT ca
Professor, Psychology Dept.
York University  Voice: 416 736-5115 x66249 Fax: 416 736-5814
4700 Keele StreetWeb:   http://www.datavis.ca
Toronto, ONT  M3J 1P3 CANADA

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

2011-09-09 Thread William Dunlap
I thought that the main advantage of subset() over [()
is that you only mention the name of the data.frame once,
in the first argument, not in the second:
   x - data.frame(xin=c(1, 8, 16, 1, 8, 16), xout=c(14, 5, 884, 14, 5, 884))
   subset(x, xin  7, select = xout) # not x$xin  7
xout
  25
  3  884
  55
  6  884

A secondary advantage of subset is that treats NA's in
the subset= argument the same as FALSE's.

I think subset is handy for one-off usage, but in general
purpose functions the [ function is better: it uses standard
argument evaluation and is faster.

Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com 

 -Original Message-
 From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
 Behalf Of Rainer
 Schuermann
 Sent: Friday, September 09, 2011 5:44 AM
 To: r-help@r-project.org; stat.kk
 Subject: Re: [R] Subset function
 
 Does that help:
 
  x
   xin xout
 1   1   14
 2   85
 3  16  884
 4   1   14
 5   85
 6  16  884
 
  subset( x, x$xin  7, select = xout )
   xout
 25
 3  884
 55
 6  884
 
 Rgds,
 Rainer
 
 
 On Friday 09 September 2011 04:38:44 stat.kk wrote:
  Hi,
 
  can anyone help me how to use 'subset' function on my data frame?
  I have created data frame 'data' with a few variables and with row names.
  Now I would like to subset rows with concrete row names.
  Using data[] I know how to do it. But I dont know how to formulate the
  subset condition:
  subset(data, subset = ?, select = c(var1, var2))
 
  Thank you very much,
  stat.kk
 
  --
  View this message in context:
  http://r.789695.n4.nabble.com/Subset-function-tp3801397p3801397.html Sent
  from the R help mailing list archive at Nabble.com.
 
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] Read a list of files into named R data.frames

2011-09-09 Thread David Winsemius


On Sep 9, 2011, at 10:39 AM, Michael Friendly wrote:

I have a collection of .csv files in a directory, and want to read  
them into R data.frames whose names

are the same as the file names, without the .csv extension

e.g., from
 (files - list.files(pattern=*.csv))
[1] Allstar.csv AllstarFull.csv
[3] Appearances.csv AwardsManagers.csv
[5] AwardsPlayers.csv   AwardsShareManagers.csv
[7] AwardsSharePlayers.csv  Batting.csv
[9] BattingPost.csv Fielding.csv
[11] FieldingOF.csv  FieldingPost.csv
[13] HallOfFame.csv  HOFold.csv
[15] Managers.csvManagersHalf.csv
[17] Master.csv  Pitching.csv
[19] PitchingPost.csvSalaries.csv
[21] Schools.csv SchoolsPlayers.csv
[23] SeriesPost.csv  Teams.csv
[25] TeamsFranchises.csv TeamsHalf.csv

 Allstar - read.csv(Allstar.csv, header=TRUE)
 ...
 TeamsHalf - read.csv(TeamsHalf.csv, header=TRUE)

Below is what I tried, which reads all the files, but doesn't create  
the R objects in the global environment.

What is missing here?

for (i in 1:length(files)) {
   inp - read.csv(file=files[i], header=TRUE)
   name - sub(.csv, , files[i])
   cat(Read , files[i], \trows: , nrow(inp),  cols: ,  
ncol(inp), \n)


Generally the assign function is used to create objects with a  
particular name. If you wanted to use eval then the text needs to be  
passed through parse() before being given to eval, but that is not the  
preferred method.


Perhaps:

 assign( files[i], inp)

Inside a for loop I think that gets done in the calling environment  
but if you were in a function  you would need to use the environment  
argument to get it to stick.




   eval(paste(name, - inp))
}

Read  Allstar.csv   rows:  4475  cols:  3
Read  AllstarFull.csv   rows:  4676  cols:  8
Read  Appearances.csv   rows:  94157  cols:  20
Read  AwardsManagers.csvrows:  57  cols:  6
Read  AwardsPlayers.csv rows:  2679  cols:  6
Read  AwardsShareManagers.csv   rows:  344  cols:  7
Read  AwardsSharePlayers.csvrows:  6354  cols:  7
Read  Batting.csv   rows:  93955  cols:  24
Read  BattingPost.csv   rows:  9840  cols:  22
Read  Fielding.csv  rows:  160710  cols:  18
Read  FieldingOF.csvrows:  12028  cols:  6
Read  FieldingPost.csv  rows:  10458  cols:  17
Read  HallOfFame.csvrows:  3913  cols:  8
Read  HOFold.csvrows:  289  cols:  7
Read  Managers.csv  rows:  3238  cols:  10
Read  ManagersHalf.csv  rows:  93  cols:  10
Read  Master.csvrows:  17674  cols:  33
Read  Pitching.csv  rows:  40432  cols:  30
Read  PitchingPost.csv  rows:  4284  cols:  30
Read  Salaries.csv  rows:  21464  cols:  5
Read  Schools.csv   rows:  749  cols:  5
Read  SchoolsPlayers.csvrows:  6147  cols:  4
Read  SeriesPost.csvrows:  256  cols:  9
Read  Teams.csv rows:  2655  cols:  48
Read  TeamsFranchises.csv   rows:  120  cols:  4
Read  TeamsHalf.csv rows:  52  cols:  10
Read  Xref_Stats.csvrows:  2753  cols:  3
 ls()
[1] files i inp   name


--
Michael Friendly Email: friendly AT yorku DOT ca
Professor, Psychology Dept.
York University  Voice: 416 736-5115 x66249 Fax: 416 736-5814
4700 Keele StreetWeb:   http://www.datavis.ca
Toronto, ONT  M3J 1P3 CANADA

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Read a list of files into named R data.frames

2011-09-09 Thread Jean V Adams
Michael Friendly  wrote on 09/09/2011 09:39:27 AM:
 
 I have a collection of .csv files in a directory, and want to read them 
 into R data.frames whose names
 are the same as the file names, without the .csv extension
 
 e.g., from
   (files - list.files(pattern=*.csv))
   [1] Allstar.csv AllstarFull.csv
   [3] Appearances.csv AwardsManagers.csv
   [5] AwardsPlayers.csv   AwardsShareManagers.csv
   [7] AwardsSharePlayers.csv  Batting.csv
   [9] BattingPost.csv Fielding.csv
 [11] FieldingOF.csv  FieldingPost.csv
 [13] HallOfFame.csv  HOFold.csv
 [15] Managers.csvManagersHalf.csv
 [17] Master.csv  Pitching.csv
 [19] PitchingPost.csvSalaries.csv
 [21] Schools.csv SchoolsPlayers.csv
 [23] SeriesPost.csv  Teams.csv
 [25] TeamsFranchises.csv TeamsHalf.csv
 
   Allstar - read.csv(Allstar.csv, header=TRUE)
...
   TeamsHalf - read.csv(TeamsHalf.csv, header=TRUE)
 
 Below is what I tried, which reads all the files, but doesn't create the 

 R objects in the global environment.
 What is missing here?
 
 for (i in 1:length(files)) {
  inp - read.csv(file=files[i], header=TRUE)
  name - sub(.csv, , files[i])
  cat(Read , files[i], \trows: , nrow(inp),  cols: , ncol(inp), 

 \n)
  eval(paste(name, - inp))
 }


Check out the assign() function.

?assign

Something like this should work:

for (i in 1:length(files)) {
 inp - read.csv(file=files[i], header=TRUE)
 name - sub(.csv, , files[i])
assign(name, inp)
}

Jean


 
 Read  Allstar.csv   rows:  4475  cols:  3
 Read  AllstarFull.csv   rows:  4676  cols:  8
 Read  Appearances.csv   rows:  94157  cols:  20
 Read  AwardsManagers.csvrows:  57  cols:  6
 Read  AwardsPlayers.csv rows:  2679  cols:  6
 Read  AwardsShareManagers.csv   rows:  344  cols:  7
 Read  AwardsSharePlayers.csvrows:  6354  cols:  7
 Read  Batting.csv   rows:  93955  cols:  24
 Read  BattingPost.csv   rows:  9840  cols:  22
 Read  Fielding.csv  rows:  160710  cols:  18
 Read  FieldingOF.csvrows:  12028  cols:  6
 Read  FieldingPost.csv  rows:  10458  cols:  17
 Read  HallOfFame.csvrows:  3913  cols:  8
 Read  HOFold.csvrows:  289  cols:  7
 Read  Managers.csv  rows:  3238  cols:  10
 Read  ManagersHalf.csv  rows:  93  cols:  10
 Read  Master.csvrows:  17674  cols:  33
 Read  Pitching.csv  rows:  40432  cols:  30
 Read  PitchingPost.csv  rows:  4284  cols:  30
 Read  Salaries.csv  rows:  21464  cols:  5
 Read  Schools.csv   rows:  749  cols:  5
 Read  SchoolsPlayers.csvrows:  6147  cols:  4
 Read  SeriesPost.csvrows:  256  cols:  9
 Read  Teams.csv rows:  2655  cols:  48
 Read  TeamsFranchises.csv   rows:  120  cols:  4
 Read  TeamsHalf.csv rows:  52  cols:  10
 Read  Xref_Stats.csvrows:  2753  cols:  3
   ls()
 [1] files i inp   name
  
 
 -- 
 Michael Friendly Email: friendly AT yorku DOT ca
 Professor, Psychology Dept.
 York University  Voice: 416 736-5115 x66249 Fax: 416 736-5814
 4700 Keele StreetWeb:   http://www.datavis.ca
 Toronto, ONT  M3J 1P3 CANADA

[[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] Very simple question about list components

2011-09-09 Thread Clemontina Alexander
I have a list 'ans' from the following code:

tt - rnorm(50)
rr - rnorm(50)
ans - lm(rr~tt)

ans[1] is $coefficients, ans[2] is $residuals, ans[3] is
$effects, ... and so on up to ans[12]. Is there an easy way to
display just these names and not the data they contain? I thought I
saw my advisor type ans$ and they were displayed, but when I tried
it, it just goes to the next line waiting for input like this:

 ans$
+

I know the answer is really simple, but I can't find it.
Thanks!
Tina

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

2011-09-09 Thread R. Michael Weylandt
Try names(ans)

Hope this helps,

Michael Weylandt

On Fri, Sep 9, 2011 at 10:01 AM, Clemontina Alexander ckale...@ncsu.eduwrote:

 I have a list 'ans' from the following code:

 tt - rnorm(50)
 rr - rnorm(50)
 ans - lm(rr~tt)

 ans[1] is $coefficients, ans[2] is $residuals, ans[3] is
 $effects, ... and so on up to ans[12]. Is there an easy way to
 display just these names and not the data they contain? I thought I
 saw my advisor type ans$ and they were displayed, but when I tried
 it, it just goes to the next line waiting for input like this:

  ans$
 +

 I know the answer is really simple, but I can't find it.
 Thanks!
 Tina

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

2011-09-09 Thread Jorge I Velez
Try

 names(ans)
 [1] coefficients  residuals effects   rank
 [5] fitted.values assignqrdf.residual
 [9] xlevels   call  terms model

HTH,
Jorge
*
*

On Fri, Sep 9, 2011 at 11:01 AM, Clemontina Alexander  wrote:

 I have a list 'ans' from the following code:

 tt - rnorm(50)
 rr - rnorm(50)
 ans - lm(rr~tt)

 ans[1] is $coefficients, ans[2] is $residuals, ans[3] is
 $effects, ... and so on up to ans[12]. Is there an easy way to
 display just these names and not the data they contain? I thought I
 saw my advisor type ans$ and they were displayed, but when I tried
 it, it just goes to the next line waiting for input like this:

  ans$
 +

 I know the answer is really simple, but I can't find it.
 Thanks!
 Tina

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

2011-09-09 Thread Sarah Goslee
 names(ans)
 [1] coefficients  residuals effects   rank
 [5] fitted.values assignqrdf.residual
 [9] xlevels   call  terms model


And thank you for providing a simple reproducible example.

Sarah

On Fri, Sep 9, 2011 at 11:01 AM, Clemontina Alexander ckale...@ncsu.edu wrote:
 I have a list 'ans' from the following code:

 tt - rnorm(50)
 rr - rnorm(50)
 ans - lm(rr~tt)

 ans[1] is $coefficients, ans[2] is $residuals, ans[3] is
 $effects, ... and so on up to ans[12]. Is there an easy way to
 display just these names and not the data they contain? I thought I
 saw my advisor type ans$ and they were displayed, but when I tried
 it, it just goes to the next line waiting for input like this:

 ans$
 +

 I know the answer is really simple, but I can't find it.
 Thanks!
 Tina

-- 
Sarah Goslee
http://www.functionaldiversity.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] Very simple question about list components

2011-09-09 Thread Steve Lianoglou
Hi Clemontina,

On Fri, Sep 9, 2011 at 11:01 AM, Clemontina Alexander ckale...@ncsu.edu wrote:
 I have a list 'ans' from the following code:

 tt - rnorm(50)
 rr - rnorm(50)
 ans - lm(rr~tt)

 ans[1] is $coefficients, ans[2] is $residuals, ans[3] is
 $effects, ... and so on up to ans[12]. Is there an easy way to
 display just these names and not the data they contain?

In this case, you can simply do:

R names(ans)
 [1] coefficients  residuals effects   rank
 [5] fitted.values assignqrdf.residual
 [9] xlevels   call  terms model

 saw my advisor type ans$ and they were displayed, but when I tried
 it, it just goes to the next line waiting for input like this:

 ans$

If you are in the R prompt and you hit tab twice after `ans$`, it
should show you the names of the elements in the `ans` list. If you
use Rstudio[1], the doube tab hit will give you a handy dropdown
menu list you can choose from.

HTH,
-steve

[1] For what it's worth, I'm not sure how far along in your R-learning
you are, but if I were starting out learning R today, I might pick
Rstudio as my IDE/interface of choice.

-- 
Steve Lianoglou
Graduate Student: Computational Systems Biology
 | Memorial Sloan-Kettering Cancer Center
 | Weill Medical College of Cornell University
Contact Info: http://cbio.mskcc.org/~lianos/contact

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] The elegant way to test if a number is a whole number

2011-09-09 Thread William Dunlap
Which of these functions is best really depends on
why you are interested in knowing whether a number
is integral or not.  Why are people interested in this?

I can think of a few

1) I have a C routine with the prototype
   void func(int *n)
   and to call it with .C() I need to make sure
   the data is stored as a C 'int' [4-byte signed
   integer].  The appropriate test is is.integer(n).

2) I have computed something that ought to be the
   length of an output vector.  E.g.,
n - 7
len - 1/6 * n * (n+1) * (n+2)
len
   [1] 84
length(complex(len)) # what?!
   [1] 83
   Here I wanted the test that len==trunc(len) or len==as.integer(len)
   because the constructor complex(len) converts len to an integer using
   truncation.

3) I am computing an integral that I know should give
   an integral result for certain inputs (e.g, the gamma
   function for nonnegative integral inputs) and I want to
   do some sanity checks on the output.  Then all.equal(x, round(x))
   or abs(x-round(x))tolerance would be the right check.

Why are people interested in this question?

Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com 

 -Original Message-
 From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
 Behalf Of Martin Maechler
 Sent: Friday, September 09, 2011 1:41 AM
 To: Marc Schwartz
 Cc: r-help; Alexander Engelhardt
 Subject: Re: [R] The elegant way to test if a number is a whole number
... [ many lines elided ] ...
  is.whole - function(x)
  is.numeric(x)  floor(x)==x
 
  Are you sure? I thought the test would have been all.equal(x,
  round(x,0) )
 
 Yes, I was *very* sure that  the is.whole() above should be
 exactly what it was:
 
 The purpose was really different.
 There, the discussion started from the fact that
   as.integer(1) gives FALSE in R
 and that some people where making a fuzz about the allegedly
 poor design of S and R in that regard.
 
 The purpose of the above  is.whole()
 was solely to be a
 concise  *alternative* to is.integer(),
 -- which even would work with complex numbers, as Marc
 mentioned...
 and BTW, instead of floor(), trunc(), or ceiling() or even round()
 would have been equally valid.  I had chosen floor() because  in
 mathematics, it is the most frequently used of those four functions.
 
 Your topic below ... about accidental rounding or not with
 floating point arithmetic is an entirely different issue,
 *and* I agree that a solution to your problem should be based on
 the same calculations as all.equal.numeric() .. and hence will
 be considerably more involved if it should work in all border
 cases.
 
 Martin Maechler, ETH Zurich
... [ many more lines elided ] ...

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

2011-09-09 Thread Clemontina Alexander
Thank you all. And thanks for pointing out I have an lm object, not a list.

@ Steve, I knew I wasn't crazy when I saw him type ans$ in the
prompt. I didn't see him hit tab though. So a special thanks to you
for helping me validate my sanity! :-)

Tina



On Fri, Sep 9, 2011 at 11:08 AM, Steve Lianoglou
mailinglist.honey...@gmail.com wrote:
 Hi Clemontina,

 On Fri, Sep 9, 2011 at 11:01 AM, Clemontina Alexander ckale...@ncsu.edu 
 wrote:
 I have a list 'ans' from the following code:

 tt - rnorm(50)
 rr - rnorm(50)
 ans - lm(rr~tt)

 ans[1] is $coefficients, ans[2] is $residuals, ans[3] is
 $effects, ... and so on up to ans[12]. Is there an easy way to
 display just these names and not the data they contain?

 In this case, you can simply do:

 R names(ans)
  [1] coefficients  residuals     effects       rank
  [5] fitted.values assign        qr            df.residual
  [9] xlevels       call          terms         model

 saw my advisor type ans$ and they were displayed, but when I tried
 it, it just goes to the next line waiting for input like this:

 ans$

 If you are in the R prompt and you hit tab twice after `ans$`, it
 should show you the names of the elements in the `ans` list. If you
 use Rstudio[1], the doube tab hit will give you a handy dropdown
 menu list you can choose from.

 HTH,
 -steve

 [1] For what it's worth, I'm not sure how far along in your R-learning
 you are, but if I were starting out learning R today, I might pick
 Rstudio as my IDE/interface of choice.

 --
 Steve Lianoglou
 Graduate Student: Computational Systems Biology
  | Memorial Sloan-Kettering Cancer Center
  | Weill Medical College of Cornell University
 Contact Info: http://cbio.mskcc.org/~lianos/contact




-- 
Clemontina Alexander
Ph.D Student
Department of Statistics
NC State University
Raleigh, NC 27695
Phone: (850) 322-6878
Email: ckale...@ncsu.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] Exception in NeweyWest - Pre-Whitening necessary?

2011-09-09 Thread Simon Zehnder
Hi guyz,

I have run my algorithm in R (see http://pastebin.com/q84Tujfg) and got the 
following error: 

Error in ar.ols(x, aic = aic, order.max = order.max, na.action = na.action,  :
  'order.max' must be  'n.used'

I am pretty sure, that the error comes from the NeweyWest function in line 45, 
as the NeweyWest function uses the ar.ols() function for pre whitening. Does 
anyone has an idea how to circumvent this? Can I just shut off pre whitening?

I am thankful for every suggestion. 

Simon


[[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] Rgl and plotmath symbols (via sprites): a trial

2011-09-09 Thread Duncan Murdoch

On 09/09/2011 11:10 AM, Marius Hofert wrote:

Dear all,

Below is some code where I try to get plotmath symbols in an rgl plot. Duncan
Murdoch kindly suggested to use a sprite for this. As you can see, on can get
it to work, but my knowledge about grid and rgl is too limited to perfectly
solve the problem.
1) As you can see (please rotate the plot a little bit so that (0,0,0) is in 
front),
the quality of the .png seems poor. Can this be improved?


Generally plotting much larger and shrinking is the way to do this.


2) When I print the file to a .pdf, the label is just a black square... I assume
one then has to print to png again... hmmm... not perfect.


rgl.postscript has a lot of limitations.  This one may be unavoidable.

3) Some parts of the z-axis are covered by the label. How can this be fixed?
The problem seems to be that sprites have a radius. Specifying a rectangle would
be better here, but I am not sure how this works (or even if it does).



I think I don't see this in Windows, which shows the black text for the 
label on a transparent background.  Do you see a solid background?  Or 
is it just the positioning?  The tick marks go out at an angle, whereas 
your sprite seems to be in line with the zy plane.  I would increase the 
x value to something like 1.4 to place it properly.  (Just tried that, 
and I think I see the problem you had:  the new sprite is obscured by 
the old one.
You may need to play with the depth_test or depth_mask material 
properties.  Not sure if those are visible on CRAN; if not, try the 
R-forge version.)


The other thing I'd recommend is that before plotting the sprite, you 
run par3d(ignoreExtent=TRUE) (and then set it back to the original value 
afterwards).  This will stop the label from affecting the bounding box 
of the plot.


Duncan Murdoch



Cheers,

Marius





require(rgl)
require(grid)

s- seq(0, 1, length.out=21)
M- function(u) apply(u, 1, min)
u- s
v- s
z- outer(u, v, function(u,v) M(cbind(u,v)))

## create z-axis label
png(zlabel.png, bg=transparent, width=70, height=70)
grid.newpage()
pushViewport(plotViewport(c(1,1,1,1)))
print(grid.text(expression(W(u[1],u[2])==c),rot=90))
popViewport()
dev.off()

## plot
persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=,
 ylab=, zlab=)
axes3d(edges=c('x--','y--','z+-'))
sprites3d(1,-0.4,0.65, radius=0.5, lit=FALSE, textype=alpha, 
texture=zlabel.png)
rgl.postscript(W.pdf, fmt=pdf) # print to file

# par3d(c(userMatrix, zoom, FOV)):
#
# $userMatrix
#[,1]   [,2]   [,3] [,4]
# [1,]  0.7262191 -0.6867201 0.031957750
# [2,]  0.2750087  0.3328032 0.902004660
# [3,] -0.6300603 -0.6462642 0.430542020
# [4,]  0.000  0.000 0.1
#
# $zoom
# [1] 1
#
# $FOV
# [1] 30






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

2011-09-09 Thread Martin Batholdy
Hi,

is there a function or an easy way to convert a variable with continuous values 
into a categorial variable (with x levels)?

here is what I mean:


I want to transform x:

x - c(3.2,  1.5,  6.8,  6.9,  8.5,  9.6,  1.1,  0.6)

into a 'categorial'-variable with four levels so that I get:

[1] 2 2 3 3 4 4 1 1

so each element is converted into its rank- value / categorial-value
(in this example four levels are created).



thanks for any suggestions!

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

2011-09-09 Thread Ben Tupper
Hi,

On Sep 9, 2011, at 11:34 AM, Martin Batholdy wrote:

 Hi,
 
 is there a function or an easy way to convert a variable with continuous 
 values into a categorial variable (with x levels)?
 
 here is what I mean:
 
 
 I want to transform x:
 
 x - c(3.2,  1.5,  6.8,  6.9,  8.5,  9.6,  1.1,  0.6)
 
 into a 'categorial'-variable with four levels so that I get:
 
 [1] 2 2 3 3 4 4 1 1
 
 so each element is converted into its rank- value / categorial-value
 (in this example four levels are created).
 
 

Try findInterval()

 x - c(3.2,  1.5,  6.8,  6.9,  8.5,  9.6,  1.1,  0.6)
 findInterval(x, c(0, 1.5, 6, 8))
[1] 2 2 3 3 4 4 1 1


CHeers,
Ben






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

Ben Tupper
Bigelow Laboratory for Ocean Sciences
180 McKown Point Rd. P.O. Box 475
West Boothbay Harbor, Maine   04575-0475 
http://www.bigelow.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] split variable / create categories

2011-09-09 Thread Jean-Christophe BOUËTTÉ
Hi there,
library(lattice)
equal.count(x,number=4,overlap=0)


JC

2011/9/9 Martin Batholdy batho...@googlemail.com:
 Hi,

 is there a function or an easy way to convert a variable with continuous 
 values into a categorial variable (with x levels)?

 here is what I mean:


 I want to transform x:

 x - c(3.2,  1.5,  6.8,  6.9,  8.5,  9.6,  1.1,  0.6)

 into a 'categorial'-variable with four levels so that I get:

 [1] 2 2 3 3 4 4 1 1

 so each element is converted into its rank- value / categorial-value
 (in this example four levels are created).



 thanks for any suggestions!

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


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


[R] Problems using quantile regression (rq) to model GLD random variables in R

2011-09-09 Thread James Shaw
Everyone:

I am working on a simulation of the efficiencies of regression
estimators when applied to model a specific form of highly skewed
data.  The outcome variable (y) is being simulated from a generalized
lambda distribution (GLD) to reflect the characteristics (mean,
variance, skewness, kurtosis) of an observed variable.  The regressor
of interest (x) is simply a binary indicator of group membership.  The
relevant R code is as follows:

x[1:200]-0
x[201:400]-1

params - c(2.684864244,.0144182,26.01913,711.0554)
y0-rgl(200, params, param=rs)

params - c(-0.113741088,.0523072,15.98282,426.4456)
y1-rgl(200, params, param=rs)

y-c(y0, y1)

I have verified that the GLD parameters in each case are valid using
gl.check.lambda (in the GLD package).

While I experienced no difficulty when using OLS to fit models to y,
the quantile regression estimator and robust (e.g., M) regression
estimators yielded minute (or missing) variance estimates and
infinitely large t statistics for the coefficient for x.  The problem
appears to be related to the number of duplicate observations in my
simulated data.  As I understand it, the GLD is a transformation of
the uniform distribution.  Given the parameters specified above, I end
up with many duplicate observations that happen to be equal to the
true median value.  This lack of variation around the median appears
to be causing problems for the quantile regression estimator (as
implemented using rq) and robust regression estimator.

I am unaware of a viable alternative to the GLD that can be readily
implemented in R.  In the absence of an alternative distribution, I am
wondering whether jittering could be used as a practical (and
hopefully valid) solution to my dilemma.  That is, add a small
residual drawn from U(-.5,.5) to each GLD observation and model the
composite variable as a function of x.  This would be expected to
preserve the mean and median of y over repeated simulations, and the
added variance would be expected to be negligible.  When using this
procedure, I derive reasonable variance estimates and get results that
make intuitive sense (i.e., the efficiency of the M estimator =
quantile regression  OLS).  I have seen a similar jittering procedure
applied in a paper on the modeling of quantiles of count data (Machado
and Santa Silva. JASA. 2005; 100: 1226).

I would appreciate others' thoughts regarding the validity of the
proposed jittering procedure or suggestions for alternative approaches
I could use to deal with my problem.  Many thanks!

Regards,

Jim


P.S.:  Although I do not think it has any bearing on my problem, here
is the quantile regression code I am using:

fit1-summary(rq(y~x, tau = .5, ci=FALSE),se=ker)

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Exception in NeweyWest - Pre-Whitening necessary?

2011-09-09 Thread Achim Zeileis

On Fri, 9 Sep 2011, Simon Zehnder wrote:


Hi guyz,

I have run my algorithm in R (see http://pastebin.com/q84Tujfg) and got 
the following error:


This is not reproducible. There is only the script but no information 
about how it is called. Also, within your script you call 
NeweyWest(regression) -- however, regression is not defined within the 
script!



Error in ar.ols(x, aic = aic, order.max = order.max, na.action = na.action,  :
 'order.max' must be  'n.used'


NeweyWest call ar() essentially with order.max = as.integer(prewhite) 
which defaults to 1. Thus, the error must be caused by a model with
n.used = 1 which suggests that it is not a very sensible model to start 
with.


I am pretty sure, that the error comes from the NeweyWest function in 
line 45, as the NeweyWest function uses the ar.ols() function for pre 
whitening. Does anyone has an idea how to circumvent this? Can I just 
shut off pre whitening?


Why don't you read the documentation in ?NeweyWest or in 
vignette(sandwich, package = sandwich)?

Z


I am thankful for every suggestion.

Simon


[[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] get mean from cdf

2011-09-09 Thread annie Zhang
Thanks all for your replies. Your responses are very helpful. But it is
indeed not my homework! I didn't know these functions, and wanted to avoid
the loop.

Annie

On Fri, Sep 9, 2011 at 1:38 AM, peter dalgaard pda...@gmail.com wrote:


 On Sep 9, 2011, at 03:51 , Jeff Newmiller wrote:

  The diff function would be helpful.


 ...or not:

  sum(ppois(0:1000, lambda=2.345, lower=F))
 [1] 2.345

 (if this was homework: the hard bit is to figure out _why_ this works.)

  annie Zhang annie.zhang2...@gmail.com wrote:
 
  Hi All,
 
  How can I get the expected value from a discrete cdf? Is there any R
  function that can do this?
 
  Thanks,
  Annie

 --
 Peter Dalgaard, Professor,
 Center for Statistics, Copenhagen Business School
 Solbjerg Plads 3, 2000 Frederiksberg, Denmark
 Phone: (+45)38153501
 Email: pd@cbs.dk  Priv: pda...@gmail.com
 Døden skal tape! --- Nordahl Grieg









[[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] split variable / create categories

2011-09-09 Thread Martin Batholdy
Thanks for the suggestions!

However all these functions don't produce exactly what I want
(at least with my actual data).


I need a split-algorithm that converts the values of my vectors into four 
factors.
And the crucial part is, that I need exactly the same number of elements in 
each factor-level
and no overlapping.



cut() seems to find equal intervals – but that leads to different numbers of 
elements in each interval.


library(lattice)
equal.count(x,number=4,overlap=0)

seems to do the job, but strangely enough, it seems to ignore the argument 
'overlap = 0' in my actual vector –
I get factor-borders that overlap.
And I really have to prevent this.




On 09.09.2011, at 17:49, Andrea Spano wrote:

 cut ( x , c(0, 1.4 ,6, 8, Inf ), labels = 1:4, include.lowest = T)
 
 On 9 September 2011 17:34, Martin Batholdy batho...@googlemail.com wrote:
 Hi,
 
 is there a function or an easy way to convert a variable with continuous 
 values into a categorial variable (with x levels)?
 
 here is what I mean:
 
 
 I want to transform x:
 
 x - c(3.2,  1.5,  6.8,  6.9,  8.5,  9.6,  1.1,  0.6)
 
 into a 'categorial'-variable with four levels so that I get:
 
 [1] 2 2 3 3 4 4 1 1
 
 so each element is converted into its rank- value / categorial-value
 (in this example four levels are created).
 
 
 
 thanks for any suggestions!
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.
 
 
 
 -- 
 Andrea Spano'
 
 Quantide s.r.l
 +39 347 747 04 92.
 andrea.sp...@quantide.com
 
 
 This electronic mail transmission may contain confidential information
 addressed only to the person (s) named. Any use, distribution, copyng or
 disclosure by any other person and/or entities other than the intended
 recipient is prohibited. If you received this transmission in error,
 please inform the sender immediately and delete the material.
 


[[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] How does one start R within Emacs/ESS with root privileges?

2011-09-09 Thread Karl Brand

Thanks a lot for sharing your experience and suggestions.

I see now the logic of making use of my user library which doesn't 
require root privs. Apparently linux is actually set up to work well - 
something i'm not used to from the other OS i've used till now.


And for the rare times i do need to install into the system-wide 
library, i think i'll survive in the terminal :P


cheers,

Karl



On 2011-09-07 17:08, Spencer Graves wrote:

Under Vista and Windows 7, I install R in a local directory pgms so I
never have to worry about permissions.


Under Linux, I use su not sudo. Then I call R and install.packages.
Then I quit R and su and continue with what I want.


hope this helps.
Spencer


On 9/7/2011 5:11 AM, Karl Brand wrote:

Cheers Paul.

Its a very good point. Although i am curious how badly i can damage my
R install by running as root. I always ran R in windows with admin.
privileges without problems (touch wood). Probably best to never find
out by sticking with user privileges.

However, even for taking care of R install/maint. i'd prefer to do
this interactively within Emacs rather than the terminal. Motivated by
this, i'd still like to find out how to invoke R with root privileges.

I've also reposted the original email on perhaps a more appropriate
forum at: ess-h...@stat.math.ethz.ch

Karl


On 2011-09-07 11:02, Paul Hiemstra wrote:

On 09/07/2011 08:54 AM, Karl Brand wrote:

Esteemed UseRs and DevelopeRs,

Apologies if this question belongs else where, but it does concern R's
package installation/maintenance.

How does one start R within Emacs/ESS with root privileges?

I tried without success:


M-x sudo R


Why i'm motivated to do so:

It seems logical to me, as the only user of the PC, to keep my R
library consolidated in the universal library rather than splitting
into universal and user libraries. Hence the desire to run R as root.


Hi Karl,

Why the need to install packages in root? As you are the only user there
is not reason to install them system wide (to make them available to all
users, which is just you). Installing the packages in your homedir
solves your problem much easier, without the need to run R as root
continuously. I think you should not run anything as root if it is not
absolutely needed, which could potentially damage your system
(accidentally overwriting something).

hope this helps,
Paul



In addition, it's nice to be able to install packages 'on the fly'
when and as needed and not need to launch a separate R session (as
root) in the terminal just to install a package.

Migrating from windows, i'm completey new to linux (ubuntu) and am
seeing for myself if Emacs/ESS is as good as its purported to be. So
maybe my motivation is nonsensical to expereinced ESS/R users. If so
i'd really appreciate tips on efficient package
installation/maintenance using Emacs/ESS.

TIA,

karl







--
Karl Brand k.br...@erasmusmc.nl
Department of Genetics
Erasmus MC
Dr Molewaterplein 50
3015 GE Rotterdam
P +31 (0)10 704 3455 | F +31 (0)10 704 4743 | M +31 (0)642 777 268

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Comments inside a line? Not possible?

2011-09-09 Thread Dimitri Liakhovitski
Hello!
In the guidelines I've read:

Comment your code. Entire commented lines should begin with # and one space.
Short comments can be placed after code preceded by two spaces, #, and
then one space. 

Just wondering if there is no way to comment something out in the
middle of a line? Example:
Original code:
mystrings-c(a,b,c)

Desired commented code:
mystrings-c(a, # want b to be commented out, but not c # c)
So that it is read as: mystrings-c(a,c)

Not possible?
Thanks!

-- 
Dimitri Liakhovitski
marketfusionanalytics.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] Rgl and plotmath symbols (via sprites): a trial

2011-09-09 Thread Marius Hofert
Dear Duncan,

attached my next trial. It seems to work :-) Unfortunately, not the printing to 
pdf. Not getting a high-quality pdf is certainly one of the major drawbacks. If 
anybody knows how to fix that, please let me know. Even when generating the 
labels by hand (with grid), the code is less than setting up a nice lattice 
wireframe plot... so that makes rgl a great tool. 

Cheers,

Marius


require(rgl)
require(grid)

s - seq(0, 1, length.out=21)
M - function(u) apply(u, 1, min)
u - s
v - s
z - outer(u, v, function(u,v) M(cbind(u,v)))

## create z-axis label
png(zlabel.png, bg=transparent, width=500, height=500)
grid.newpage()
pushViewport(plotViewport(c(1,1,1,1)))
print(grid.text(expression(W(u[1],u[2])==c), rot=90, gp=gpar(fontsize=85)))
popViewport()
dev.off()

## plot
persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=, 
ylab=, zlab=) 
axes3d(edges=c('x--','y--','z+-')) 
par3d(ignoreExtent=TRUE)
sprites3d(1.4,-0.2,0.5, radius=0.5, lit=FALSE, textype=alpha, 
texture=zlabel.png)
par3d(windowRect=c(0,0,600,600), zoom=1.6) 
rgl.snapshot(W.png, fmt=png) # print to file



On 2011-09-09, at 17:26 , Duncan Murdoch wrote:

 On 09/09/2011 11:10 AM, Marius Hofert wrote:
 Dear all,
 
 Below is some code where I try to get plotmath symbols in an rgl plot. Duncan
 Murdoch kindly suggested to use a sprite for this. As you can see, on can 
 get
 it to work, but my knowledge about grid and rgl is too limited to perfectly
 solve the problem.
 1) As you can see (please rotate the plot a little bit so that (0,0,0) is 
 in front),
 the quality of the .png seems poor. Can this be improved?
 
 Generally plotting much larger and shrinking is the way to do this.

 2) When I print the file to a .pdf, the label is just a black square... I 
 assume
 one then has to print to png again... hmmm... not perfect.
 
 rgl.postscript has a lot of limitations.  This one may be unavoidable.
 3) Some parts of the z-axis are covered by the label. How can this be fixed?
 The problem seems to be that sprites have a radius. Specifying a rectangle 
 would
 be better here, but I am not sure how this works (or even if it does).
 
 
 I think I don't see this in Windows, which shows the black text for the label 
 on a transparent background.  Do you see a solid background?  Or is it just 
 the positioning?  The tick marks go out at an angle, whereas your sprite 
 seems to be in line with the zy plane.  I would increase the x value to 
 something like 1.4 to place it properly.  (Just tried that, and I think I see 
 the problem you had:  the new sprite is obscured by the old one.
 You may need to play with the depth_test or depth_mask material properties.  
 Not sure if those are visible on CRAN; if not, try the R-forge version.)
 
 The other thing I'd recommend is that before plotting the sprite, you run 
 par3d(ignoreExtent=TRUE) (and then set it back to the original value 
 afterwards).  This will stop the label from affecting the bounding box of the 
 plot.
 
 Duncan Murdoch
 
 
 Cheers,
 
 Marius
 
 
 
 
 
 require(rgl)
 require(grid)
 
 s- seq(0, 1, length.out=21)
 M- function(u) apply(u, 1, min)
 u- s
 v- s
 z- outer(u, v, function(u,v) M(cbind(u,v)))
 
 ## create z-axis label
 png(zlabel.png, bg=transparent, width=70, height=70)
 grid.newpage()
 pushViewport(plotViewport(c(1,1,1,1)))
 print(grid.text(expression(W(u[1],u[2])==c),rot=90))
 popViewport()
 dev.off()
 
 ## plot
 persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=,
 ylab=, zlab=)
 axes3d(edges=c('x--','y--','z+-'))
 sprites3d(1,-0.4,0.65, radius=0.5, lit=FALSE, textype=alpha, 
 texture=zlabel.png)
 rgl.postscript(W.pdf, fmt=pdf) # print to file
 
 # par3d(c(userMatrix, zoom, FOV)):
 #
 # $userMatrix
 #[,1]   [,2]   [,3] [,4]
 # [1,]  0.7262191 -0.6867201 0.031957750
 # [2,]  0.2750087  0.3328032 0.902004660
 # [3,] -0.6300603 -0.6462642 0.430542020
 # [4,]  0.000  0.000 0.1
 #
 # $zoom
 # [1] 1
 #
 # $FOV
 # [1] 30
 
 
 
 
 

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

2011-09-09 Thread Jean-Christophe BOUËTTÉ
so what about :
library(ggplot2)
z - breaks(x,equal=NULL,nbins=4)
a - cut(x,z,include.lowest = TRUE)
it works for me !
JC

2011/9/9 Martin Batholdy batho...@googlemail.com:
 Thanks for the suggestions!

 However all these functions don't produce exactly what I want
 (at least with my actual data).


 I need a split-algorithm that converts the values of my vectors into four 
 factors.
 And the crucial part is, that I need exactly the same number of elements in 
 each factor-level
 and no overlapping.



 cut() seems to find equal intervals – but that leads to different numbers of 
 elements in each interval.


 library(lattice)
 equal.count(x,number=4,overlap=0)

 seems to do the job, but strangely enough, it seems to ignore the argument 
 'overlap = 0' in my actual vector –
 I get factor-borders that overlap.
 And I really have to prevent this.




 On 09.09.2011, at 17:49, Andrea Spano wrote:

 cut ( x , c(0, 1.4 ,6, 8, Inf ), labels = 1:4, include.lowest = T)

 On 9 September 2011 17:34, Martin Batholdy batho...@googlemail.com wrote:
 Hi,

 is there a function or an easy way to convert a variable with continuous 
 values into a categorial variable (with x levels)?

 here is what I mean:


 I want to transform x:

 x - c(3.2,  1.5,  6.8,  6.9,  8.5,  9.6,  1.1,  0.6)

 into a 'categorial'-variable with four levels so that I get:

 [1] 2 2 3 3 4 4 1 1

 so each element is converted into its rank- value / categorial-value
 (in this example four levels are created).



 thanks for any suggestions!

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



 --
 Andrea Spano'

 Quantide s.r.l
 +39 347 747 04 92.
 andrea.sp...@quantide.com


 This electronic mail transmission may contain confidential information
 addressed only to the person (s) named. Any use, distribution, copyng or
 disclosure by any other person and/or entities other than the intended
 recipient is prohibited. If you received this transmission in error,
 please inform the sender immediately and delete the material.



        [[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] Comments inside a line? Not possible?

2011-09-09 Thread David Winsemius


On Sep 9, 2011, at 12:21 PM, Dimitri Liakhovitski wrote:


Hello!
In the guidelines I've read:

Comment your code. Entire commented lines should begin with # and  
one space.

Short comments can be placed after code preceded by two spaces, #, and
then one space. 

Just wondering if there is no way to comment something out in the
middle of a line? Example:
Original code:
mystrings-c(a,b,c)

Desired commented code:
mystrings-c(a, # want b to be commented out, but not c # c)
So that it is read as: mystrings-c(a,c)

Not possible?


Not possible. Comments are ended by newlines. And why would you not  
want:


mystrings-c(a, # want b to be commented out, but not c
 c)

--

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] Comments inside a line? Not possible?

2011-09-09 Thread Dimitri Liakhovitski
Thank you, David.
The only reason why I would prefer the in-line commenting is that it
would save space.
Having to start a new line stretches the code downwards.
Dimitri

On Fri, Sep 9, 2011 at 12:36 PM, David Winsemius dwinsem...@comcast.net wrote:

 On Sep 9, 2011, at 12:21 PM, Dimitri Liakhovitski wrote:

 Hello!
 In the guidelines I've read:

 Comment your code. Entire commented lines should begin with # and one
 space.
 Short comments can be placed after code preceded by two spaces, #, and
 then one space. 

 Just wondering if there is no way to comment something out in the
 middle of a line? Example:
 Original code:
 mystrings-c(a,b,c)

 Desired commented code:
 mystrings-c(a, # want b to be commented out, but not c # c)
 So that it is read as: mystrings-c(a,c)

 Not possible?

 Not possible. Comments are ended by newlines. And why would you not want:

 mystrings-c(a, # want b to be commented out, but not c
             c)

 --

 David Winsemius, MD
 West Hartford, CT





-- 
Dimitri Liakhovitski
marketfusionanalytics.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] Global optimisation with inequality constraints

2011-09-09 Thread Ravi Varadhan
Hi,

Take a look at the multiStart() function in the BB package.  This allows you 
to specify random multiple starting values, which need not satisfy constraints. 
 Then you need to specify the constraints using the `projectLinear' argument.

Ravi.

---
Ravi Varadhan, Ph.D.
Assistant Professor,
Division of Geriatric Medicine and Gerontology School of Medicine Johns Hopkins 
University

Ph. (410) 502-2619
email: rvarad...@jhmi.edumailto:rvarad...@jhmi.edu


[[alternative HTML version deleted]]

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


Re: [R] split variable / create categories

2011-09-09 Thread David Winsemius


On Sep 9, 2011, at 12:10 PM, Martin Batholdy wrote:


Thanks for the suggestions!

However all these functions don't produce exactly what I want
(at least with my actual data).


I need a split-algorithm that converts the values of my vectors into  
four factors.
And the crucial part is, that I need exactly the same number of  
elements in each factor-level

and no overlapping.



cut() seems to find equal intervals – but that leads to different  
numbers of elements in each interval.




You may want to look at 'cut2' with the 'g' argument in the Hmisc  
package. It's defaults are to include.lowest and it give the option of  
specifying equal sized groups. It fits with my assumptions about how a  
factor-ing variable _should_ be constructed, so apparently Harrell and  
I think alike.


(See code and output below. You displayed output from a numeric vector  
rather than a factor.)


--
David.


library(lattice)
equal.count(x,number=4,overlap=0)

seems to do the job, but strangely enough, it seems to ignore the  
argument 'overlap = 0' in my actual vector –

I get factor-borders that overlap.
And I really have to prevent this.




On 09.09.2011, at 17:49, Andrea Spano wrote:


cut ( x , c(0, 1.4 ,6, 8, Inf ), labels = 1:4, include.lowest = T)

On 9 September 2011 17:34, Martin Batholdy  
batho...@googlemail.com wrote:

Hi,

is there a function or an easy way to convert a variable with  
continuous values into a categorial variable (with x levels)?


here is what I mean:


I want to transform x:

x - c(3.2,  1.5,  6.8,  6.9,  8.5,  9.6,  1.1,  0.6)

into a 'categorial'-variable with four levels so that I get:

[1] 2 2 3 3 4 4 1 1


 as.numeric(cut2 (x - c(3.2,  1.5,  6.8,  6.9,  8.5,  9.6,  1.1,   
0.6) ,g=4) )

[1] 2 2 3 3 4 4 1 1




so each element is converted into its rank- value / categorial-value
(in this example four levels are created).




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] split variable / create categories

2011-09-09 Thread David L Carlson
This will work as long as the number of values is evenly divisible by four
and no two observations are identical. But if either of those assumptions is
false, you cannot be certain to get the even split you are looking for

x - c(3.2,  1.5,  6.8,  6.9,  8.5,  9.6,  1.1,  0.6)

ceiling(rank(x)/(length(x)/4))

[1] 2 2 3 3 4 4 1 1

--
David L Carlson
Associate Professor of Anthropology
Texas AM University
College Station, TX 77843-4352



-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
Behalf Of Jean-Christophe BOUËTTÉ
Sent: Friday, September 09, 2011 11:21 AM
To: Martin Batholdy
Cc: R Help
Subject: Re: [R] split variable / create categories

so what about :
library(ggplot2)
z - breaks(x,equal=NULL,nbins=4)
a - cut(x,z,include.lowest = TRUE)
it works for me !
JC

2011/9/9 Martin Batholdy batho...@googlemail.com:
 Thanks for the suggestions!

 However all these functions don't produce exactly what I want
 (at least with my actual data).


 I need a split-algorithm that converts the values of my vectors into four
factors.
 And the crucial part is, that I need exactly the same number of elements
in each factor-level
 and no overlapping.



 cut() seems to find equal intervals – but that leads to different numbers
of elements in each interval.


 library(lattice)
 equal.count(x,number=4,overlap=0)

 seems to do the job, but strangely enough, it seems to ignore the argument
'overlap = 0' in my actual vector –
 I get factor-borders that overlap.
 And I really have to prevent this.




 On 09.09.2011, at 17:49, Andrea Spano wrote:

 cut ( x , c(0, 1.4 ,6, 8, Inf ), labels = 1:4, include.lowest = T)

 On 9 September 2011 17:34, Martin Batholdy batho...@googlemail.com
wrote:
 Hi,

 is there a function or an easy way to convert a variable with continuous
values into a categorial variable (with x levels)?

 here is what I mean:


 I want to transform x:

 x - c(3.2,  1.5,  6.8,  6.9,  8.5,  9.6,  1.1,  0.6)

 into a 'categorial'-variable with four levels so that I get:

 [1] 2 2 3 3 4 4 1 1

 so each element is converted into its rank- value / categorial-value
 (in this example four levels are created).



 thanks for any suggestions!

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



 --
 Andrea Spano'

 Quantide s.r.l
 +39 347 747 04 92.
 andrea.sp...@quantide.com


 This electronic mail transmission may contain confidential information
 addressed only to the person (s) named. Any use, distribution, copyng or
 disclosure by any other person and/or entities other than the intended
 recipient is prohibited. If you received this transmission in error,
 please inform the sender immediately and delete the material.



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

2011-09-09 Thread Ben Bolker
Eric Tiger footballinlsu at yahoo.com writes:

 
 Hello,
 
 I got one question on the Latin hypercube sampling.
 
 suppose
  there are three variables a, b, c, all of them follow the normal 
 distribution. the mean value and standard deviation for each are  a(32, 
 2), b(35,5), c(37,3). I would like to use Latin hypercube sampling to 
 random generate 1000 samples. but it needs to satisfy the condition that
  abc.
 
 How can I implement this conditional sampling?
 

   Inefficient, but perhaps you could just sample a lot more than
you need and retain the ones that satisfy your criterion?  (That
would be efficient in terms of your time, and I can't immediately
think of any reason that the result would be biased ...) If all
the variables had the same mean you'd need about 8000 points total
(since the ordering relationship you have here would apply 1/(2^3)
of the time), it should be even easier in your case since the means
line up in the order you want.

  Just a thought.

  Ben Bolker

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] reshape data from long to wide format

2011-09-09 Thread Dennis Murphy
Hi:

There are two ways to go about this, depending on what you want.

# 1: Jean's output format:

# I'm using the reshape package here with cast(),
# but it should work similarly with dcast() in reshape2:
 cast(example, DATE ~ SENSOR, value = 'VALUE')
   DATE   A B  C  D E   F
1  01/01/2010 1 270NA NA NANA  NA
2  01/01/2010 2  NA 292.5 NA NANA  NA
3  01/01/2010 3  NANA  0 NANA  NA
4  01/01/2010 4  NANA NA 45NA  NA
5  01/01/2010 5  NANA NA NA 247.5  NA
6  01/01/2010 6  NANA NA NANA 315

# Method 2: One line per day - requires redefinition of date
example$date - as.Date('01/01/2010', format = '%m/%d/%Y')
cast(example, date ~ SENSOR, value = 'VALUE')
date   A B C  D E   F
1 2010-01-01 270 292.5 0 45 247.5 315

In both cases, note the use of value = 'VALUE' in the cast() call. It
needs value_var = 'VALUE' instead of value = 'VALUE' if you're using
dcast() from reshape2 instead:

library('reshape2')
dcast(example, DATE ~ SENSOR, value_var = 'VALUE')
dcast(example, date ~ SENSOR, value_var = 'VALUE')

HTH,
Dennis

On Fri, Sep 9, 2011 at 4:28 AM, maxbre mbres...@arpa.veneto.it wrote:
 This is my reproducible example:

 example-structure(list(SENSOR = structure(1:6, .Label = c(A, B, C,
 D, E, F), class = factor), VALUE = c(270, 292.5, 0, 45,
 247.5, 315), DATE = structure(1:6, .Label = c( 01/01/2010 1,
  01/01/2010 2,  01/01/2010 3,  01/01/2010 4,  01/01/2010 5,
  01/01/2010 6), class = factor)), .Names = c(SENSOR, VALUE,
 DATE), class = data.frame, row.names = c(1, 2, 3, 4,
 5, 6))

 I need to resahpe example  in a wide format so that “SENSOR” appear as
 columns and “DATE” as rows with corresponding “VALUE” as value;

 I thought it was very simple so that I’ve been trying this:

 dcast(example,DATE~SENSOR)

 But I've got this message:

 Using DATE as value column.  Use the value argument to cast to override this
 choice
 Errore in `[.data.frame`(data, , variables, drop = FALSE) :
  undefined columns selected

 and even if by using the value argument sorted out any good result…

 sorry for the very trivial question, I've been looking at documentation and
 forums anywhere but I was not successful at all and now I’m somehow in the
 right middle of nowhere…

 any help or hint for this?

 thank you

 max

 --
 View this message in context: 
 http://r.789695.n4.nabble.com/reshape-data-from-long-to-wide-format-tp3801381p3801381.html
 Sent from the R help mailing list archive at Nabble.com.

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


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


Re: [R] suggestion for proportions

2011-09-09 Thread array chip
Thanks all again for the suggestions. I agree with Wolfgang that mcnemar.test() 
is what I am looking for. The accuracy is the proportion of correct diagnosis 
compared to a gold standard, and I am interested in which diagnosis test is 
better, not particular interested in assessing the agreement between the 2 
tests.

Thanks again.


John


- Original Message -
From: Viechtbauer Wolfgang (STAT) wolfgang.viechtba...@maastrichtuniversity.nl
To: r-help@r-project.org r-help@r-project.org
Cc: csrabak cesar.ra...@gmail.com; array chip arrayprof...@yahoo.com
Sent: Thursday, September 8, 2011 1:24 AM
Subject: RE: [R] suggestion for proportions

I assume you mean Cohen's kappa. This is not what the OP is asking about. The 
OP wants to know how to test for a difference in the proportions of 1's. 
Cohen's kappa will tell you what the level of agreement is between the two 
tests. This is something different.

Also, the OP has now clarified that the data are paired. Therefore, prop.test() 
and binom.test() are not appropriate. So, to answer the OPs question: yes, 
mcnemar.test() is what you should be using.

Wolfgang

 -Original Message-
 From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org]
 On Behalf Of csrabak
 Sent: Thursday, September 08, 2011 02:31
 To: array chip
 Cc: r-help@r-project.org
 Subject: Re: [R] suggestion for proportions
 
 Em 7/9/2011 16:53, array chip escreveu:
  Hi all, thanks very much for sharing your thoughts. and sorry for my
 describing the problem not clearly, my fault.
 
  My data is paired, that is 2 different diagnostic tests were performed
 on the same individuals. Each individual will have a test results from
 each of the 2 tests. Then in the end, 2 accuracy rates were calculated for
 the 2 tests. And I want to test if there is a significant difference in
 the accuracy (proportion) between the 2 tests. My understanding is that
 prop.test() is appropriate for 2 independent proportions,  whereas in my
 situation, the 2 proportions are not independent calculated from paired
 data, right?
 
  the data would look like:
 
  pid   test1    test2
  p1      1         0
  p2      1         1
  p3      0         1
  :
  :
 
  1=test is correct; 0=not correct
 
  from the data above, we can calculate accuracy for test1 and test2, then
 to compare
 
 
  So mcnemar.test() is good for that, right?
 
  Thanks
 
 John,
 
  From above clarifying I suggest you consider the use of kappa test. For
 a list of possible ways of doing it in R try:
 RSiteSearch(kappa,restrict=functions)
 
 HTH
 
 --
 Cesar Rabak

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

2011-09-09 Thread Weinfeld, Jillian
To whom it may concern:

My name is Jillian Weinfeld. I am currently and undergraduate student at New 
York University and working at Mount Sinai School of Medicine doing research 
with epilepsy patients.

At the moment I am creating a manhattan plot with my data set. After reading 
many forums and such, I have appropriately plotted my data, however, I wanted 
to see how I can change the colors of the data points for my plot. The normal 
function that I use for color is not working. I have tried many different 
arguments, but none seem to work. How should I go about changing the colors of 
my data points?

Thank you very much for your assistance.

Best,
Jillian Weinfeld



Jillian Weinfeld

Department of Psychiatry

One Gustave L. Levy Place Box 1230

New York, NY 10029-6574

Phone (212) 659-5629

Email: jillian.weinf...@mssm.edu

[[alternative HTML version deleted]]

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


[R] simple t test for multiple colimns

2011-09-09 Thread Mihovil Pletikos
hi
i am new to r and i would like to how could i do something like this in
spss:

T-TEST GROUPS=who('1' '2')
 /MISSING=ANALYSIS
 /VARIABLES= X3881686 X2672712 X2842570 X3526544 X2902531 X2402942 X3382216
X3771800 X2427469 X2392945 X2453006 X2562821 X3416651 X3552083 X3721851
X2477438 X2926969 X3026969 X3442176 X3796335 X2622772 X3636470 X2952497
X2367537 X3332276 X2427500 X2647315
/CRITERIA=CI(.95).

so i wouldn't have to enter the name of columns by hand, but automaticly
... please help me

Mihovil

[[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] Manhattan Plot

2011-09-09 Thread R. Michael Weylandt
You need to supply some information about how you actually are plotting your
data so we can figure out what you've tried. Example code (with some toy
data) would be greatly appreciated.

Michael Weylandt

On Fri, Sep 9, 2011 at 11:40 AM, Weinfeld, Jillian 
jillian.weinf...@mssm.edu wrote:

 To whom it may concern:

 My name is Jillian Weinfeld. I am currently and undergraduate student at
 New York University and working at Mount Sinai School of Medicine doing
 research with epilepsy patients.

 At the moment I am creating a manhattan plot with my data set. After
 reading many forums and such, I have appropriately plotted my data, however,
 I wanted to see how I can change the colors of the data points for my plot.
 The normal function that I use for color is not working. I have tried many
 different arguments, but none seem to work. How should I go about changing
 the colors of my data points?

 Thank you very much for your assistance.

 Best,
 Jillian Weinfeld

 

 Jillian Weinfeld

 Department of Psychiatry

 One Gustave L. Levy Place Box 1230

 New York, NY 10029-6574

 Phone (212) 659-5629

 Email: jillian.weinf...@mssm.edu

[[alternative HTML version deleted]]

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


[[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] Comments inside a line? Not possible?

2011-09-09 Thread Duncan Murdoch

On 09/09/2011 12:21 PM, Dimitri Liakhovitski wrote:

Hello!
In the guidelines I've read:

Comment your code. Entire commented lines should begin with # and one space.
Short comments can be placed after code preceded by two spaces, #, and
then one space. 


Which guidelines are those?  They seem reasonable, but they put more 
restrictions on comments than the parser enforces.  In particular, the 
spaces are irrelevant.  If you have a # symbol in the line, not in a 
quoted string, it starts a comment that runs to the end of the line.


Duncan Murdoch

Just wondering if there is no way to comment something out in the
middle of a line? Example:
Original code:
mystrings-c(a,b,c)

Desired commented code:
mystrings-c(a, # want b to be commented out, but not c # c)
So that it is read as: mystrings-c(a,c)

Not possible?
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] simple t test for multiple colimns

2011-09-09 Thread Jean-Christophe BOUËTTÉ
The examples in ?t.test and ?formula will help you.
Don't forget An Introduction to R too to import your data in a data
frame first...

2011/9/9 Mihovil Pletikos mihovil.pleti...@gmail.com:
 hi
 i am new to r and i would like to how could i do something like this in
 spss:

 T-TEST GROUPS=who('1' '2')
  /MISSING=ANALYSIS
  /VARIABLES= X3881686 X2672712 X2842570 X3526544 X2902531 X2402942 X3382216
 X3771800 X2427469 X2392945 X2453006 X2562821 X3416651 X3552083 X3721851
 X2477438 X2926969 X3026969 X3442176 X3796335 X2622772 X3636470 X2952497
 X2367537 X3332276 X2427500 X2647315
 /CRITERIA=CI(.95).

 so i wouldn't have to enter the name of columns by hand, but automaticly
 ... please help me

 Mihovil

        [[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] Manhattan Plot

2011-09-09 Thread Jorge I Velez
Hi Jillian,

Take a look at http://www.oga-lab.net/RGM2/func.php?rd_id=gap:mhtplot.

Best,
Jorge
*
*

On Fri, Sep 9, 2011 at 12:40 PM, Weinfeld, Jillian  wrote:

 To whom it may concern:

 My name is Jillian Weinfeld. I am currently and undergraduate student at
 New York University and working at Mount Sinai School of Medicine doing
 research with epilepsy patients.

 At the moment I am creating a manhattan plot with my data set. After
 reading many forums and such, I have appropriately plotted my data, however,
 I wanted to see how I can change the colors of the data points for my plot.
 The normal function that I use for color is not working. I have tried many
 different arguments, but none seem to work. How should I go about changing
 the colors of my data points?

 Thank you very much for your assistance.

 Best,
 Jillian Weinfeld

 

 Jillian Weinfeld

 Department of Psychiatry

 One Gustave L. Levy Place Box 1230

 New York, NY 10029-6574

 Phone (212) 659-5629

 Email: jillian.weinf...@mssm.edu

[[alternative HTML version deleted]]

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


[[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] Rgl and plotmath symbols (via sprites): a trial

2011-09-09 Thread Duncan Murdoch

On 09/09/2011 12:19 PM, Marius Hofert wrote:

Dear Duncan,

attached my next trial. It seems to work :-) Unfortunately, not the printing to pdf. Not 
getting a high-quality pdf is certainly one of the major drawbacks. If anybody knows how 
to fix that, please let me know. Even when generating the labels by hand 
(with grid), the code is less than setting up a nice lattice wireframe plot... so that 
makes rgl a great tool.


Two suggestions:

1. remember the par3d(ignoreExtent=FALSE) call to restore the default 
behaviour, or you'll be unhappy when you draw your next plot.


2. Take a look at the mouseOneAxis mouse callback in 
demo(mouseCallbacks).  It would make sense in your example to keep the z 
axis always pointing up; that callback will do it.


And one todo item for me:  it appears that the sprites aren't being 
recognized as partly transparent, so are being drawn in the wrong 
order.  Got to fix that.


Duncan Murdoch

Cheers,

Marius


require(rgl)
require(grid)

s- seq(0, 1, length.out=21)
M- function(u) apply(u, 1, min)
u- s
v- s
z- outer(u, v, function(u,v) M(cbind(u,v)))

## create z-axis label
png(zlabel.png, bg=transparent, width=500, height=500)
grid.newpage()
pushViewport(plotViewport(c(1,1,1,1)))
print(grid.text(expression(W(u[1],u[2])==c), rot=90, gp=gpar(fontsize=85)))
popViewport()
dev.off()

## plot
persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=,
 ylab=, zlab=)
axes3d(edges=c('x--','y--','z+-'))
par3d(ignoreExtent=TRUE)


sprites3d(1.4,-0.2,0.5, radius=0.5, lit=FALSE, textype=alpha, 
texture=zlabel.png)


sprites3d(1.4,-0.2,0.5, radius=0.5, lit=FALSE, textype=alpha, 
texture=zlabel.png)
par3d(windowRect=c(0,0,600,600), zoom=1.6)
rgl.snapshot(W.png, fmt=png) # print to file



On 2011-09-09, at 17:26 , Duncan Murdoch wrote:

  On 09/09/2011 11:10 AM, Marius Hofert wrote:
  Dear all,

  Below is some code where I try to get plotmath symbols in an rgl plot. 
Duncan
  Murdoch kindly suggested to use a sprite for this. As you can see, on can 
get
  it to work, but my knowledge about grid and rgl is too limited to perfectly
  solve the problem.
  1) As you can see (please rotate the plot a little bit so that (0,0,0) is in 
front),
  the quality of the .png seems poor. Can this be improved?

  Generally plotting much larger and shrinking is the way to do this.

  2) When I print the file to a .pdf, the label is just a black square... I 
assume
  one then has to print to png again... hmmm... not perfect.

  rgl.postscript has a lot of limitations.  This one may be unavoidable.
  3) Some parts of the z-axis are covered by the label. How can this be fixed?
  The problem seems to be that sprites have a radius. Specifying a rectangle 
would
  be better here, but I am not sure how this works (or even if it does).


  I think I don't see this in Windows, which shows the black text for the 
label on a transparent background.  Do you see a solid background?  Or is it just 
the positioning?  The tick marks go out at an angle, whereas your sprite seems to 
be in line with the zy plane.  I would increase the x value to something like 1.4 
to place it properly.  (Just tried that, and I think I see the problem you had:  
the new sprite is obscured by the old one.
  You may need to play with the depth_test or depth_mask material properties.  
Not sure if those are visible on CRAN; if not, try the R-forge version.)

  The other thing I'd recommend is that before plotting the sprite, you run 
par3d(ignoreExtent=TRUE) (and then set it back to the original value afterwards).  
This will stop the label from affecting the bounding box of the plot.

  Duncan Murdoch


  Cheers,

  Marius





  require(rgl)
  require(grid)

  s- seq(0, 1, length.out=21)
  M- function(u) apply(u, 1, min)
  u- s
  v- s
  z- outer(u, v, function(u,v) M(cbind(u,v)))

  ## create z-axis label
  png(zlabel.png, bg=transparent, width=70, height=70)
  grid.newpage()
  pushViewport(plotViewport(c(1,1,1,1)))
  print(grid.text(expression(W(u[1],u[2])==c),rot=90))
  popViewport()
  dev.off()

  ## plot
  persp3d(u, v, z, aspect=iso, front=line, lit=FALSE, axes=FALSE, xlab=,
  ylab=, zlab=)
  axes3d(edges=c('x--','y--','z+-'))
  sprites3d(1,-0.4,0.65, radius=0.5, lit=FALSE, textype=alpha, 
texture=zlabel.png)
  rgl.postscript(W.pdf, fmt=pdf) # print to file

  # par3d(c(userMatrix, zoom, FOV)):
  #
  # $userMatrix
  #[,1]   [,2]   [,3] [,4]
  # [1,]  0.7262191 -0.6867201 0.031957750
  # [2,]  0.2750087  0.3328032 0.902004660
  # [3,] -0.6300603 -0.6462642 0.430542020
  # [4,]  0.000  0.000 0.1
  #
  # $zoom
  # [1] 1
  #
  # $FOV
  # [1] 30








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


Re: [R] Read a list of files into named R data.frames

2011-09-09 Thread Dennis Murphy
Hi:

Hmm, these files look familiar...Lahman database? :)

I have the 2010 set, so here's how I got this to work on my system:

files - list.files(pattern = '*.csv')
 files
 [1] Allstar.csv AllstarFull.csv
 [3] Appearances.csv AwardsManagers.csv
 [5] AwardsPlayers.csv   AwardsShareManagers.csv
 [7] AwardsSharePlayers.csv  Batting.csv
 [9] BattingPost.csv Fielding.csv
[11] FieldingOF.csv  FieldingPost.csv
[13] HallOfFame.csv  HOFold.csv
[15] lahman58-csv.zipManagers.csv
[17] ManagersHalf.csvMaster.csv
[19] Pitching.csvPitchingPost.csv
[21] Salaries.csvSchools.csv
[23] SchoolsPlayers.csv  SeriesPost.csv
[25] Teams.csv   TeamsFranchises.csv
[27] TeamsHalf.csv   Xref_Stats.csv

# Have to get rid of #15, then we can proceed:
files - files[-15]

# Get the root names out and create a similar vector
# of output files with extension .Rdata
strp - unlist(lapply(strsplit(files, '\\.'), '[', 1))
outfiles - paste(strp, 'Rdata', sep = '.')

# Function to read csv and write out an .Rdata file
f - function(x, y) {
   d - read.csv(x, header = TRUE, stringsAsFactors = FALSE)
   save(d, file = y)
  }

# Load and fire:
mapply(f, files, outfiles)

Worked for me on all but the Master file, so check your results. I
loaded in the Master file separately and saved it without problem.

HTH,
Dennis

On Fri, Sep 9, 2011 at 7:39 AM, Michael Friendly frien...@yorku.ca wrote:
 I have a collection of .csv files in a directory, and want to read them into
 R data.frames whose names
 are the same as the file names, without the .csv extension

 e.g., from
 (files - list.files(pattern=*.csv))
  [1] Allstar.csv             AllstarFull.csv
  [3] Appearances.csv         AwardsManagers.csv
  [5] AwardsPlayers.csv       AwardsShareManagers.csv
  [7] AwardsSharePlayers.csv  Batting.csv
  [9] BattingPost.csv         Fielding.csv
 [11] FieldingOF.csv          FieldingPost.csv
 [13] HallOfFame.csv          HOFold.csv
 [15] Managers.csv            ManagersHalf.csv
 [17] Master.csv              Pitching.csv
 [19] PitchingPost.csv        Salaries.csv
 [21] Schools.csv             SchoolsPlayers.csv
 [23] SeriesPost.csv          Teams.csv
 [25] TeamsFranchises.csv     TeamsHalf.csv

 Allstar - read.csv(Allstar.csv, header=TRUE)
  ...
 TeamsHalf - read.csv(TeamsHalf.csv, header=TRUE)

 Below is what I tried, which reads all the files, but doesn't create the R
 objects in the global environment.
 What is missing here?

 for (i in 1:length(files)) {
    inp - read.csv(file=files[i], header=TRUE)
    name - sub(.csv, , files[i])
    cat(Read , files[i], \trows: , nrow(inp),  cols: , ncol(inp), \n)
    eval(paste(name, - inp))
 }

 Read  Allstar.csv       rows:  4475  cols:  3
 Read  AllstarFull.csv   rows:  4676  cols:  8
 Read  Appearances.csv   rows:  94157  cols:  20
 Read  AwardsManagers.csv        rows:  57  cols:  6
 Read  AwardsPlayers.csv         rows:  2679  cols:  6
 Read  AwardsShareManagers.csv   rows:  344  cols:  7
 Read  AwardsSharePlayers.csv    rows:  6354  cols:  7
 Read  Batting.csv       rows:  93955  cols:  24
 Read  BattingPost.csv   rows:  9840  cols:  22
 Read  Fielding.csv      rows:  160710  cols:  18
 Read  FieldingOF.csv    rows:  12028  cols:  6
 Read  FieldingPost.csv  rows:  10458  cols:  17
 Read  HallOfFame.csv    rows:  3913  cols:  8
 Read  HOFold.csv        rows:  289  cols:  7
 Read  Managers.csv      rows:  3238  cols:  10
 Read  ManagersHalf.csv  rows:  93  cols:  10
 Read  Master.csv        rows:  17674  cols:  33
 Read  Pitching.csv      rows:  40432  cols:  30
 Read  PitchingPost.csv  rows:  4284  cols:  30
 Read  Salaries.csv      rows:  21464  cols:  5
 Read  Schools.csv       rows:  749  cols:  5
 Read  SchoolsPlayers.csv        rows:  6147  cols:  4
 Read  SeriesPost.csv    rows:  256  cols:  9
 Read  Teams.csv         rows:  2655  cols:  48
 Read  TeamsFranchises.csv       rows:  120  cols:  4
 Read  TeamsHalf.csv     rows:  52  cols:  10
 Read  Xref_Stats.csv    rows:  2753  cols:  3
 ls()
 [1] files i     inp   name


 --
 Michael Friendly     Email: friendly AT yorku DOT ca
 Professor, Psychology Dept.
 York University      Voice: 416 736-5115 x66249 Fax: 416 736-5814
 4700 Keele Street    Web:   http://www.datavis.ca
 Toronto, ONT  M3J 1P3 CANADA

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

2011-09-09 Thread Halabi, Anan

Hi,
Below is the code I write,
I am trying to create a metric of h and t while the values is out of function R.
First I have message warning
Second, the metric is not created


 h - seq(0.1, 0.9, by=0.1)
 t - seq(0,11000, by=100)
 z - cbind(t)
 eta=1
 beta=2
 R - array (1:1100, dim= c(110,10))
 R= exp(-(z/eta*(1-h))^(beta*(1-h)))
Warning messages:
1: In z/eta * (1 - h) :
  longer object length is not a multiple of shorter object length
2: In (z/eta * (1 - h))^(beta * (1 - h)) :
  longer object length is not a multiple of shorter object length
 R
   t
  [1,] 1.000
  [2,] 0.9995586
  [3,] 0.9974647
  [4,] 0.9919725
  [5,] 0.9801987
  [6,] 0.9572081
  [7,] 0.9141353
  [8,] 0.8341665
  [9,] 0.6833618
 [10,] 0.9892125
 [11,] 0.9825766
 [12,] 0.9727664
 [13,] 0.9583520
 [14,] 0.9370675
 [15,] 0.9051390
 [16,] 0.8559267
 [17,] 0.7769472
 [18,] 0.6423128
 [19,] 0.9629361
 [20,] 0.9521001
 [21,] 0.9382259
 [22,] 0.9201104
 [23,] 0.8958341
 [24,] 0.8622051
 [25,] 0.8136275
 [26,] 0.7395534
 [27,] 0.6175830
 [28,] 0.9246320
 [29,] 0.9127580
 [30,] 0.8982795
 [31,] 0.8800821
 [32,] 0.8564152
 [33,] 0.8244048
 [34,] 0.7790512
 [35,] 0.7109198
 [36,] 0.5996181
 [37,] 0.8767646
 [38,] 0.8671159
 [39,] 0.8550329
 [40,] 0.8394499
 [41,] 0.8187308
 [42,] 0.7902253
 [43,] 0.7493487
 [44,] 0.6874269
 [45,] 0.5854263
 [46,] 0.8215803
 [47,] 0.8170930
 [48,] 0.8098560
 [49,] 0.7988930
 [50,] 0.7827045
 [51,] 0.7588540
 [52,] 0.7231046
 [53,] 0.6673758
 [54,] 0.5736597
 [55,] 0.7611968
 [56,] 0.7642511
 [57,] 0.7637415
 [58,] 0.7588504
 [59,] 0.7482636
 [60,] 0.7297788
 [61,] 0.6994875
 [62,] 0.6498159
 [63,] 0.5635906
 [64,] 0.6975913
 [65,] 0.7098972
 [66,] 0.7174446
 [67,] 0.7196221
 [68,] 0.7153381
 [69,] 0.7026455
 [70,] 0.6779572
 [71,] 0.6341555
 [72,] 0.5547793
 [73,] 0.6325678
 [74,] 0.6551296
 [75,] 0.6715521
 [76,] 0.6814181
 [77,] 0.6838614
 [78,] 0.6771931
 [79,] 0.6581383
 [80,] 0.6199984
 [81,] 0.5469392
 [82,] 0.5677232
 [83,] 0.6008650
 [84,] 0.6265229
 [85,] 0.6443861
 [86,] 0.6537698
 [87,] 0.6532207
 [88,] 0.6397566
 [89,] 0.6070646
 [90,] 0.5398727
 [91,] 0.5044222
 [92,] 0.5478565
 [93,] 0.5827146
 [94,] 0.6086282
 [95,] 0.6250023
 [96,] 0.6305684
 [97,] 0.6226043
 [98,] 0.5951486
 [99,] 0.5334377
[100,] 0.4437834
[101,] 0.4967086
[102,] 0.5404019
[103,] 0.5742126
[104,] 0.5975006
[105,] 0.6091056
[106,] 0.6065193
[107,] 0.5840942
[108,] 0.5275280
[109,] 0.3866770
[110,] 0.4478912
[111,] 0.4997913
 length(R)
[1] 111

Anan Halabi
Reliability Eng, RD
HP Scitex
Tel: 972-9-8924648
mobil: 972-52-6624231


[[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] Reliability metric

2011-09-09 Thread R. Michael Weylandt
Your code works fine -- the warning message is just a heads up that the
object z/eta (which has length 111) does not quite match up with 1-h which
has length 10: usually this is a sign that the user defined something
improperly so R gives you a nice warning, but it doesn't halt the
calculation.

As to what you mean about your metric not being created -- if it's R (odd
choice for a variable name given your platform) it sure seems to
existespecially since you just gave us 111 values for it.

Though I think you don't know what you mean by [the] function R. Perhaps
you should read about how user defined functions are set up in R (the
fantastically useful insightful warning message giving language):
http://www.statmethods.net/management/userfunctions.html

Michael Weylandt

On Fri, Sep 9, 2011 at 2:14 PM, Halabi, Anan anan.hal...@hp.com wrote:


 Hi,
 Below is the code I write,
 I am trying to create a metric of h and t while the values is out of
 function R.
 First I have message warning
 Second, the metric is not created


  h - seq(0.1, 0.9, by=0.1)
  t - seq(0,11000, by=100)
  z - cbind(t)
  eta=1
  beta=2
  R - array (1:1100, dim= c(110,10))
  R= exp(-(z/eta*(1-h))^(beta*(1-h)))
 Warning messages:
 1: In z/eta * (1 - h) :
  longer object length is not a multiple of shorter object length
 2: In (z/eta * (1 - h))^(beta * (1 - h)) :
  longer object length is not a multiple of shorter object length
  R
   t
  [1,] 1.000
  [2,] 0.9995586
  [3,] 0.9974647
  [4,] 0.9919725
  [5,] 0.9801987
  [6,] 0.9572081
  [7,] 0.9141353
  [8,] 0.8341665
  [9,] 0.6833618
  [10,] 0.9892125
  [11,] 0.9825766
  [12,] 0.9727664
  [13,] 0.9583520
  [14,] 0.9370675
  [15,] 0.9051390
  [16,] 0.8559267
  [17,] 0.7769472
  [18,] 0.6423128
  [19,] 0.9629361
  [20,] 0.9521001
  [21,] 0.9382259
  [22,] 0.9201104
  [23,] 0.8958341
  [24,] 0.8622051
  [25,] 0.8136275
  [26,] 0.7395534
  [27,] 0.6175830
  [28,] 0.9246320
  [29,] 0.9127580
  [30,] 0.8982795
  [31,] 0.8800821
  [32,] 0.8564152
  [33,] 0.8244048
  [34,] 0.7790512
  [35,] 0.7109198
  [36,] 0.5996181
  [37,] 0.8767646
  [38,] 0.8671159
  [39,] 0.8550329
  [40,] 0.8394499
  [41,] 0.8187308
  [42,] 0.7902253
  [43,] 0.7493487
  [44,] 0.6874269
  [45,] 0.5854263
  [46,] 0.8215803
  [47,] 0.8170930
  [48,] 0.8098560
  [49,] 0.7988930
  [50,] 0.7827045
  [51,] 0.7588540
  [52,] 0.7231046
  [53,] 0.6673758
  [54,] 0.5736597
  [55,] 0.7611968
  [56,] 0.7642511
  [57,] 0.7637415
  [58,] 0.7588504
  [59,] 0.7482636
  [60,] 0.7297788
  [61,] 0.6994875
  [62,] 0.6498159
  [63,] 0.5635906
  [64,] 0.6975913
  [65,] 0.7098972
  [66,] 0.7174446
  [67,] 0.7196221
  [68,] 0.7153381
  [69,] 0.7026455
  [70,] 0.6779572
  [71,] 0.6341555
  [72,] 0.5547793
  [73,] 0.6325678
  [74,] 0.6551296
  [75,] 0.6715521
  [76,] 0.6814181
  [77,] 0.6838614
  [78,] 0.6771931
  [79,] 0.6581383
  [80,] 0.6199984
  [81,] 0.5469392
  [82,] 0.5677232
  [83,] 0.6008650
  [84,] 0.6265229
  [85,] 0.6443861
  [86,] 0.6537698
  [87,] 0.6532207
  [88,] 0.6397566
  [89,] 0.6070646
  [90,] 0.5398727
  [91,] 0.5044222
  [92,] 0.5478565
  [93,] 0.5827146
  [94,] 0.6086282
  [95,] 0.6250023
  [96,] 0.6305684
  [97,] 0.6226043
  [98,] 0.5951486
  [99,] 0.5334377
 [100,] 0.4437834
 [101,] 0.4967086
 [102,] 0.5404019
 [103,] 0.5742126
 [104,] 0.5975006
 [105,] 0.6091056
 [106,] 0.6065193
 [107,] 0.5840942
 [108,] 0.5275280
 [109,] 0.3866770
 [110,] 0.4478912
 [111,] 0.4997913
  length(R)
 [1] 111

 Anan Halabi
 Reliability Eng, RD
 HP Scitex
 Tel: 972-9-8924648
 mobil: 972-52-6624231


[[alternative HTML version deleted]]

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


[[alternative HTML version deleted]]

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


Re: [R] Read a list of files into named R data.frames

2011-09-09 Thread Michael Friendly

Thanks, David  Jean

assign(), as you suggested works to create the named objects in the 
global environment.


In my application, I also want to :
(a) save each object to a similarly-named .RData file with that object
(b) create a corresponding .Rd documentation file.

(b) works fine in my loop below, but for (a) I'm unable to find a way to 
use save() so that the

saved .RData file contains the one data.frame of the same name.

The code below
does produce all the separate .RData files, but each one contains object 
called 'inp'
rather than the name of the object created by assign().  Using 
save(name, file=...)

below produces similar .RData files with one object called 'name'.

#setwd(c:/sasuser/data/lahman)

# Read the Lahman MLB .csv files and create .Rdata and .Rd files
(files - list.files(pattern=*.csv))

for (i in 1:length(files)) {
inp - read.csv(file=files[i], header=TRUE)
cname - name - sub(.csv, , files[i])
cat(Read:, files[i], \trows: , nrow(inp),  cols: , ncol(inp), 
\n)

assign( name, inp)
save(inp, file=paste(cname, .RData, sep=))
promptData(inp, name=cname)
}

To get what I want for (a), I can save() each separately, as in

save(Allstar, file=Allstar.RData)
save(AllstarFull, file=AllstarFull.RData)
save(Appearances, file=Appearances.RData)
...

Is there some special incantation for the name of the saved object in a 
loop that I've missed?


TIA
-Michael


On 9/9/2011 10:51 AM, David Winsemius wrote:


On Sep 9, 2011, at 10:39 AM, Michael Friendly wrote:

I have a collection of .csv files in a directory, and want to read 
them into R data.frames whose names

are the same as the file names, without the .csv extension

e.g., from
 (files - list.files(pattern=*.csv))
[1] Allstar.csv AllstarFull.csv
[3] Appearances.csv AwardsManagers.csv
[5] AwardsPlayers.csv   AwardsShareManagers.csv
[7] AwardsSharePlayers.csv  Batting.csv
[9] BattingPost.csv Fielding.csv
[11] FieldingOF.csv  FieldingPost.csv
[13] HallOfFame.csv  HOFold.csv
[15] Managers.csvManagersHalf.csv
[17] Master.csv  Pitching.csv
[19] PitchingPost.csvSalaries.csv
[21] Schools.csv SchoolsPlayers.csv
[23] SeriesPost.csv  Teams.csv
[25] TeamsFranchises.csv TeamsHalf.csv

 Allstar - read.csv(Allstar.csv, header=TRUE)
 ...
 TeamsHalf - read.csv(TeamsHalf.csv, header=TRUE)

Below is what I tried, which reads all the files, but doesn't create 
the R objects in the global environment.

What is missing here?

for (i in 1:length(files)) {
   inp - read.csv(file=files[i], header=TRUE)
   name - sub(.csv, , files[i])
   cat(Read , files[i], \trows: , nrow(inp),  cols: , 
ncol(inp), \n)


Generally the assign function is used to create objects with a 
particular name. If you wanted to use eval then the text needs to be 
passed through parse() before being given to eval, but that is not the 
preferred method.


Perhaps:

 assign( files[i], inp)

Inside a for loop I think that gets done in the calling environment 
but if you were in a function  you would need to use the environment 
argument to get it to stick.




   eval(paste(name, - inp))
}

Read  Allstar.csv   rows:  4475  cols:  3
Read  AllstarFull.csv   rows:  4676  cols:  8
Read  Appearances.csv   rows:  94157  cols:  20
Read  AwardsManagers.csvrows:  57  cols:  6
Read  AwardsPlayers.csv rows:  2679  cols:  6
Read  AwardsShareManagers.csv   rows:  344  cols:  7
Read  AwardsSharePlayers.csvrows:  6354  cols:  7
Read  Batting.csv   rows:  93955  cols:  24
Read  BattingPost.csv   rows:  9840  cols:  22
Read  Fielding.csv  rows:  160710  cols:  18
Read  FieldingOF.csvrows:  12028  cols:  6
Read  FieldingPost.csv  rows:  10458  cols:  17
Read  HallOfFame.csvrows:  3913  cols:  8
Read  HOFold.csvrows:  289  cols:  7
Read  Managers.csv  rows:  3238  cols:  10
Read  ManagersHalf.csv  rows:  93  cols:  10
Read  Master.csvrows:  17674  cols:  33
Read  Pitching.csv  rows:  40432  cols:  30
Read  PitchingPost.csv  rows:  4284  cols:  30
Read  Salaries.csv  rows:  21464  cols:  5
Read  Schools.csv   rows:  749  cols:  5
Read  SchoolsPlayers.csvrows:  6147  cols:  4
Read  SeriesPost.csvrows:  256  cols:  9
Read  Teams.csv rows:  2655  cols:  48
Read  TeamsFranchises.csv   rows:  120  cols:  4
Read  TeamsHalf.csv rows:  52  cols:  10
Read  Xref_Stats.csvrows:  2753  cols:  3
 ls()
[1] files i inp   name


--
Michael Friendly Email: friendly AT yorku DOT ca
Professor, Psychology Dept.
York University  Voice: 416 736-5115 x66249 Fax: 416 736-5814
4700 Keele StreetWeb:   http://www.datavis.ca
Toronto, ONT  M3J 1P3 CANADA

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

and provide 

Re: [R] Reliability metric

2011-09-09 Thread R. Michael Weylandt
Well, I would point out to you that you only have 9 values of h which
certainly is a problem in calculating 10 colmans(?) -- more generally, I'd
suggest you learn more about the difference between vectorized operations
and operations on vectors: things like `%*%` and outer() will help you
greatly. If you are perhaps coming from a Matlab background, keep in mind
that R's default behavior is to include the . before all operations. E.g.,


R x = 1:3
R x^2
[1] 1 4 9

R x * x
[1] 1 4 9

R x %*% x
  [,1]
[1,] 14

R outer(x,x)
 [,1] [,2] [,3]
[1,]123
[2,]246
[3,]369

Michael


On Fri, Sep 9, 2011 at 2:32 PM, Halabi, Anan anan.hal...@hp.com wrote:

 What I mean by the metric not created,

 Desired metric is of 1100 values by 110 rows (t), and 10 Colman (h)

 Thank, 

 ** **

 *From:* R. Michael Weylandt [mailto:michael.weyla...@gmail.com]
 *Sent:* Friday, September 09, 2011 10:25 PM
 *To:* Halabi, Anan
 *Cc:* r-help@R-project.org
 *Subject:* Re: [R] Reliability metric

 ** **

 Your code works fine -- the warning message is just a heads up that the
 object z/eta (which has length 111) does not quite match up with 1-h which
 has length 10: usually this is a sign that the user defined something
 improperly so R gives you a nice warning, but it doesn't halt the
 calculation.

 As to what you mean about your metric not being created -- if it's R (odd
 choice for a variable name given your platform) it sure seems to
 existespecially since you just gave us 111 values for it.

 Though I think you don't know what you mean by [the] function R. Perhaps
 you should read about how user defined functions are set up in R (the
 fantastically useful insightful warning message giving language):
 http://www.statmethods.net/management/userfunctions.html

 Michael Weylandt

 On Fri, Sep 9, 2011 at 2:14 PM, Halabi, Anan anan.hal...@hp.com wrote:**
 **


 Hi,
 Below is the code I write,
 I am trying to create a metric of h and t while the values is out of
 function R.
 First I have message warning
 Second, the metric is not created


  h - seq(0.1, 0.9, by=0.1)
  t - seq(0,11000, by=100)
  z - cbind(t)
  eta=1
  beta=2
  R - array (1:1100, dim= c(110,10))
  R= exp(-(z/eta*(1-h))^(beta*(1-h)))
 Warning messages:
 1: In z/eta * (1 - h) :
  longer object length is not a multiple of shorter object length
 2: In (z/eta * (1 - h))^(beta * (1 - h)) :
  longer object length is not a multiple of shorter object length
  R
   t
  [1,] 1.000
  [2,] 0.9995586
  [3,] 0.9974647
  [4,] 0.9919725
  [5,] 0.9801987
  [6,] 0.9572081
  [7,] 0.9141353
  [8,] 0.8341665
  [9,] 0.6833618
  [10,] 0.9892125
  [11,] 0.9825766
  [12,] 0.9727664
  [13,] 0.9583520
  [14,] 0.9370675
  [15,] 0.9051390
  [16,] 0.8559267
  [17,] 0.7769472
  [18,] 0.6423128
  [19,] 0.9629361
  [20,] 0.9521001
  [21,] 0.9382259
  [22,] 0.9201104
  [23,] 0.8958341
  [24,] 0.8622051
  [25,] 0.8136275
  [26,] 0.7395534
  [27,] 0.6175830
  [28,] 0.9246320
  [29,] 0.9127580
  [30,] 0.8982795
  [31,] 0.8800821
  [32,] 0.8564152
  [33,] 0.8244048
  [34,] 0.7790512
  [35,] 0.7109198
  [36,] 0.5996181
  [37,] 0.8767646
  [38,] 0.8671159
  [39,] 0.8550329
  [40,] 0.8394499
  [41,] 0.8187308
  [42,] 0.7902253
  [43,] 0.7493487
  [44,] 0.6874269
  [45,] 0.5854263
  [46,] 0.8215803
  [47,] 0.8170930
  [48,] 0.8098560
  [49,] 0.7988930
  [50,] 0.7827045
  [51,] 0.7588540
  [52,] 0.7231046
  [53,] 0.6673758
  [54,] 0.5736597
  [55,] 0.7611968
  [56,] 0.7642511
  [57,] 0.7637415
  [58,] 0.7588504
  [59,] 0.7482636
  [60,] 0.7297788
  [61,] 0.6994875
  [62,] 0.6498159
  [63,] 0.5635906
  [64,] 0.6975913
  [65,] 0.7098972
  [66,] 0.7174446
  [67,] 0.7196221
  [68,] 0.7153381
  [69,] 0.7026455
  [70,] 0.6779572
  [71,] 0.6341555
  [72,] 0.5547793
  [73,] 0.6325678
  [74,] 0.6551296
  [75,] 0.6715521
  [76,] 0.6814181
  [77,] 0.6838614
  [78,] 0.6771931
  [79,] 0.6581383
  [80,] 0.6199984
  [81,] 0.5469392
  [82,] 0.5677232
  [83,] 0.6008650
  [84,] 0.6265229
  [85,] 0.6443861
  [86,] 0.6537698
  [87,] 0.6532207
  [88,] 0.6397566
  [89,] 0.6070646
  [90,] 0.5398727
  [91,] 0.5044222
  [92,] 0.5478565
  [93,] 0.5827146
  [94,] 0.6086282
  [95,] 0.6250023
  [96,] 0.6305684
  [97,] 0.6226043
  [98,] 0.5951486
  [99,] 0.5334377
 [100,] 0.4437834
 [101,] 0.4967086
 [102,] 0.5404019
 [103,] 0.5742126
 [104,] 0.5975006
 [105,] 0.6091056
 [106,] 0.6065193
 [107,] 0.5840942
 [108,] 0.5275280
 [109,] 0.3866770
 [110,] 0.4478912
 [111,] 0.4997913
  length(R)
 [1] 111

 Anan Halabi
 Reliability Eng, RD
 HP Scitex
 Tel: 972-9-8924648
 mobil: 972-52-6624231


[[alternative HTML version deleted]]

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

 ** **


[[alternative HTML version deleted]]


Re: [R] Reliability metric

2011-09-09 Thread David Winsemius


On Sep 9, 2011, at 3:14 PM, Halabi, Anan wrote:



Hi,
Below is the code I write,
I am trying to create a metric of h and t while the values is out of  
function R.

First I have message warning


It is only a warning.


Second, the metric is not created


What not created means has already been question by Weylandt. What  
you seem to have missed in your introductory material is the concept  
of recycling of arguments. Look at this code and its results:


 (1:10) * (1:2)
# [1]  1  4  3  8  5 12  7 16  9 20

 (1:10) * (1:3)
# [1]  1  4  9  4 10 18  7 16 27 10
Warning message:
In (1:10) * (1:3) :
  longer object length is not a multiple of shorter object length

--
David.




h - seq(0.1, 0.9, by=0.1)
t - seq(0,11000, by=100)
z - cbind(t)
eta=1
beta=2
R - array (1:1100, dim= c(110,10))
R= exp(-(z/eta*(1-h))^(beta*(1-h)))

Warning messages:
1: In z/eta * (1 - h) :
 longer object length is not a multiple of shorter object length
2: In (z/eta * (1 - h))^(beta * (1 - h)) :
 longer object length is not a multiple of shorter object length

R

  t
 [1,] 1.000
 [2,] 0.9995586
 [3,] 0.9974647
 [4,] 0.9919725
 [5,] 0.9801987
 [6,] 0.9572081
 [7,] 0.9141353
 [8,] 0.8341665
 [9,] 0.6833618
[10,] 0.9892125
[11,] 0.9825766
[12,] 0.9727664
[13,] 0.9583520
[14,] 0.9370675
[15,] 0.9051390
[16,] 0.8559267
[17,] 0.7769472
[18,] 0.6423128
[19,] 0.9629361
[20,] 0.9521001
[21,] 0.9382259
[22,] 0.9201104
[23,] 0.8958341
[24,] 0.8622051
[25,] 0.8136275
[26,] 0.7395534
[27,] 0.6175830
[28,] 0.9246320
[29,] 0.9127580
[30,] 0.8982795
[31,] 0.8800821
[32,] 0.8564152
[33,] 0.8244048
[34,] 0.7790512
[35,] 0.7109198
[36,] 0.5996181
[37,] 0.8767646
[38,] 0.8671159
[39,] 0.8550329
[40,] 0.8394499
[41,] 0.8187308
[42,] 0.7902253
[43,] 0.7493487
[44,] 0.6874269
[45,] 0.5854263
[46,] 0.8215803
[47,] 0.8170930
[48,] 0.8098560
[49,] 0.7988930
[50,] 0.7827045
[51,] 0.7588540
[52,] 0.7231046
[53,] 0.6673758
[54,] 0.5736597
[55,] 0.7611968
[56,] 0.7642511
[57,] 0.7637415
[58,] 0.7588504
[59,] 0.7482636
[60,] 0.7297788
[61,] 0.6994875
[62,] 0.6498159
[63,] 0.5635906
[64,] 0.6975913
[65,] 0.7098972
[66,] 0.7174446
[67,] 0.7196221
[68,] 0.7153381
[69,] 0.7026455
[70,] 0.6779572
[71,] 0.6341555
[72,] 0.5547793
[73,] 0.6325678
[74,] 0.6551296
[75,] 0.6715521
[76,] 0.6814181
[77,] 0.6838614
[78,] 0.6771931
[79,] 0.6581383
[80,] 0.6199984
[81,] 0.5469392
[82,] 0.5677232
[83,] 0.6008650
[84,] 0.6265229
[85,] 0.6443861
[86,] 0.6537698
[87,] 0.6532207
[88,] 0.6397566
[89,] 0.6070646
[90,] 0.5398727
[91,] 0.5044222
[92,] 0.5478565
[93,] 0.5827146
[94,] 0.6086282
[95,] 0.6250023
[96,] 0.6305684
[97,] 0.6226043
[98,] 0.5951486
[99,] 0.5334377
[100,] 0.4437834
[101,] 0.4967086
[102,] 0.5404019
[103,] 0.5742126
[104,] 0.5975006
[105,] 0.6091056
[106,] 0.6065193
[107,] 0.5840942
[108,] 0.5275280
[109,] 0.3866770
[110,] 0.4478912
[111,] 0.4997913

length(R)

[1] 111

Anan Halabi
Reliability Eng, RD
HP Scitex
Tel: 972-9-8924648
mobil: 972-52-6624231


[[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] Reliability metric

2011-09-09 Thread William Dunlap
I may be wrong, but I thought the original poster
was trying to make a matrix ('metric' sp.) of values,
R, so that
   R[i,j] == f(t[i], h[i])
for all i and j where
   f - function(t, h)exp(-(t/eta*(1-h))^(beta*(1-h))
outer() would do that:
   R - outer(t, h, function(t, h)exp(-(t/eta*(1-h))^(beta*(1-h
It dimensions are length(t) by length(h), 111 by 9, not the
110 by 10 mentioned in the post.

Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com 

 -Original Message-
 From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
 Behalf Of David Winsemius
 Sent: Friday, September 09, 2011 1:25 PM
 To: Halabi, Anan
 Cc: r-help@R-project.org
 Subject: Re: [R] Reliability metric
 Importance: High
 
 
 On Sep 9, 2011, at 3:14 PM, Halabi, Anan wrote:
 
 
  Hi,
  Below is the code I write,
  I am trying to create a metric of h and t while the values is out of
  function R.
  First I have message warning
 
 It is only a warning.
 
  Second, the metric is not created
 
 What not created means has already been question by Weylandt. What
 you seem to have missed in your introductory material is the concept
 of recycling of arguments. Look at this code and its results:
 
   (1:10) * (1:2)
 # [1]  1  4  3  8  5 12  7 16  9 20
 
   (1:10) * (1:3)
 # [1]  1  4  9  4 10 18  7 16 27 10
 Warning message:
 In (1:10) * (1:3) :
longer object length is not a multiple of shorter object length
 
 --
 David.
 
 
  h - seq(0.1, 0.9, by=0.1)
  t - seq(0,11000, by=100)
  z - cbind(t)
  eta=1
  beta=2
  R - array (1:1100, dim= c(110,10))
  R= exp(-(z/eta*(1-h))^(beta*(1-h)))
  Warning messages:
  1: In z/eta * (1 - h) :
   longer object length is not a multiple of shorter object length
  2: In (z/eta * (1 - h))^(beta * (1 - h)) :
   longer object length is not a multiple of shorter object length
  R
t
   [1,] 1.000
   [2,] 0.9995586
   [3,] 0.9974647
   [4,] 0.9919725
   [5,] 0.9801987
   [6,] 0.9572081
   [7,] 0.9141353
   [8,] 0.8341665
   [9,] 0.6833618
  [10,] 0.9892125
  [11,] 0.9825766
  [12,] 0.9727664
  [13,] 0.9583520
  [14,] 0.9370675
  [15,] 0.9051390
  [16,] 0.8559267
  [17,] 0.7769472
  [18,] 0.6423128
  [19,] 0.9629361
  [20,] 0.9521001
  [21,] 0.9382259
  [22,] 0.9201104
  [23,] 0.8958341
  [24,] 0.8622051
  [25,] 0.8136275
  [26,] 0.7395534
  [27,] 0.6175830
  [28,] 0.9246320
  [29,] 0.9127580
  [30,] 0.8982795
  [31,] 0.8800821
  [32,] 0.8564152
  [33,] 0.8244048
  [34,] 0.7790512
  [35,] 0.7109198
  [36,] 0.5996181
  [37,] 0.8767646
  [38,] 0.8671159
  [39,] 0.8550329
  [40,] 0.8394499
  [41,] 0.8187308
  [42,] 0.7902253
  [43,] 0.7493487
  [44,] 0.6874269
  [45,] 0.5854263
  [46,] 0.8215803
  [47,] 0.8170930
  [48,] 0.8098560
  [49,] 0.7988930
  [50,] 0.7827045
  [51,] 0.7588540
  [52,] 0.7231046
  [53,] 0.6673758
  [54,] 0.5736597
  [55,] 0.7611968
  [56,] 0.7642511
  [57,] 0.7637415
  [58,] 0.7588504
  [59,] 0.7482636
  [60,] 0.7297788
  [61,] 0.6994875
  [62,] 0.6498159
  [63,] 0.5635906
  [64,] 0.6975913
  [65,] 0.7098972
  [66,] 0.7174446
  [67,] 0.7196221
  [68,] 0.7153381
  [69,] 0.7026455
  [70,] 0.6779572
  [71,] 0.6341555
  [72,] 0.5547793
  [73,] 0.6325678
  [74,] 0.6551296
  [75,] 0.6715521
  [76,] 0.6814181
  [77,] 0.6838614
  [78,] 0.6771931
  [79,] 0.6581383
  [80,] 0.6199984
  [81,] 0.5469392
  [82,] 0.5677232
  [83,] 0.6008650
  [84,] 0.6265229
  [85,] 0.6443861
  [86,] 0.6537698
  [87,] 0.6532207
  [88,] 0.6397566
  [89,] 0.6070646
  [90,] 0.5398727
  [91,] 0.5044222
  [92,] 0.5478565
  [93,] 0.5827146
  [94,] 0.6086282
  [95,] 0.6250023
  [96,] 0.6305684
  [97,] 0.6226043
  [98,] 0.5951486
  [99,] 0.5334377
  [100,] 0.4437834
  [101,] 0.4967086
  [102,] 0.5404019
  [103,] 0.5742126
  [104,] 0.5975006
  [105,] 0.6091056
  [106,] 0.6065193
  [107,] 0.5840942
  [108,] 0.5275280
  [109,] 0.3866770
  [110,] 0.4478912
  [111,] 0.4997913
  length(R)
  [1] 111
 
  Anan Halabi
  Reliability Eng, RD
  HP Scitex
  Tel: 972-9-8924648
  mobil: 972-52-6624231
 
 
  [[alternative HTML version deleted]]
 
  __
  R-help@r-project.org mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
 
 David Winsemius, MD
 West Hartford, CT
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] function censReg in panel data setting

2011-09-09 Thread Igors
Dear Arne,


Thank you for your answer. However it doesn't solve my problem fully. 

The problem is that I have much bigger data set than I sent to you (it was
only a small part : 874 obs.). My full data set is 546718 obs.

If I try to use censReg on full data set, then it still gives me the same
already mentioned error about Na's in the initial gradient.


I have sent you an e-mail with full dataset and the code. I would really
appreciate if you could check how it works and suggest me how to solve this
problem. 

I have noticed that you use iterlim argumet to specify maximum number of
iterations.  How this argument could affect possibility of obtaining
estimates? 



/Igors

--
View this message in context: 
http://r.789695.n4.nabble.com/function-censReg-in-panel-data-setting-tp3792227p3802648.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] Looping through multiple datasets already in memory

2011-09-09 Thread MatAra
Gents,

I have the following loop:

for (i in (seq(along=Draw.1[,1]))) {print(i)} # from 1 - 5 (counter)
for (i in (seq(along=Draw.1[,1]))) {print(Draw[i,4])} # from the 1 - column
(passes parameter per column) 

Now: I have multiple Draws [Draw.1 - Draw.100] in memory with different
orders
Question: How do I loop through  each Draw.[j]? OR in other way, how do I
pick up all the datasets from memory dynamically?

Something like: 
for (i in (seq(along=Draw.[J][,1]))) {print(Draw[i,4])}

paste() would not allow me to pick the [row,col] combination I need. It just
passes the list of all the elements.
paste(as.matrix(Draw), 1, sep=.)[10,5] DOES NOT WORK 

Thanks in advanced for the ideas!

An example of how the draw (control files) looks like:

Draw.1
  ID Name Tab   Folder  L/S 
 [1,]   38 Stoxx50_24_08_11   38  Stoxx -1
 [2,]   47 Stoxx50_24_08_11   47  Stoxx  1
 [3,]  153 DAX_31_08_11   29  DAX1
 [4,]  256 FT100_31_08_11 12  UK 1
 [5,]  303 FT100_31_08_11 59  UK 1


--
View this message in context: 
http://r.789695.n4.nabble.com/Looping-through-multiple-datasets-already-in-memory-tp3802453p3802453.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] Question about plot.mona {cluster}

2011-09-09 Thread trekvana
Hello all, I what to print the banner plot that is output from the mona
method in the cluster package but the problem is I dont want to print all
that red ink. Here is an example:

data(animals)
ma - mona(animals)
ma
## Plot similar to Figure 10 in Struyf et al (1996)
plot(ma)

I can change the bar color by using the argument col=c(0,0) -
plot(ma,col=c(0,0)) - but then the variable labels also get colored white.
Is there a way to remove the bar colors but leave the variables colored?

Thank you
George Skountrianos

--
View this message in context: 
http://r.789695.n4.nabble.com/Question-about-plot-mona-cluster-tp3802672p3802672.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] GLM Question

2011-09-09 Thread Jim Silverton
The minimum achievable level of significance is defined asthe minimum of
Prob(Y=y) over all y's. If I have GLM with a treatment and replicate and I
would like to find out how to compute the minimum achievable level of
significance for that GLM in R

For example, how do I do this for the following data:
Treat_Rep1Treat_Rep2   Control_Rep1 Control_Rep2
  4  10 2  8




-- 
Thanks,
Jim.

[[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] Looping through multiple datasets already in memory

2011-09-09 Thread Greg Snow
The best approach is to put all those datasets into a list, then you can loop 
through the elements of the list, or use lapply or sapply on the list.

If you insist on keeping them outside of a list, or need a loop to put them 
into the list, then look at the 'get' function.

-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
Behalf Of MatAra
Sent: Friday, September 09, 2011 12:26 PM
To: r-help@r-project.org
Subject: [R] Looping through multiple datasets already in memory

Gents,

I have the following loop:

for (i in (seq(along=Draw.1[,1]))) {print(i)} # from 1 - 5 (counter)
for (i in (seq(along=Draw.1[,1]))) {print(Draw[i,4])} # from the 1 - column
(passes parameter per column) 

Now: I have multiple Draws [Draw.1 - Draw.100] in memory with different
orders
Question: How do I loop through  each Draw.[j]? OR in other way, how do I
pick up all the datasets from memory dynamically?

Something like: 
for (i in (seq(along=Draw.[J][,1]))) {print(Draw[i,4])}

paste() would not allow me to pick the [row,col] combination I need. It just
passes the list of all the elements.
paste(as.matrix(Draw), 1, sep=.)[10,5] DOES NOT WORK 

Thanks in advanced for the ideas!

An example of how the draw (control files) looks like:

Draw.1
  ID Name Tab   Folder  L/S 
 [1,]   38 Stoxx50_24_08_11   38  Stoxx -1
 [2,]   47 Stoxx50_24_08_11   47  Stoxx  1
 [3,]  153 DAX_31_08_11   29  DAX1
 [4,]  256 FT100_31_08_11 12  UK 1
 [5,]  303 FT100_31_08_11 59  UK 1


--
View this message in context: 
http://r.789695.n4.nabble.com/Looping-through-multiple-datasets-already-in-memory-tp3802453p3802453.html
Sent from the R help mailing list archive at Nabble.com.

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

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


[R] Different results with arima in R 2.12.2 and R 2.11.1

2011-09-09 Thread Luis Felipe Parra
Hello , I have estimated the following model, a sarima:

p=9
d=1
q=2
P=0
D=1
Q=1
S=12


In R 2.12.2
Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q),
period = S),
optim.control = list(reltol = tol))

Coefficients:
 ar1 ar2  ar3 ar4 ar5 ar6  ar7  ar8
ar9
  0.3152  0.8762  -0.4413  0.0152  0.1500  0.0001  -0.0413  -0.1811
 0.0646
s.e.  0.0865  0.0885   0.1141  0.1181  0.1196  0.1220   0.1120   0.0908
 0.0865
  ma1  ma2 sma1
  -0.0221  -0.9779  -0.7635
s.e.   0.0539   0.0534   0.0834

sigma^2 estimated as 1.965e+17:  log likelihood = -3316.07,  aic = 6658.13


and in In R 2.11.1
Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q),
period = S),
optim.control = list(reltol = tol))

Coefficients:
 ar1 ar2  ar3 ar4 ar5ar6  ar7  ar8
ar9
  0.3152  0.8761  -0.4413  0.0153  0.1500  0.000  -0.0413  -0.1810
 0.0646
s.e.  0.0865  0.0885   0.1141  0.1181  0.1196  0.122   0.1120   0.0908
 0.0865
  ma1  ma2 sma1
  -0.0221  -0.9779  -0.7635
s.e.   0.0539   0.0534   0.0834

sigma^2 estimated as 1.965e+17:  log likelihood = -3316.07,  aic = 6658.13

and as you can see in the results some coefficients (for example ar2 and
ar8) are different in the different R versions. does anybody know what might
be going on. Was there any change in the arima function between the two
versions?

Thank you

Felipe Parra

-- 

Este mensaje de correo electrónico es enviado por Quantil S.A.S y puede
contener información confidencial o privilegiada.

This e-mail is sent by Quantil S.A.S and may contain confidential or
privileged information

[[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] Using SVG in R

2011-09-09 Thread Fong Chun Chan
Hi,

I've been trying to use SVG in R using the function svg, but all my my
graphics end up with the text being very small to be point you can't even
see them.  I've attached a screenshot of what the output for a histogram
looks like.  Here is my code:

CairoSVG('file/location')
hist()
dev.off()

Am I missing something here?  Thanks,

Fong
attachment: Screen shot 2011-09-09 at 3.16.50 PM.png__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.


Re: [R] Using SVG in R

2011-09-09 Thread R. Michael Weylandt
Perhaps the problem could be solved by playing with the pointsize argument
to CairoSVG or the various plotting parameters of hist().

More generally, this code:

R library(CairoSVG)
R CairoSVG()
R hist(rpois(500,3))
R dev.off()

produced an entirely acceptable SVG image on my computer so I can't
replicate the problem...

As a general life question, why in a question about SVG files would you send
us a PNG screenshot rather than the SVG itself? Can SVG's not be sent over
the mailing list?

Michael


On Fri, Sep 9, 2011 at 5:19 PM, Fong Chun Chan fongchunc...@gmail.comwrote:

 Hi,

 I've been trying to use SVG in R using the function svg, but all my my
 graphics end up with the text being very small to be point you can't even
 see them.  I've attached a screenshot of what the output for a histogram
 looks like.  Here is my code:

 CairoSVG('file/location')
 hist()
 dev.off()

 Am I missing something here?  Thanks,

 Fong

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


  1   2   >