[R] why does this simple example NOT work?

2012-07-20 Thread Martin Ivanov
 Hello,

I want to create and save objects in a loop, but this is precluded by the 
following obstacle:
this part of the script fails to work:

assign(x=paste(a, 1, sep=), value=1);
save(paste(a, 1, sep=), file=paste(paste(a, 1, sep=), .RData, sep=))

Do you know any workaround? I am already out of ideas.

Best regards,

Martin

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] why does this simple example NOT work?

2012-07-20 Thread Berend Hasselman

Martin Ivanov wrote
 
 Hello,
 
 I want to create and save objects in a loop, but this is precluded by the
 following obstacle:
 this part of the script fails to work:
 
 assign(x=paste(a, 1, sep=), value=1);
 save(paste(a, 1, sep=), file=paste(paste(a, 1, sep=), .RData,
 sep=))
 
 Do you know any workaround? I am already out of ideas.
 

You haven't given the error message R provides.
Carefully read the documentation of save.

save(list=paste(a, 1, sep=), file=paste(paste(a, 1, sep=), .RData,
sep=))

This works for for me.

Berend




--
View this message in context: 
http://r.789695.n4.nabble.com/why-does-this-simple-example-NOT-work-tp4637158p4637159.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] convert deldir$delsgs to a X3D IndexedTriangleSet

2012-07-20 Thread Erdal Karaca
I have managed to use the tripack package which is a bit simpler to use
(for me) and performs better (timely)...

##
# inputFile format: x y z
inputFile - d:\\inputfile.xyz.txt

# output file will contain a x3d scene embedded in a html file to be viewed
in a browser supporting WebGL/X3DOM
outputFile - d:\\scene.x3d.html

library(tripack)
points - read.table( inputFile )
mesh - tri.mesh(x, y)
tris - triangles(mesh)
nodes - cbind(tris[,1], tris[,2], tris[,3])

write(htmlheadlink rel='stylesheet' type='text/css' href='
http://www.x3dom.org/download/x3dom.css'/link/headbody
X3D profile='Interchange'
style='width:100%;height100%'SceneTransformShape DEF='river'
AppearanceMaterial ambientIntensity='0.0' diffuseColor='0.7 0.4 0'
shininess='1.0'//Appearance
IndexedTriangleSet index=', file=outputFile, append=FALSE)
write.table(nodes, file=outputFile, row.names=FALSE, col.names=FALSE,
append = TRUE)
write('Coordinate point=', file=outputFile, append=TRUE)
write.table(cbind(x,y,z), file=outputFile, row.names=FALSE,
col.names=FALSE, append = TRUE)
write('//IndexedTriangleSet/Shape/Transform/Scene/X3D
script type='text/javascript'
src='http://www.x3dom.org/download/x3dom.js'/script/body/html,
file=outputFile, append=TRUE)

##

2012/7/20 Rolf Turner rolf.tur...@xtra.co.nz

 On 19/07/12 20:40, Erdal Karaca wrote:

 Oh, I dont want anyone to do anything for me. I am just asking for
 hints/infos to be able to do it myself...
 I was hoping that someone has already done this...
 But anyways, I will post my solution once I have succeeded...


 If you haven't already looked at it, you may find it helpful to experiment
 with the triang.list() function and the z argument to deldir().

 If I understand correctly what you are after, the output of triang.list()
 may go some way toward providing it.

 cheers,

 Rolf Turner


[[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] median comparison tests

2012-07-20 Thread Tal Galili
You can do a transformation on your Y, and than use the normal tools.
Using something like a median test is not advisable due to low power, and a
wilcox.test is considered a better choice (although the meaning is more
complex when observations are not symmetrically distributed).

If you wish to perform a one way anova, which is non parametric, you can
use:
?kruskal.test
And after that you can use the tools suggested here by David
Winsemiushttp://r.789695.n4.nabble.com/template/NamlServlet.jtp?macro=user_nodesuser=4999
:
http://r.789695.n4.nabble.com/post-hoc-test-with-kruskal-test-td898661.html

With regards,
Tal



Contact
Details:---
Contact me: tal.gal...@gmail.com |  972-52-7275845
Read me: www.talgalili.com (Hebrew) | www.biostatistics.co.il (Hebrew) |
www.r-statistics.com (English)
--




On Thu, Jul 19, 2012 at 11:14 PM, Data Analytics Corp. 
w...@dataanalyticscorp.com wrote:

 Hi,

 A client has a consumption measure on each of four products.  The sample
 size is 75.  The consumption distributions are highly skewed for each
 product.  He would like a pairwise comparison test of the products, much
 like Tukey's HSD but using medians rather than means.  Is there such a
 median comparison test in R?

 Thanks,

 Walt

 

 Walter R. Paczkowski, Ph.D.
 Data Analytics Corp.
 44 Hamilton Lane
 Plainsboro, NJ 08536
 
 (V) 609-936-8999
 (F) 609-936-3733
 w...@dataanalyticscorp.com
 www.dataanalyticscorp.com

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


[[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] UTK mirror

2012-07-20 Thread Jason Curole
Hi all,

The past couple of days I have been trying to update my R installation
(Debian) and I am getting the following:


Get:9 http://http.us.debian.org squeeze/main amd64 Packages [6,542 kB]
Ign http://mirrors.nics.utk.edu squeeze-cran/ Release
Ign http://mirrors.nics.utk.edu squeeze-cran/ Sources/DiffIndex
Ign http://mirrors.nics.utk.edu squeeze-cran/ Packages/DiffIndex
Ign http://mirrors.nics.utk.edu squeeze-cran/ Sources
Ign http://mirrors.nics.utk.edu squeeze-cran/ Packages
Err http://mirrors.nics.utk.edu squeeze-cran/ Sources
  Undetermined Error
Err http://mirrors.nics.utk.edu squeeze-cran/ Packages
  Undetermined Error
Fetched 6,766 kB in 11min 26s (9,854 B/s)
W: Failed to fetch
http://mirrors.nics.utk.edu/cran/bin/linux/debian/squeeze-cran/Sources.gz
 Undetermined Error

W: Failed to fetch
http://mirrors.nics.utk.edu/cran/bin/linux/debian/squeeze-cran/Packages.gz
 Undetermined Error

E: Some index files failed to download, they have been ignored, or old
ones used instead.


I notice that http://mirros.nics.utk.edu now takes me to
https://magic.nics.utk.edu/home and that the link to the CRAN UTK
archive gives a 404 error (page not found).  Is this mirror dead?

I searched the archive and can't find any threads on the UTK mirror
for June or July.

Thanks, Jason

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Removing values from a string

2012-07-20 Thread arun
Hi,
This should make much more sense

onenew-data.frame(keyword=gsub((NA){0,1}\\|, ,one$keyword))
onenew1-data.frame(keyword=gsub((NA){0,1},,onenew$keyword))
onenew1
   keyword
1 auto
2 auto insurance quote
3   auto insurance
4    insurance
5   auto insurance
6 

A.K.



- Original Message -
From: Abraham Mathew abmathe...@gmail.com
To: r-help@r-project.org
Cc: 
Sent: Thursday, July 19, 2012 3:21 PM
Subject: [R] Removing values from a string

So I have the following data frame and I want to know how I can remove all
NA values from each string, and also
remove all | values from the START of the string. So they should
something like auto|insurance or auto|insurance|quote

one = data.frame(keyword=c(|auto, NA|auto|insurance|quote,
NA|auto|insurance,
                           NA|insurance, NA|auto|insurance, NA))

one


Can anyone point me in the right direction? I'm still not too familiar with
regex or gsub to find a solution, and there doesn't seem
to be anything helpful in the stringr package for this task.


Thanks

-- 
*Abraham Mathew
Statistical Analyst
www.amathew.com
720-648-0108
@abmathewks*

    [[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] A graph with a positive y-axis intercept

2012-07-20 Thread arun
HI,

I guess you wanted to change the scale of Y-axis:

 plot(X,Y,ylim=c(0,25))

A.K.



- Original Message -
From: Lekgatlhamang, lexi Setlhare lexisetlh...@yahoo.com
To: r-h...@stat.math.ethz.ch r-h...@stat.math.ethz.ch
Cc: 
Sent: Thursday, July 19, 2012 5:31 PM
Subject: [R] A graph with a positive y-axis intercept

Dear all,
 
I have a challenge with a supposedly simple graph (a scatter). I wanted to use 
R to create a plot/graph with a positive y-axis intercept but it does not seem 
to give me what I expect to see. As an example, I used the following values
 
 X- c(0, 4, 8, 11)
 Y- c(8, 12, 15, 22)

then I used the plot function in the base package
 
 plot(X, Y)
 
This gave me a graph which might be mistaken to be starting from the origin, 
and yet, using the same numbers in Excel, and plotting them using the scatter 
option, I get a graph which clearly has a positive y-axis intercept.
Please help.
 
Lexi
    [[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] generic functions question in building a new package

2012-07-20 Thread zgu9
Hi all, 

I'm trying to build a package 'Thiele', but run into problems with generic
functions. 
I have a class Benefit, and defined the function print.Benefit. 
When I try R CMD check Thiele in X11, this warning came out
---
Functions/methods with usage in documentation object 'Benefit' but not in
code:
  print
---
Do I need to write a new print Rd document into R file in that Package? If
so, how can I write it?
I read the section 7 Generic functions and methods in the 'Writing R
Extensions' manual. But I'm still not clear how to figure out my problem.

Thanks!
Eric

--
View this message in context: 
http://r.789695.n4.nabble.com/generic-functions-question-in-building-a-new-package-tp4637129.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] convert deldir$delsgs to a X3D IndexedTriangleSet

2012-07-20 Thread Rolf Turner

On 19/07/12 20:40, Erdal Karaca wrote:
Oh, I dont want anyone to do anything for me. I am just asking for 
hints/infos to be able to do it myself...

I was hoping that someone has already done this...
But anyways, I will post my solution once I have succeeded...


If you haven't already looked at it, you may find it helpful to experiment
with the triang.list() function and the z argument to deldir().

If I understand correctly what you are after, the output of triang.list()
may go some way toward providing it.

cheers,

Rolf Turner

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

2012-07-20 Thread arun
Hi,

Check this link.

http://r.789695.n4.nabble.com/help-comparing-two-median-with-R-td822998.html
A.K.



- Original Message -
From: Data Analytics Corp. w...@dataanalyticscorp.com
To: r-help@r-project.org
Cc: 
Sent: Thursday, July 19, 2012 4:14 PM
Subject: [R] median comparison tests

Hi,

A client has a consumption measure on each of four products.  The sample size 
is 75.  The consumption distributions are highly skewed for each product.  He 
would like a pairwise comparison test of the products, much like Tukey's HSD 
but using medians rather than means.  Is there such a median comparison test in 
R?

Thanks,

Walt



Walter R. Paczkowski, Ph.D.
Data Analytics Corp.
44 Hamilton Lane
Plainsboro, NJ 08536

(V) 609-936-8999
(F) 609-936-3733
w...@dataanalyticscorp.com
www.dataanalyticscorp.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] change file name from file0.1_data.RData to file1_data.Rdata

2012-07-20 Thread Yasir Kaheil
use function system in R
for (x in 1:100){
system (paste(rename file0.,x,_data.RData file,x,_data.Rdata,sep=))
}

-
Yasir Kaheil

--
View this message in context: 
http://r.789695.n4.nabble.com/change-file-name-from-file0-1-data-RData-to-file1-data-Rdata-tp4637135p4637136.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] Creating Multiple Repeating samples and Cross Correlating them.

2012-07-20 Thread Rbc1356
I'm having a lot of trouble getting this done and nothing I've written has
been remotely successful. Basically I have about 64 binary data sets stored
as vectors, with 72900 entries in each. I am not very familiar with cross
correlations, but I was advised that I should create 1000 randomized date
sets (used sample() for that) to use as a control, then compare the cross
correlation of my real data sets to the controls. First of all I am having
trouble creating 1000 different sample()'s of the same vector, and then I am
not all too sure how to cross correlate them Nothing seems to produce what I
want, all I get are copies of the same sample() of the vector repeated over
and over. Please Help?

--
View this message in context: 
http://r.789695.n4.nabble.com/Creating-Multiple-Repeating-samples-and-Cross-Correlating-them-tp4637131.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] Truncating (rounding down) to nearest half hour.

2012-07-20 Thread arun
Hi,

Just a variant of creating rounded down minutes vector:


minutes-ifelse(ceiling(as.numeric(substr(t[,1],15,16))30),30,00)
#or
minutes1-ifelse(floor(as.numeric(substr(dat1[,1],15,16))30),30,00)
identical(minutes,minutes1)
#[1] TRUE


minutes3-floor( as.numeric( substr( t[,1], 15, 16 ) ) / 30 ) * 30
 minutes3- ifelse( minutes == 30, 30, 00 )
 identical(minutes,minutes3)
[1] TRUE


A.K.

- Original Message -
From: Rainer Schuermann rainer.schuerm...@gmx.net
To: r-help@r-project.org; APOCooter mikeedinge...@gmail.com
Cc: 
Sent: Friday, July 20, 2012 1:07 AM
Subject: Re: [R] Truncating (rounding down) to nearest half hour.

My amateur approach:

I put your data in a dataframe called t:
 head( t )
                 Date Score
1 2008-05-01 08:58:00    80
2 2008-05-01 13:31:00    11
3 2008-05-01 16:35:00    81
4 2008-05-01 23:20:00   152
5 2008-05-02 01:01:00   130
6 2008-05-02 03:35:00   122

Then I created a vector with rounded down minutes:
 minutes - floor( as.numeric( substr( t[,1], 15, 16 ) ) / 30 ) * 30
 minutes - ifelse( minutes == 30, 30, 00 )
 head( minutes )
[1] 30 30 30 00 00 30

Next the replacement:
 tc - t
 tc[,1] - as.character( t[,1] )
 substr( tc[,1] , 15, 16 ) - minutes
 tc[,1] - as.POSIXct( tc[,1] )
 head( tc )
                 Date Score
1 2008-05-01 08:30:00    80
2 2008-05-01 13:30:00    11
3 2008-05-01 16:30:00    81
4 2008-05-01 23:00:00   152
5 2008-05-02 01:00:00   130
6 2008-05-02 03:30:00   122

Is that what you were looking for?

Rgds,
Rainer



 Original-Nachricht 
 Datum: Thu, 19 Jul 2012 10:03:23 -0700 (PDT)
 Von: APOCooter mikeedinge...@gmail.com
 An: r-help@r-project.org
 Betreff: [R] Truncating (rounding down) to nearest half hour.

 I couldn't find anything in the chron or timeDate packages, and a good
 search
 yielded rounding to the nearest half hour, which I don't want.
 
 The data:
 
 structure(list(Date = structure(c(1209625080, 1209641460, 1209652500, 
 1209676800, 1209682860, 1209692100, 1209706980, 1209722580, 1209726300, 
 1209739620, 1209762780, 1209765720, 1209770520, 1209791040, 1209812580, 
 1209829920, 1209837180, 1209848160, 1209854640, 1209859440, 1209870780, 
 1209887760, 1209901080, 1209921660, 1209929280, 1209945600, 1209957240, 
 1209980280, 1210001760, 1210017000, 1210021140, 1210034820, 1210042800, 
 1210048980, 1210061520, 1210074480, 1210081200, 1210089300, 1210095960, 
 1210104120, 1210110900, 1210110900, 1210118400, 1210126980, 1210134180, 
 1210142640, 1210156080, 1210164180, 1210176840, 1210183740), class =
 c(POSIXct, 
 POSIXt), tzone = ), Score = c(80L, 11L, 81L, 152L, 130L, 
 122L, 142L, 20L, 1L, 31L, 93L, 136L, 128L, 112L, 48L, 57L, 92L, 
 108L, 100L, 107L, 81L, 37L, 47L, 70L, 114L, 125L, 99L, 46L, 108L, 
 106L, 111L, 75L, 75L, 136L, 36L, 13L, 35L, 71L, 105L, 113L, 116L, 
 116L, 94L, 130L, 102L, 19L, 1L, 33L, 78L, 89L)), .Names = c(Date, 
 Score), row.names = c(NA, 50L), class = data.frame)
 
 I'm trying to round the time down to the nearest half hour, without losing
 the date.  For example, 01/01/2009 1:29 would round to 01/01/2009 1:00,
 while 01/01/2009 1:31 would round to 01/01/2009 1:30.
 
 Any help is greately appreciated.
 
 --
 View this message in context:
 http://r.789695.n4.nabble.com/Truncating-rounding-down-to-nearest-half-hour-tp4637083.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.

-- 
---

Gentoo Linux with KDE

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] A graph with a positive y-axis intercept

2012-07-20 Thread Lekgatlhamang, lexi Setlhare
Dear Sarah and A.K.,
Yes, both your suggestions give me exactly what I wanted. Thank you kindly.
 
Lexi
 


 From: Sarah Goslee sarah.gos...@gmail.com

Cc: r-h...@stat.math.ethz.ch r-h...@stat.math.ethz.ch 
Sent: Thursday, July 19, 2012 11:36 PM
Subject: Re: [R] A graph with a positive y-axis intercept
  
Your question is not entirely clear to me, but I think you might want:

plot(X, Y, xlim=c(0, 25), ylim=c(0, 25))

Sarah

On Thu, Jul 19, 2012 at 5:31 PM, Lekgatlhamang, lexi Setlhare

 Dear all,

 I have a challenge with a supposedly simple graph (a scatter). I wanted to 
 use R to create a plot/graph with a positive y-axis intercept but it does not 
 seem to give me what I expect to see. As an example, I used the following 
 values

 X- c(0, 4, 8, 11)
 Y- c(8, 12, 15, 22)

 then I used the plot function in the base package

 plot(X, Y)

 This gave me a graph which might be mistaken to be starting from the origin, 
 and yet, using the same numbers in Excel, and plotting them using the scatter 
 option, I get a graph which clearly has a positive y-axis intercept.
 Please help.

 Lexi

-- 
Sarah Goslee
http://www.functionaldiversity.org
[[alternative HTML version deleted]]

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


Re: [R] Subsetting problem data, 2

2012-07-20 Thread Lib Gray
I'm getting this error message:

nms-names(data)[grep(vars,names(data))]
Warning message:
In grep(vars, names(data)) :
  argument 'pattern' has length  1 and only the first element will be used

Is there a way around this?


On Thu, Jul 19, 2012 at 6:17 PM, Rui Barradas ruipbarra...@sapo.pt wrote:

 Hello,

 I guess so, and I can save you some typing.

 vars - sort(apply(expand.grid(L, 1:8, 1:2), 1, paste, collapse=))


 Then use it and see the result.

 Rui Barradas

 Em 20-07-2012 00:00, Lib Gray escreveu:

 The variables are actually L11, L12, L21, L22, ... , L81, L82. Would just
 creating a vector c(L11,... ,L82) be fine? (I'm about to try it, but I
 wanted to check to see if that was going to be a big issue).

 On Thu, Jul 19, 2012 at 3:27 PM, Rui Barradas ruipbarra...@sapo.pt
 wrote:

  Hello,

 Try the following. The data is your example of Patient A through E, but
 from the output of dput().

 dat - structure(list(Patient = structure(c(1L, 1L, 1L, 1L, 1L, 2L,
 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 5L, 5L, 5L), .Label = c(A,
 B, C, D, E), class = factor), Cycle = c(1L, 2L, 3L,
 4L, 5L, 1L, 2L, 1L, 3L, 4L, 5L, 1L, 2L, 4L, 5L, 1L, 2L, 3L),
  V1 = c(0.4, 0.3, 0.3, 0.4, 0.5, 0.4, 0.4, 0.9, 0.3, NA, 0.4,
  0.2, 0.5, 0.6, 0.5, 0.1, 0.5, 0.4), V2 = c(0.1, 0.2, NA,
  NA, 0.2, NA, NA, 0.9, 0.5, NA, NA, 0.5, 0.7, 0.4, 0.5, NA,
  0.3, 0.3), V3 = c(0.5, 0.5, 0.6, 0.4, 0.5, NA, NA, 0.9, 0.6,
  NA, NA, NA, NA, NA, NA, NA, NA, NA), V4 = c(1.5, 1.6, 1.7,
  1.8, 1.5, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
  NA), V5 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
  NA, NA, NA, NA, NA, NA)), .Names = c(Patient, Cycle,
 V1, V2, V3, V4, V5), class = data.frame, row.names = c(NA,
 -18L))

 dat

 nms - names(dat)[grep(^V[1-9]$, names(dat))]
 dd - split(dat, dat$Patient)
 fun - function(x) any(is.na(x))  any(!is.na(x))
 ix - sapply(dd, function(x) Reduce(`|`, lapply(x[, nms], fun)))

 dd[ix]
 do.call(rbind, dd[ix])


 I'm assuming that the variables names are as posted, V followed by one
 single digit 1-9. To keep the Patients with complete cases just negate
 the
 index 'ix', it's a logical index.
 Note also that dput() is the best way of posting a data example.

 Hope this helps,

 Rui Barradas

 Em 19-07-2012 15:15, Lib Gray escreveu:

  Hello,

 I didn't give enough information when I sent an query before, so I'm
 trying
 again with a more detailed explanation:

 In this data set, each patient has a different number of measured
 variables
 (they represent tumors, so some people had 2 tumors, some had 5, etc).
 The
 problem I have is that often in later cycles for a patient, tumors that
 were originally measured are now missing (or a new tumor showed up).
 We
 assume there are many different reasons for why a tumor would be
 measured
 in one cycle and not another, and so I want to subset OUT the problem
 patients to better study these patterns.

 An example:

 Patient  Cycle  V1  V2  V3  V4  V5
 A  1  0.4  0.1  0.5  1.5  NA
 A  2  0.3  0.2  0.5  1.6  NA
 A  3  0.3  NA  0.6  1.7  NA
 A  4  0.4  NA  0.4  1.8  NA
 A  5  0.5  0.2  0.5  1.5  NA

 I want to keep patient A; they have 4 measured tumors, but tumor 2 is
 missing data for cycles 3 and 4

 B  1  0.4  NA  NA  NA  NA
 B  2  0.4  NA  NA  NA  NA

 I do not want to keep patient B; they have 1 tumor that is measure
 consistently in both cycles

 C  1  0.9  0.9  0.9  NA  NA
 C  3  0.3  0.5  0.6  NA  NA
 C  4  NA  NA  NA  NA  NA
 C  5  0.4  NA  NA  NA  NA

 I do want to keep patient C; all their data is missing for cycle 4 and
 cycle 5 only measured one tumor

 D  1  0.2  0.5  NA  NA  NA
 D  2  0.5  0.7  NA  NA  NA
 D  4  0.6  0.4  NA  NA  NA
 D  5  0.5  0.5  NA  NA  NA

 I do not want patient D, their two tumors were measured each cycle

 E  1  0.1  NA  NA  NA  NA
 E  2  0.5  0.3  NA  NA  NA
 E  3  0.4  0.3  NA  NA  NA

 I DO want patient E; they only had one tumor register in Cycle 1, but
 cycles 2 and 3 had two tumors.


 Thanks for any help!

  [[alternative HTML version deleted]]

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-helphttps://stat.ethz.ch/mailman/**listinfo/r-help
 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.htmlhttp://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] R code for to check outliers

2012-07-20 Thread Angus Wallace
Really appreciate the discussion on outliers.

I come from an engineering signal processing background, and my thinking
has generally been that an outlier is outside a threshold of

   - distance from the mean
   - rarity

that we don't need/want to capture in whatever model we're building.

In my recent work (bioinformatics), I've seen that it's common to Winsorize
the data. I am a bit uncomfortable with this, though it seems to be
standard practice. Do people have thoughts here?

Cheers,
-Gus

[[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] Subsetting problem data, 2

2012-07-20 Thread Lib Gray
I'm still getting the message (if this is what you were suggesting I try).
The data set I'm using has many more columns other than these variables;
could that be a problem? I didn't think it would affect it.

pattern - L[1-8][12]
 nms-names(data)[grep(vars,names(data))]
Warning message:
In grep(vars, names(data)) :
  argument 'pattern' has length  1 and only the first element will be used


On Thu, Jul 19, 2012 at 6:55 PM, Rui Barradas ruipbarra...@sapo.pt wrote:

 Hello,

 Sorry, forgot about that. It's trickier to write code without a dataset to
 test it.

 Try

 pattern - L[1-8][12]

 and after the grep print nms to see if it's right.

 Rui Barradas

 Em 20-07-2012 00:33, Lib Gray escreveu:

 I'm getting this error message:

 nms-names(data)[grep(vars,**names(data))]
 Warning message:
 In grep(vars, names(data)) :
argument 'pattern' has length  1 and only the first element will be
 used

 Is there a way around this?


 On Thu, Jul 19, 2012 at 6:17 PM, Rui Barradas ruipbarra...@sapo.pt
 wrote:

  Hello,

 I guess so, and I can save you some typing.

 vars - sort(apply(expand.grid(L, 1:8, 1:2), 1, paste, collapse=))


 Then use it and see the result.

 Rui Barradas

 Em 20-07-2012 00:00, Lib Gray escreveu:

  The variables are actually L11, L12, L21, L22, ... , L81, L82. Would
 just
 creating a vector c(L11,... ,L82) be fine? (I'm about to try it, but I
 wanted to check to see if that was going to be a big issue).

 On Thu, Jul 19, 2012 at 3:27 PM, Rui Barradas ruipbarra...@sapo.pt
 wrote:

   Hello,

 Try the following. The data is your example of Patient A through E, but
 from the output of dput().

 dat - structure(list(Patient = structure(c(1L, 1L, 1L, 1L, 1L, 2L,
 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 5L, 5L, 5L), .Label = c(A,
 B, C, D, E), class = factor), Cycle = c(1L, 2L, 3L,
 4L, 5L, 1L, 2L, 1L, 3L, 4L, 5L, 1L, 2L, 4L, 5L, 1L, 2L, 3L),
   V1 = c(0.4, 0.3, 0.3, 0.4, 0.5, 0.4, 0.4, 0.9, 0.3, NA, 0.4,
   0.2, 0.5, 0.6, 0.5, 0.1, 0.5, 0.4), V2 = c(0.1, 0.2, NA,
   NA, 0.2, NA, NA, 0.9, 0.5, NA, NA, 0.5, 0.7, 0.4, 0.5, NA,
   0.3, 0.3), V3 = c(0.5, 0.5, 0.6, 0.4, 0.5, NA, NA, 0.9, 0.6,
   NA, NA, NA, NA, NA, NA, NA, NA, NA), V4 = c(1.5, 1.6, 1.7,
   1.8, 1.5, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
   NA), V5 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
   NA, NA, NA, NA, NA, NA)), .Names = c(Patient, Cycle,
 V1, V2, V3, V4, V5), class = data.frame, row.names = c(NA,
 -18L))

 dat

 nms - names(dat)[grep(^V[1-9]$, names(dat))]
 dd - split(dat, dat$Patient)
 fun - function(x) any(is.na(x))  any(!is.na(x))
 ix - sapply(dd, function(x) Reduce(`|`, lapply(x[, nms], fun)))

 dd[ix]
 do.call(rbind, dd[ix])


 I'm assuming that the variables names are as posted, V followed by one
 single digit 1-9. To keep the Patients with complete cases just negate
 the
 index 'ix', it's a logical index.
 Note also that dput() is the best way of posting a data example.

 Hope this helps,

 Rui Barradas

 Em 19-07-2012 15:15, Lib Gray escreveu:

   Hello,

 I didn't give enough information when I sent an query before, so I'm
 trying
 again with a more detailed explanation:

 In this data set, each patient has a different number of measured
 variables
 (they represent tumors, so some people had 2 tumors, some had 5, etc).
 The
 problem I have is that often in later cycles for a patient, tumors
 that
 were originally measured are now missing (or a new tumor showed up).
 We
 assume there are many different reasons for why a tumor would be
 measured
 in one cycle and not another, and so I want to subset OUT the
 problem
 patients to better study these patterns.

 An example:

 Patient  Cycle  V1  V2  V3  V4  V5
 A  1  0.4  0.1  0.5  1.5  NA
 A  2  0.3  0.2  0.5  1.6  NA
 A  3  0.3  NA  0.6  1.7  NA
 A  4  0.4  NA  0.4  1.8  NA
 A  5  0.5  0.2  0.5  1.5  NA

 I want to keep patient A; they have 4 measured tumors, but tumor 2 is
 missing data for cycles 3 and 4

 B  1  0.4  NA  NA  NA  NA
 B  2  0.4  NA  NA  NA  NA

 I do not want to keep patient B; they have 1 tumor that is measure
 consistently in both cycles

 C  1  0.9  0.9  0.9  NA  NA
 C  3  0.3  0.5  0.6  NA  NA
 C  4  NA  NA  NA  NA  NA
 C  5  0.4  NA  NA  NA  NA

 I do want to keep patient C; all their data is missing for cycle 4 and
 cycle 5 only measured one tumor

 D  1  0.2  0.5  NA  NA  NA
 D  2  0.5  0.7  NA  NA  NA
 D  4  0.6  0.4  NA  NA  NA
 D  5  0.5  0.5  NA  NA  NA

 I do not want patient D, their two tumors were measured each cycle

 E  1  0.1  NA  NA  NA  NA
 E  2  0.5  0.3  NA  NA  NA
 E  3  0.4  0.3  NA  NA  NA

 I DO want patient E; they only had one tumor register in Cycle 1, but
 cycles 2 and 3 had two tumors.


 Thanks for any help!

   [[alternative HTML version deleted]]

 __**
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/**listinfo/r-helphttps://stat.ethz.ch/mailman/listinfo/r-help
 

Re: [R] Changing ungrouped cases to grouped cases

2012-07-20 Thread arun
Hi,

Try this:

dat1-read.table(text=
y    A  B  C
0    1    1  2
0    1    2  1
1    1    1  2
0    1    1  2
1    1    1  2
1    1    2  1
0    1    2  2
,sep=,header=TRUE)
 dat2-aggregate(y~A+B+C,data=dat1,sum)
 dat2-dat2[,c(4,1:3)]
dat3-dat2[with(dat2,rev(order(y,A,B,C))),]
 dat3
  y A B C
2 2 1 1 2
1 1 1 2 1
3 0 1 2 2

A.K.





- Original Message -
From: Christopher Desjardins cddesjard...@gmail.com
To: R help r-help@r-project.org
Cc: 
Sent: Thursday, July 19, 2012 8:34 PM
Subject: [R] Changing ungrouped cases to grouped cases

Hi,
I have my data the following way:

y     A   B   C
0     1    1   2
0     1    2   1
1     1    1   2
0     1    1   2
1     1    1   2
1     1    2   1
0     1    2   2
.
.
.
And so on.  How can I make my data look like the following:
y   A  B  C
2   1   1  2
1   1   2  1
0   1   2   2
.
.
.

In other words how can I change my ungrouped cases into grouped cases?
Thanks!
Chris

    [[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] (no subject)

2012-07-20 Thread Suman Pramanik
Dear Sir,
I am using timereg package for analysing competing risk data.
There are two options in method: proportional and additive.
How to choose among the methods for analysis and how to check which of the
methods are appropiate for the specific data.
Tanks a lot.
Dr Suman Kumar

[[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] Last answer

2012-07-20 Thread arun
Hi,

Try this:

 ans-function(){.Last.value}
 2+3
[1] 5
 ans()
[1] 5

A.K.



- Original Message -
From: darnold dwarnol...@suddenlink.net
To: r-help@r-project.org
Cc: 
Sent: Friday, July 20, 2012 1:12 AM
Subject: [R] Last answer

Hi,

In Matlab, I can access the last computation as follows:

 2+3

ans =

     5

 ans

ans =

     5

Anything similar in R?

David



--
View this message in context: 
http://r.789695.n4.nabble.com/Last-answer-tp4637151.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] Subsetting problem data, 2

2012-07-20 Thread Chris Campbell
Hi!

# toy data   

toyData - data.frame(x = 1:4, y = 5:8, xy = 9:12, z = 13:16)
vars - c(x, z)  

# pattern is an argument of grep  

args(grep)  

# pattern must only consist of a single element 
# otherwise only the first element is used  

grep(pattern = vars, x = names(toyData))   

# one way to do this - a loop 
# create a vector to collect the output of each call
 
toyColIndexList - vector(length = length(vars), mode = list)

# grep each element in turn 

for (i in seq_along(vars)) {  
toyColIndexList[[i]] - grep(pattern = vars[i], x = names(toyData)) 
}  
 
# combine all of the answers 
 
toyColIndex - unlist(toyColIndexList) 

# remove duplicated columns if present

toyColIndex - toyColIndex[!duplicated(toyColIndex)] 
 
# select the elements we want

toyData[, toyColIndex] 

  
# alternatively we could use regular expressions   
 
grep(pattern = (x|z), x = names(toyData))
 
# hope this helps

Best wishes

Chris

Chris Campbell
Mango Solutions
Data Analysis that Delivers
http://www.mango-solutions.com
+44 (0) 1249 705 450  


-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
Behalf Of Lib Gray
Sent: 20 July 2012 01:17
To: Rui Barradas
Cc: r-help
Subject: Re: [R] Subsetting problem data, 2

I'm still getting the message (if this is what you were suggesting I try).
The data set I'm using has many more columns other than these variables; could 
that be a problem? I didn't think it would affect it.

pattern - L[1-8][12]
 nms-names(data)[grep(vars,names(data))]
Warning message:
In grep(vars, names(data)) :
  argument 'pattern' has length  1 and only the first element will be used


On Thu, Jul 19, 2012 at 6:55 PM, Rui Barradas ruipbarra...@sapo.pt wrote:

 Hello,

 Sorry, forgot about that. It's trickier to write code without a 
 dataset to test it.

 Try

 pattern - L[1-8][12]

 and after the grep print nms to see if it's right.

 Rui Barradas

 Em 20-07-2012 00:33, Lib Gray escreveu:

 I'm getting this error message:

 nms-names(data)[grep(vars,**names(data))]
 Warning message:
 In grep(vars, names(data)) :
argument 'pattern' has length  1 and only the first element will 
 be used

 Is there a way around this?


 On Thu, Jul 19, 2012 at 6:17 PM, Rui Barradas ruipbarra...@sapo.pt
 wrote:

  Hello,

 I guess so, and I can save you some typing.

 vars - sort(apply(expand.grid(L, 1:8, 1:2), 1, paste, 
 collapse=))


 Then use it and see the result.

 Rui Barradas

 Em 20-07-2012 00:00, Lib Gray escreveu:

  The variables are actually L11, L12, L21, L22, ... , L81, L82. 
 Would
 just
 creating a vector c(L11,... ,L82) be fine? (I'm about to try it, 
 but I wanted to check to see if that was going to be a big issue).

 On Thu, Jul 19, 2012 at 3:27 PM, Rui Barradas 
 ruipbarra...@sapo.pt
 wrote:

   Hello,

 Try the following. The data is your example of Patient A through 
 E, but from the output of dput().

 dat - structure(list(Patient = structure(c(1L, 1L, 1L, 1L, 1L, 
 2L, 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 5L, 5L, 5L), .Label = 
 c(A, B, C, D, E), class = factor), Cycle = c(1L, 2L, 
 3L, 4L, 5L, 1L, 2L, 1L, 3L, 4L, 5L, 1L, 2L, 4L, 5L, 1L, 2L, 3L),
   V1 = c(0.4, 0.3, 0.3, 0.4, 0.5, 0.4, 0.4, 0.9, 0.3, NA, 0.4,
   0.2, 0.5, 0.6, 0.5, 0.1, 0.5, 0.4), V2 = c(0.1, 0.2, NA,
   NA, 0.2, NA, NA, 0.9, 0.5, NA, NA, 0.5, 0.7, 0.4, 0.5, NA,
   0.3, 0.3), V3 = c(0.5, 0.5, 0.6, 0.4, 0.5, NA, NA, 0.9, 0.6,
   NA, NA, NA, NA, NA, NA, NA, NA, NA), V4 = c(1.5, 1.6, 1.7,
   1.8, 1.5, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
   NA), V5 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
   NA, NA, NA, NA, NA, NA)), .Names = c(Patient, Cycle, 
 V1, V2, V3, V4, V5), class = data.frame, row.names = 
 c(NA,
 -18L))

 dat

 nms - names(dat)[grep(^V[1-9]$, names(dat))] dd - split(dat, 
 dat$Patient) fun - function(x) any(is.na(x))  any(!is.na(x)) ix 
 - sapply(dd, function(x) Reduce(`|`, lapply(x[, nms], fun)))

 dd[ix]
 do.call(rbind, dd[ix])


 I'm assuming that the variables names are as posted, V followed by 
 one single digit 1-9. To keep the Patients with complete cases 
 just negate the index 'ix', it's a logical index.
 Note also that dput() is the best way of posting a data example.

 Hope this helps,

 Rui Barradas

 Em 19-07-2012 15:15, Lib Gray escreveu:

   Hello,

 I didn't give enough information when I sent an query before, so 
 I'm trying again with a more detailed explanation:

 In this data set, each patient has a different number of measured 
 variables (they represent tumors, so some people had 2 tumors, 
 some had 5, etc).
 The
 problem I have is that often in later cycles for a patient, 
 tumors that were originally measured are now missing (or a new 
 tumor showed up).
 We
 assume there are many different reasons for why a tumor would be 
 

Re: [R] generic functions question in building a new package

2012-07-20 Thread Uwe Ligges



On 19.07.2012 22:44, zgu9 wrote:

Hi all,

I'm trying to build a package 'Thiele', but run into problems with generic
functions.
I have a class Benefit, and defined the function print.Benefit.
When I try R CMD check Thiele in X11, this warning came out
---
Functions/methods with usage in documentation object 'Benefit' but not in
code:
   print



You need

\S3method{print}{Benefit}(...)

in the usage line rather than

print(...)

Best,
Uwe Ligges



---
Do I need to write a new print Rd document into R file in that Package? If
so, how can I write it?
I read the section 7 Generic functions and methods in the 'Writing R
Extensions' manual. But I'm still not clear how to figure out my problem.

Thanks!
Eric

--
View this message in context: 
http://r.789695.n4.nabble.com/generic-functions-question-in-building-a-new-package-tp4637129.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] function for inverse normal transformation

2012-07-20 Thread carol white
Hi,
What is the function for inverse normal transformation?

Thanks,

Carol

[[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] Creating a pdf with layers?

2012-07-20 Thread Christoph Scherber
Dear all,

Is it possible to create a pdf file with layers using the pdf() device in R?

Many thanks for your help!
Christoph


(using R 2.15.1 on Windows 7 64-Bit)



-- 
PD Dr Christoph Scherber
Georg-August University Goettingen
Department of Crop Science
Agroecology
Grisebachstrasse 6
D-37077 Goettingen
Germany
phone 0049 (0)551 39 8807
fax 0049 (0)551 39 8806
http://www.gwdg.de/~cscherb1
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] why does this simple example NOT work?

2012-07-20 Thread Duncan Murdoch

On 12-07-20 2:50 AM, Martin Ivanov wrote:

  Hello,

I want to create and save objects in a loop, but this is precluded by the 
following obstacle:
this part of the script fails to work:

assign(x=paste(a, 1, sep=), value=1);
save(paste(a, 1, sep=), file=paste(paste(a, 1, sep=), .RData, sep=))


paste(a, 1, sep=) is equivalent to a1.  So you are saving the 
string a1 to the file a1.RData.  You presumably want to save the 
value of the variable a1, not that string.  (Though you didn't say that. 
 It's helpful to describe what you don't like, don't just say it fails 
to work.)


If you have the name of a variable in a string, you can use get() to get 
the value.  So I think this code would do what you want:


name - paste0(a, 1)
assign(x=name, value=1)
save(get(name), file=paste0(name, .RData))

(If you don't have the paste0 function, you should upgrade to the 
current R version.)


Duncan Murdoch



Do you know any workaround? I am already out of ideas.

Best regards,

Martin

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

2012-07-20 Thread Duncan Murdoch

On 12-07-20 6:21 AM, carol white wrote:

Hi,
What is the function for inverse normal transformation?


qnorm

Duncan Murdoch


Thanks,

Carol

[[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 for inverse normal transformation

2012-07-20 Thread carol white
Thanks for your reply.

So to derive it from a given data set, is the following correct to do?

my_data.p =2*pnorm(abs(my_data),lower.tail=FALSE)

my_data.q = qnorm(my_data.p)

Cheers,



 From: Duncan Murdoch murdoch.dun...@gmail.com

Cc: r-h...@stat.math.ethz.ch r-h...@stat.math.ethz.ch 
Sent: Friday, July 20, 2012 1:23 PM
Subject: Re: [R] function for inverse normal transformation

On 12-07-20 6:21 AM, carol white wrote:
 Hi,
 What is the function for inverse normal transformation?

qnorm

Duncan Murdoch

 Thanks,

 Carol

     [[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] conditional subset and reorder dataframe rows

2012-07-20 Thread Iain Gallagher
Hi List

I have a dataframe (~1,200,000 rows deep) and I'd like to conditionally reorder 
groups of rows in this dataframe. 

I would like to reorder any rows where the Chr.Strand column contains a '-' but 
reorder within subsets delineated by the Probe.Set.Name column.

# toy example 

library(plyr)

negStrandGene - data.frame(Probe.Set.Name = rep('ENSMUSG0022174_at', 6), 
Chr = rep(14,6), Chr.Strand = rep('-', 6), Chr.From = c(54873546, 54873539, 
54873533, 54873529, 54873527, 54873416), Probe.X = 
c(388,1634,2141,2305,882,960), Probe.Y = c(2112, 1773, 1045, 862, 971, 2160))

posStrandGene - data.frame(Probe.Set.Name = rep('ENSMUSG0047459_at', 6), 
Chr = rep(2, 6), Chr.Strand = rep('+', 6), Chr.From = c(155062277, 155062304, 
155062305, 155062309, 155062326, 155062531), Probe.X = c(428, 1681, 2058, 1570, 
1293, 2125), Probe.Y = c(1484, 2090, 893, 1082, 1435, 1008))

mapping - rbind (negStrandGene, posStrandGene)

# define a function to do what we want
revSort - function(df){

  if (unique(df$Chr.Strand == '-')) 
  return (df[order(df$Chr.From), ])
  else return (df)

}


# split the data with plyr, apply the function and recombine
test2 - ddply(mapping, .(Probe.Set.Name), function(df) revSort(df)) # ok, cool 
works

So here the rows with the '-' if Chr.Strand are reordered whilst those with '+' 
are not.

My initial attempt using plyr is very inefficient and I wondered if someone 
could suggest something better.

Best

Iain

[[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] function for inverse normal transformation

2012-07-20 Thread Duncan Murdoch

On 12-07-20 7:36 AM, carol white wrote:

Thanks for your reply.

So to derive it from a given data set, is the following correct to do?

my_data.p =2*pnorm(abs(my_data),lower.tail=FALSE)

my_data.q = qnorm(my_data.p)


I don't know what you're trying to do, but that doesn't look like it 
does something sensible.  It would take a value like 2, compute the p to 
be 0.045, and return the corresponding quantile of the normal 
distribution, i.e. -1.69 or so.  I don't know why you'd want to do that.


Duncan Murdoch



Cheers,



  From: Duncan Murdoch murdoch.dun...@gmail.com
To: carol white wht_...@yahoo.com
Cc: r-h...@stat.math.ethz.ch r-h...@stat.math.ethz.ch
Sent: Friday, July 20, 2012 1:23 PM
Subject: Re: [R] function for inverse normal transformation

On 12-07-20 6:21 AM, carol white wrote:

Hi,
What is the function for inverse normal transformation?


qnorm

Duncan Murdoch


Thanks,

Carol

 [[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] rcspline.problem

2012-07-20 Thread Bart Ferket
Dear Peter,

This indeed resolves the problem.

Many thanks. My apologies for not starting a new thread. I am a new R user
and not yet fully integrated into the R community.

Kind regards, 

Bart Ferket



--
View this message in context: 
http://r.789695.n4.nabble.com/rcspline-problem-tp3501627p4637181.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] Forced inclusion of varaibles in validate command as well as step

2012-07-20 Thread D. Lofaro
Dear prof. Harrell,

I'm not able to use the force option with fastbw, here an example of the error 
I've got (dataset stagec rpart package):

 fitstc - cph(Surv(stagec$pgtime,stagec$pgstat) ~ age + eet + g2 + grade + 
 gleason + ploidy, data=stagec)
 fbwstc - fastbw(fitstc,rule=aic,type=individual)
 fbwstc

 Deleted Chi-Sq d.f. P  Residual d.f. P  AIC
 eet 0.19   10.6610 0.19 10.6610 -1.81
 age 0.74   10.3889 0.93 20.6266 -3.07
 gleason 1.91   10.1673 2.84 30.4167 -3.16

Approximate Estimates after Deleting Factors

   CoefS.E. Wald Z P
g2 -0.03498 0.01843 -1.898 5.774e-02
grade   1.71963 0.35008  4.912 9.012e-07
ploidy  0.66700 0.25084  2.659 7.836e-03

Factors in Final Model

[1] g2 grade  ploidy
 fbwstc - fastbw(fitstg,rule=aic,type=individual,force=2)
Warning message:
In fastbw(fitstg, rule = aic, type = individual, force = 2) :
  force probably does not work unless type=individual
 fbwstc

Parameters forced into all models:
 eet

 Deleted Chi-Sq d.f. P  Residual d.f. P  AIC
 age 0.72   10.3960 0.72 10.3960 -1.28
 eet 0.21   10.6432 0.93 20.6266 -3.07
 gleason 1.91   10.1673 2.84 30.4167 -3.16

Approximate Estimates after Deleting Factors

   CoefS.E. Wald Z P
g2 -0.03498 0.01843 -1.898 5.774e-02
grade   1.71963 0.35008  4.912 9.012e-07
ploidy  0.66700 0.25084  2.659 7.836e-03

Factors in Final Model

[1] g2 grade  ploidy

Sorry if I bother you

best regards

Danilo Lofaro
De informatie opgenomen in dit bericht kan vertrouwelijk zijn en is uitsluitend 
bestemd voor de geadresseerde. Indien u dit bericht onterecht ontvangt, wordt u 
verzocht de inhoud niet te gebruiken en de afzender direct te informeren door 
het bericht te retourneren. Het Academisch Medisch Centrum is een 
publiekrechtelijke rechtspersoon in de zin van de W.H.W. (Wet Hoger Onderwijs 
en Wetenschappelijk Onderzoek) en staat geregistreerd bij de Kamer van 
Koophandel voor Amsterdam onder nr. 34362777. De Algemene Inkoop Voorwaarden 
van het AMC zijn van toepassing op en maken integraal onderdeel uit van alle 
rechtsbetrekkingen, daaronder mede verstaan alle inkoop opdrachten en 
overeenkomsten, tussen AMC en derden. Deze voorwaarden zijn te raadplegen op 
www.amc.nl en worden op verzoek toegezonden.



This message may contain confidential information and is...{{dropped:10}}

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

2012-07-20 Thread carla moreira

Hi,

I would like to evaluate a function, with 3 arguments, for instance, 

myfunc-function(a,b,c) { sqrt(a)-exp(b)+4*c
}

How to execute  myfunc(x,y,z), for all x, all y and all z, where x,y,z are
vectors?

Thank you very much in advance



--
View this message in context: 
http://r.789695.n4.nabble.com/Execute-a-function-tp4637182.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Creating a pdf with layers?

2012-07-20 Thread Prof Brian Ripley

On 20/07/12 12:07, Christoph Scherber wrote:

Dear all,

Is it possible to create a pdf file with layers using the pdf() device in R?


No.  Is it possible to specify layers in the R graphics language or any 
device?  (From what I understand by 'layers', no.)


The author of pdf()



Many thanks for your help!
Christoph


(using R 2.15.1 on Windows 7 64-Bit)




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

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


Re: [R] Line chart with a double matrix

2012-07-20 Thread ramonovelar
Thanks David for the solution, I paste here the code with a few changes:

par(xpd = NA, oma = c(5, 0, 0, 0))
matplot(t(buffersump[,1:6]), type=l, xaxt=n, lwd=3, ylab=Porcentajes de
respuestas afirmativas)
axis(1, 1:6, colnames(buffersump2))
legend(bottomright, legend=rownames(buffersump2), lty=1:5, col=1:5,  inset
= c(0, -0.7))



--
View this message in context: 
http://r.789695.n4.nabble.com/Line-chart-with-a-double-matrix-tp4637047p4637177.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] TukeyHSD not working

2012-07-20 Thread Michael Eisenring
Dear r-help members.
I would like to compare species numbers of moths between eight different 
forests (each sampled for six nights). I would like to do a nested anova to 
compare species numbers between forests and nights.
For more site specific details I wanted to do a Tukey test (TukeyHSD). 
Unfortunately this test is not working (message:  nicht anwendbare Methode für 
'TukeyHSD' auf Objekt der Klasse c('anova', 'data.frame ) /(message: method 
not appropriate for TukeyHSD for objects of the class anova, data.frame).
I tried it also with the aov instead of the anova function but neither the 
aov nor the anova approach is working.
Could anyone help me with this  problem?
Attached you can find the txt.file and the r-script.
I'm working with a MAC OSX snow leopard, with R-Studio 0.96 316.

Thank you very much.
Michael 
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Crosstab with Average and Count

2012-07-20 Thread vioravis
I have the following data:

x - as.factor(c(1,1,1,2,2,2,3,3,3))
y - as.factor(c(10,10,10,20,20,20,30,30,30))
z - c(100,100,NA,200,200,200,300,300,300)

I could create the cross tab of x and y with Sum of z as its elements using
the xtabs function as follows:

# X Vs. Y with Sum Z

xtabs(z ~ x + y)

   y
x10  20  30
  1 200   0   0
  2   0 600   0
  3   0   0 900

How do I replace the sum with average and count so that I can get the
following outputs??

# X Vs. Y with Average of Z
   y
x  10  20  30
  1100 0   0
  20   200 0
  30   0   300

# X Vs. Y with Count Z
  y
x10  20  30
 12   0   0
 20   3   0
 30   0   3

Would appreciate any help on these? Thank you.

Ravi





--
View this message in context: 
http://r.789695.n4.nabble.com/Crosstab-with-Average-and-Count-tp4637180.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] qplot/ggplot

2012-07-20 Thread Anamika K
Hello,
I have a file like this (just a snapshot) where I have numerical
values for various genes, I want a line plot with shading (may be
using smooth ?) using qplot or ggplot :

Gene1 10 14 12 23 11 11 33 1 ..(multiple columns)
Gene2 4 2 1 1 3 4 1 2 .
Gene3 2 5 7 5 6 89 7 3 ..
Gene4 1 9 8 3 90 8 8  .

Please help.

Ana

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

2012-07-20 Thread arun
HI,

Probably this might be helpful for you.

http://rgm2.lab.nig.ac.jp/RGM2/func.php?rd_id=fBasics:dist-nigFit
A.K.



- Original Message -
From: carol white wht_...@yahoo.com
To: r-h...@stat.math.ethz.ch r-h...@stat.math.ethz.ch
Cc: 
Sent: Friday, July 20, 2012 6:21 AM
Subject: [R] function for inverse normal transformation

Hi,
What is the function for inverse normal transformation?

Thanks,

Carol

    [[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 for inverse normal transformation

2012-07-20 Thread carol white
I have a continuous data. So to calculate the inverse normal transformation, I 
thought that I should first calculate the Z-score normalized data and then, 
calculate the p-value et the quantile transformation. Does this seem to be more 
sensible


my_data.p =2*pnorm(abs(scale(my_data)),lower.tail=FALSE)
my_data.q= qnorm(my_data.p)


The attached file shows the histogram of a small data set before 
transformation, the p-value generated from the Z-score normalized data and 
then, the qnorm-transformed data.

Thanks for your feedback,



 From: Duncan Murdoch murdoch.dun...@gmail.com
To: carol white wht_...@yahoo.com 
Cc: r-h...@stat.math.ethz.ch r-h...@stat.math.ethz.ch 
Sent: Friday, July 20, 2012 1:42 PM
Subject: Re: [R] function for inverse normal transformation
 
On 12-07-20 7:36 AM, carol white wrote:
 Thanks for your reply.

 So to derive it from a given data set, is the following correct to do?

 my_data.p =2*pnorm(abs(my_data),lower.tail=FALSE)

 my_data.q = qnorm(my_data.p)

I don't know what you're trying to do, but that doesn't look like it 
does something sensible.  It would take a value like 2, compute the p to 
be 0.045, and return the corresponding quantile of the normal 
distribution, i.e. -1.69 or so.  I don't know why you'd want to do that.

Duncan Murdoch


 Cheers,


 
   From: Duncan Murdoch murdoch.dun...@gmail.com
 To: carol white wht_...@yahoo.com
 Cc: r-h...@stat.math.ethz.ch r-h...@stat.math.ethz.ch
 Sent: Friday, July 20, 2012 1:23 PM
 Subject: Re: [R] function for inverse normal transformation

 On 12-07-20 6:21 AM, carol white wrote:
 Hi,
 What is the function for inverse normal transformation?

 qnorm

 Duncan Murdoch

 Thanks,

 Carol

      [[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.
attachment: tmp.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] function for inverse normal transformation

2012-07-20 Thread Rui Barradas
Hello,

No it's not correct, you are computing a what seems to be a two-tailed 
probabiity, so the inverse should account for it. Look closely: you take 
the absolute value, then the upper tail probability, then multiply 2 
into it. Reverse these steps to get the correct value.

# Helper function
equal - function(x, y, tol=.Machine$double.eps^0.5) all(abs(x - y)   tol)

m - rnorm(5)
p - 2*pnorm(abs(m), lower.tail=FALSE)
m2 - qnorm(p/2, lower.tail=FALSE)*sign(m)

equal(m, m2)

(The helper function is just to test floating point values computed 
differently for equality.)

Hope this helps,

Rui Barradas

Em 20-07-2012 12:36, carol white escreveu:
 Thanks for your reply.

 So to derive it from a given data set, is the following correct to do?

 my_data.p =2*pnorm(abs(my_data),lower.tail=FALSE)

 my_data.q = qnorm(my_data.p)

 Cheers,


 
   From: Duncan Murdoch murdoch.dun...@gmail.com

 Cc: r-h...@stat.math.ethz.ch r-h...@stat.math.ethz.ch
 Sent: Friday, July 20, 2012 1:23 PM
 Subject: Re: [R] function for inverse normal transformation

 On 12-07-20 6:21 AM, carol white wrote:
 Hi,
 What is the function for inverse normal transformation?
 qnorm

 Duncan Murdoch

 Thanks,

 Carol

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


[[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] Execute a function

2012-07-20 Thread Jessica Streicher
You mean executing the function for all combinations of values?
For example, if you have a-b-c-1:2
you would get back the values of

myfunc(1,1,1)
myfunc(1,1,2)
myfunc(1,2,1)
myfunc(1,2,2)
myfunc(2,1,1)
myfunc(2,1,2)
myfunc(2,2,1)
myfunc(2,2,2)

?

On 20.07.2012, at 13:05, carla moreira wrote:

 
 Hi,
 
 I would like to evaluate a function, with 3 arguments, for instance, 
 
 myfunc-function(a,b,c) { sqrt(a)-exp(b)+4*c
}
 
 How to execute  myfunc(x,y,z), for all x, all y and all z, where x,y,z are
 vectors?
 
 Thank you very much in advance
 
 
 
 --
 View this message in context: 
 http://r.789695.n4.nabble.com/Execute-a-function-tp4637182.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] Execute a function

2012-07-20 Thread Benno Pütz
Not quite sure what you are aiming at, but looking at ?mapply or ?expand.grid 
could be helpful

Benno

On Jul 20, 2012, at 1:05 PM, carla moreira wrote:

 
 Hi,
 
 I would like to evaluate a function, with 3 arguments, for instance, 
 
 myfunc-function(a,b,c) { sqrt(a)-exp(b)+4*c
}
 
 How to execute  myfunc(x,y,z), for all x, all y and all z, where x,y,z are
 vectors?
 
 Thank you very much in advance
 
 
 
 --
 View this message in context: 
 http://r.789695.n4.nabble.com/Execute-a-function-tp4637182.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.
 

Benno Pütz
Statistical Genetics
MPI of Psychiatry
Kraepelinstr. 2-10
80804 Munich, Germany
T: ++49-(0)89-306 22 222
F: ++49-(0)89-306 22 601




[[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] Defining a variable outside of optim or differential equation solver.

2012-07-20 Thread Thomas Petzoldt

Hi Tjun Kat,

you can define variables outside the ode function, but normally NOT 
state variables, because their values need to be updated by the solver 
during the simulation process.


But, if you want to block this for any debugging purposes and want to 
e.g. fix a derivative to a certain value, even this is possible. Note 
however that this is a very special case and I suspect that you don't 
want this.


Can you please tell, why you want to define states outside? I guess you 
want to emulate a feature that is already available in deSolve, e.g. 
forcings or events. In that case, please have a look into the 
documentation and one of the papers tutorial slides etc. that can be 
found on:


http://desolve.r-forge.r-project.org



Note also that your code contains 3 errors:

1) The call must be function(t, y, p), i.e. with p even if this is 
not required by the model, because ode needs this interface.


2) the closing parenthesis ) of list is missing.

3) dvdpol vs. vdpol

Hope it helps

Thomas Petzoldt


On 7/18/2012 3:59 AM, Tjun Kiat Teo wrote:
 This is applicable to either using optim  or the differential equation
 solver or any similar solver

 Suppose I want to use the differential equation solver and this is my 
code


 d-y[2]

 vdpol-function(t,y)
 {
 list(c(1,
 d,
 3,
 4
)
 }


 stiff-ode(y=rep(0,4),times=c(0,1),func=dvdpol,parms=1)


 The thing is I want d to be composed of one of state variables in the
 differential function vdopl. Can it be done ?

 tjun kiat

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

2012-07-20 Thread Jan van der Laan


See below for the complete mail to which I reply which was not sent to rhelp.


==

emptyexpandlist2-list(ne=0,l=array(NA, dim=c(1, 1000L)),len=1000L)

addexpandlist2-function(x,prev){
  if(prev$len==prev$ne){
n2-prev$len*2
prev - list(ne=prev$ne, l=array(prev$l, dim=c(1, n2)), len=n2)
  }
  prev$ne-prev$ne+1
  prev$l[prev$ne]-x
  return(prev)
}

compressexpandlist2-function(prev){
  return(prev$l[seq.int(prev$ne)])
}

h3-function(dotot){
  v-emptyexpandlist2
  for(i in 1:dotot){
v-addexpandlist2(FALSE,v)
  }
  return(compressexpandlist2(v))
}

===


The problem with your addexpandlist2 is that R in principle works with  
pass by value (certainly when you modify the objects you pass in your  
function as you do with prev). Therefore, when you pass your list to  
addexpendlist2 it makes a copy of the entire list.


You can avoid that by using environments that are passed by reference.  
The code below shows an example of this. If you would like to  
implement something like that I would recommend using reference  
classes (see ?ReferenceClasses) . Personally I don't find the 'messy  
code' that messy. You get used to it.


myvector - function(N = 1000) {
data - vector(list, N)
n- 0

append - function(d) {
n - n + 1
if (n  N) {
N - 2*N
length(data) - N
}
data[[n]] - d
}

length - function() {
return(n)
}

get - function() {
return(data[seq_len(n)])
}

return(list(append=append, length=length, get=get))
}


h4 - function(dotot){
v - myvector()
for(i in seq_len(dotot)) {
v$append(FALSE)
}
return(v$get())
}



system.time(h3(1E5))

   user  system elapsed
 22.846   0.536  23.407

system.time(h4(1E5))

   user  system elapsed
  0.700   0.000   0.702


Jan




Johan Henriksson maho...@areta.org schreef:


On Thu, Jul 19, 2012 at 5:02 PM, Jan van der Laan rh...@eoos.dds.nl wrote:


Johan,

Your 'list' and 'array doubling' code can be written much more efficient.

The following function is faster than your g and easier to read:

g2 - function(dotot) {
  v - list()
  for (i in seq_len(dotot)) {
v[[i]] - FALSE
  }
}



the reason for my highly convoluted code was to simulate a linked list - to
my knowledge a list() in R is not a linked list but a vector. I was
assuming that R would not copy the entire memory into each sublist but
rather keep a pointer. if this had worked, it would also be possible to
create fancier data structures like trees and heaps

http://en.wikipedia.org/wiki/Linked_list
http://en.wikipedia.org/wiki/Heap_%28data_structure%29





In the following line in you array doubling function

  prev$l-rbind(prev$l,matrix(**ncol=1,nrow=nextsize))

you first create a new array: the second half of your new array. Then
rbind creates a new array and has to copy the contents of both into this
new array. The following routine is much faster and almost scales linearly
(see below):

h2 - function(dotot) {
  n - 1000L
  v - array(NA, dim=c(1, n))
  for(i in seq_len(dotot)) {
if (i  n) {
  n - 2*n
  v - array(v, dim=c(1, n))
}
v[, i] - TRUE
  }
  return(v[, seq_len(i)])
}




that's blazingly fast! thanks! I've also learned some nice optimization
tricks, like L, and seq.int, and the resizing with array...

given that this works, as you see, it's rather messy to use on a day-to-day
basis. my goal was next to hide to this in a couple of convenient functions
(emptyexpandlist and addexpandlist) so that this can be reused without
cluttering otherwise fine code. but the overhead of the function call, and
the tuple, seems to kill off the major advantage. that said, I present
these convenience functions below:

==

emptyexpandlist2-list(ne=0,l=array(NA, dim=c(1, 1000L)),len=1000L)

addexpandlist2-function(x,prev){
  if(prev$len==prev$ne){
n2-prev$len*2
prev - list(ne=prev$ne, l=array(prev$l, dim=c(1, n2)), len=n2)
  }
  prev$ne-prev$ne+1
  prev$l[prev$ne]-x
  return(prev)
}

compressexpandlist2-function(prev){
  return(prev$l[seq.int(prev$ne)])
}

h3-function(dotot){
  v-emptyexpandlist2
  for(i in 1:dotot){
v-addexpandlist2(FALSE,v)
  }
  return(compressexpandlist2(v))
}

===

I haven't checked the scaling but take it works as it should. the constant
factor is really bad though:

dotot=5

system.time(f(dotot))

   user  system elapsed
  5.250   0.020   5.279

system.time(h(dotot))

   user  system elapsed
  2.650   0.060   2.713

system.time(h2(dotot))

   user  system elapsed
  0.140   0.000   0.148

system.time(h3(dotot))

   user  system elapsed
  2.480   0.020   2.495

still better than without the optimization though, and pretty much as
readable.
moral of the story: it seems possible to write fast R, but the code won't
look pretty

thanks for the answers!








Storing the data column wise makes it easier to increase the size of the
array.

As a reference for the timing I use the following routine in 

Re: [R] function for inverse normal transformation

2012-07-20 Thread carol white
Thanks Rui.
I changed my scripts to the followings and I think that it still is not 
correct. See also the attached file.

Thanks for your help,


 tmp
 [1]  2.502519  1.828576  3.755778 17.415000  3.779296  2.956850  2.379663
 [8]  1.103559  8.920316  2.744500  2.938480  7.522174 10.629200  8.552259
[15]  5.425938  4.388906  0.00  0.723887 11.337860  3.763786


 tmp.p =2*pnorm(abs(scale(tmp)),lower.tail=FALSE)
  tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)
  tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)*sign(tmp)
 equal(tmp, tmp.qnorm)
[1] FALSE
 par(mfrow = c(1,3))
 hist(tmp)
 hist(tmp.p)
 hist(tmp.qnorm)




 From: Rui Barradas ruipbarra...@sapo.pt
To: carol white wht_...@yahoo.com 
Cc: r-help r-help@r-project.org 
Sent: Friday, July 20, 2012 2:02 PM
Subject: Re: [R] function for inverse normal transformation
 

Hello,

No it's not correct, you are computing a what seems to be a
two-tailed probabiity, so the inverse should account for it. Look
closely: you take the absolute value, then the upper tail
probability, then multiply 2 into it. Reverse these steps to get the
correct value.

# Helper function
equal - function(x, y, tol=.Machine$double.eps^0.5) all(abs(x -
y)   tol)

m - rnorm(5)
p - 2*pnorm(abs(m), lower.tail=FALSE)
m2 - qnorm(p/2, lower.tail=FALSE)*sign(m)

equal(m, m2)

(The helper function is just to test floating point values computed
differently for equality.)

Hope this helps,

Rui Barradas


Em 20-07-2012 12:36, carol white escreveu:

Thanks for your reply. So to derive it from a given data set, is the following 
correct to do? my_data.p =2*pnorm(abs(my_data),lower.tail=FALSE) my_data.q = 
qnorm(my_data.p) Cheers,  From: Duncan Murdoch 
murdoch.dun...@gmail.com Cc: r-h...@stat.math.ethz.ch 
r-h...@stat.math.ethz.ch Sent: Friday, July 20, 2012 1:23 PM
Subject: Re: [R] function for inverse normal transformation On 12-07-20 6:21 
AM, carol white wrote: 
Hi,
What is the function for inverse normal transformation? 
qnorm Duncan Murdoch 
Thanks, Carol     [[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. attachment: tmp.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.


[R] Speeding up a loop

2012-07-20 Thread wwreith
General problem: I have 20 projects that can be invested in and I need to
decide which combinations meet a certain set of standards. The total
possible combinations comes out to 2^20. However I know for a fact that the
number of projects must be greater than 5 and less than 13. So far the the
code below is the best I can come up with for iteratively creating a set to
check against my set of standards.

Code
x-matrix(0,nrow=1,ncol=20)
for(i in 1:2^20)
{
x[1]-x[1]+1
  for(j in 1:20)
  {
if(x[j]1)
{
  x[j]=0
  if(j20)
  {
x[j+1]=x[j+1]+1
  }
}
  }
if(sum(x)5  sum(x)13)
{
# insert criteria here.
}
}

my code forces me to create all 2^20 x's and then use an if statement to
decide if x is within my range of projects. Is there a faster way to
increment x. Any ideas on how to kill the for loop so that it won't attempt
to process an x where the sum is greater than 12 or less than 6?



--
View this message in context: 
http://r.789695.n4.nabble.com/Speeding-up-a-loop-tp4637201.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] Execute a function

2012-07-20 Thread Carla Moreira
Yes, I do.

But I need to control how the permutations are done.

Thank you.

2012/7/20 Jessica Streicher j.streic...@micromata.de

 You mean executing the function for all combinations of values?
 For example, if you have a-b-c-1:2
 you would get back the values of

 myfunc(1,1,1)
 myfunc(1,1,2)
 myfunc(1,2,1)
 myfunc(1,2,2)
 myfunc(2,1,1)
 myfunc(2,1,2)
 myfunc(2,2,1)
 myfunc(2,2,2)

 ?

 On 20.07.2012, at 13:05, carla moreira wrote:

 
  Hi,
 
  I would like to evaluate a function, with 3 arguments, for instance,
 
  myfunc-function(a,b,c) { sqrt(a)-exp(b)+4*c
 }
 
  How to execute  myfunc(x,y,z), for all x, all y and all z, where x,y,z
 are
  vectors?
 
  Thank you very much in advance
 
 
 
  --
  View this message in context:
 http://r.789695.n4.nabble.com/Execute-a-function-tp4637182.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.




-- 
Carla Moreira

[[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] qplot/ggplot

2012-07-20 Thread John Kane
I think we need some raw data.  Have a look at ?dput for a way to supply it.

It might help if you supplied some sample code of what you have tried.  

I want a line plot is not particularly helpful.

John Kane
Kingston ON Canada


 -Original Message-
 From: anamika...@gmail.com
 Sent: Fri, 20 Jul 2012 12:57:14 +0200
 To: r-help@r-project.org
 Subject: [R] qplot/ggplot
 
 Hello,
 I have a file like this (just a snapshot) where I have numerical
 values for various genes, I want a line plot with shading (may be
 using smooth ?) using qplot or ggplot :
 
 Gene1 10 14 12 23 11 11 33 1 ..(multiple columns)
 Gene2 4 2 1 1 3 4 1 2 .
 Gene3 2 5 7 5 6 89 7 3 ..
 Gene4 1 9 8 3 90 8 8  .
 
 Please help.
 
 Ana
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


FREE 3D MARINE AQUARIUM SCREENSAVER - Watch dolphins, sharks  orcas on your 
desktop!

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

2012-07-20 Thread Jessica Streicher
Well, what do you want to control there?

Need a subset? Need an ordering?


On 20.07.2012, at 15:00, Carla Moreira wrote:

 Yes, I do.
  
 But I need to control how the permutations are done.
  
 Thank you.
 
 2012/7/20 Jessica Streicher j.streic...@micromata.de
 You mean executing the function for all combinations of values?
 For example, if you have a-b-c-1:2
 you would get back the values of
 
 myfunc(1,1,1)
 myfunc(1,1,2)
 myfunc(1,2,1)
 myfunc(1,2,2)
 myfunc(2,1,1)
 myfunc(2,1,2)
 myfunc(2,2,1)
 myfunc(2,2,2)
 
 ?
 
 On 20.07.2012, at 13:05, carla moreira wrote:
 
 
  Hi,
 
  I would like to evaluate a function, with 3 arguments, for instance,
 
  myfunc-function(a,b,c) { sqrt(a)-exp(b)+4*c
 }
 
  How to execute  myfunc(x,y,z), for all x, all y and all z, where x,y,z are
  vectors?
 
  Thank you very much in advance
 
 
 
  --
  View this message in context: 
  http://r.789695.n4.nabble.com/Execute-a-function-tp4637182.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.
 
 
 
 
 -- 
 Carla Moreira


[[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] Crosstab with Average and Count

2012-07-20 Thread John Kane
Not quite what you asked for but would this do?

library(plyr)
x - as.factor(c(1,1,1,2,2,2,3,3,3))
y - as.factor(c(10,10,10,20,20,20,30,30,30))
z - c(100,100,NA,200,200,200,300,300,300)
df1  -  data.frame(x, y, z)

(dsum - ddply(df1, .(x, y), summarise, sum = sum(z, na.rm = TRUE), 
mean = mean(z, na.rm = TRUE), count = length(z)))




John Kane
Kingston ON Canada


 -Original Message-
 From: viora...@gmail.com
 Sent: Fri, 20 Jul 2012 03:30:16 -0700 (PDT)
 To: r-help@r-project.org
 Subject: [R] Crosstab with Average and Count
 
 I have the following data:
 
 x - as.factor(c(1,1,1,2,2,2,3,3,3))
 y - as.factor(c(10,10,10,20,20,20,30,30,30))
 z - c(100,100,NA,200,200,200,300,300,300)
 
 I could create the cross tab of x and y with Sum of z as its elements
 using
 the xtabs function as follows:
 
 # X Vs. Y with Sum Z
 
 xtabs(z ~ x + y)
 
y
 x10  20  30
   1 200   0   0
   2   0 600   0
   3   0   0 900
 
 How do I replace the sum with average and count so that I can get the
 following outputs??
 
 # X Vs. Y with Average of Z
y
 x  10  20  30
   1100 0   0
   20   200 0
   30   0   300
 
 # X Vs. Y with Count Z
   y
 x10  20  30
  12   0   0
  20   3   0
  30   0   3
 
 Would appreciate any help on these? Thank you.
 
 Ravi
 
 
 
 
 
 --
 View this message in context:
 http://r.789695.n4.nabble.com/Crosstab-with-Average-and-Count-tp4637180.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.


GET FREE SMILEYS FOR YOUR IM  EMAIL - Learn more at 
http://www.inbox.com/smileys
Works with AIM®, MSN® Messenger, Yahoo!® Messenger, ICQ®, Google Talk™ and most 
webmails

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] overlay contourplot with labels over spplot with color regions?

2012-07-20 Thread Martin Ivanov
 Dear R users,

I have the following problem. I plot a SpatialPixelsDataFrame object with 
spplot and I need to add a contourplot
of another spatial variable. To be more precise, I plot the  
SpatialPixelsDataFrame of the mean precipitation over Germany with coloured 
regions and I want to overlay the isolines of the 500, 1000 and 1500 m height 
above sea level. Of course I can plot the lines as SpatialLines using the 
sp.layout argument to spplot, but so I lose the labels on the lines. And I want 
the height above sea level to be labelled just as it is labelled in a contour 
plot. Is it possible to achieve that with lattice? 

I am still new to lattice, so I beg Your excuse in case I am asking something 
trivial. I googled a lot and tried a lot of things, but none worked. I tried 
with contour, it does work with Spatial objects, but unfortunately contour does 
not seem work with lattice. And the lattice equivalent contourplot does not 
recognise Spatial objects. 

I am already out of ideas. I am looking forward to Your replies. Any 
suggestions will be appreciated.

Best regards,

Martin

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

2012-07-20 Thread Marc Schwartz
On Jul 20, 2012, at 5:30 AM, vioravis wrote:

 I have the following data:
 
 x - as.factor(c(1,1,1,2,2,2,3,3,3))
 y - as.factor(c(10,10,10,20,20,20,30,30,30))
 z - c(100,100,NA,200,200,200,300,300,300)
 
 I could create the cross tab of x and y with Sum of z as its elements using
 the xtabs function as follows:
 
 # X Vs. Y with Sum Z
 
 xtabs(z ~ x + y)
 
   y
 x10  20  30
  1 200   0   0
  2   0 600   0
  3   0   0 900
 
 How do I replace the sum with average and count so that I can get the
 following outputs??
 
 # X Vs. Y with Average of Z
   y
 x  10  20  30
  1100 0   0
  20   200 0
  30   0   300
 
 # X Vs. Y with Count Z
  y
 x10  20  30
 12   0   0
 20   3   0
 30   0   3
 
 Would appreciate any help on these? Thank you.
 
 Ravi


You can use ?tapply, albeit you will get NA's rather than 0's:

 tapply(z,  list(x, y), mean, na.rm = TRUE)
   10  20  30
1 100  NA  NA
2  NA 200  NA
3  NA  NA 300


 tapply(z,  list(x, y), function(x) sum(!is.na(x)))
  10 20 30
1  2 NA NA
2 NA  3 NA
3 NA NA  3


Regards,

Marc Schwartz

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


Re: [R] Crosstab with Average and Count

2012-07-20 Thread David L Carlson
There are lots of possibilities. Here's one using only xtabs():

dframe - na.omit(data.frame(x, y, z))
zsum - xtabs(z~x+y, dframe)
zcount - xtabs(~x+y, dframe)
zmean - ifelse(is.nan(zsum/zcount), 0, zsum/zcount)

--
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-bounces@r-
 project.org] On Behalf Of Marc Schwartz
 Sent: Friday, July 20, 2012 8:41 AM
 To: vioravis
 Cc: r-help@r-project.org
 Subject: Re: [R] Crosstab with Average and Count
 
 On Jul 20, 2012, at 5:30 AM, vioravis wrote:
 
  I have the following data:
 
  x - as.factor(c(1,1,1,2,2,2,3,3,3))
  y - as.factor(c(10,10,10,20,20,20,30,30,30))
  z - c(100,100,NA,200,200,200,300,300,300)
 
  I could create the cross tab of x and y with Sum of z as its elements
 using
  the xtabs function as follows:
 
  # X Vs. Y with Sum Z
 
  xtabs(z ~ x + y)
 
y
  x10  20  30
   1 200   0   0
   2   0 600   0
   3   0   0 900
 
  How do I replace the sum with average and count so that I can get the
  following outputs??
 
  # X Vs. Y with Average of Z
y
  x  10  20  30
   1100 0   0
   20   200 0
   30   0   300
 
  # X Vs. Y with Count Z
   y
  x10  20  30
  12   0   0
  20   3   0
  30   0   3
 
  Would appreciate any help on these? Thank you.
 
  Ravi
 
 
 You can use ?tapply, albeit you will get NA's rather than 0's:
 
  tapply(z,  list(x, y), mean, na.rm = TRUE)
10  20  30
 1 100  NA  NA
 2  NA 200  NA
 3  NA  NA 300
 
 
  tapply(z,  list(x, y), function(x) sum(!is.na(x)))
   10 20 30
 1  2 NA NA
 2 NA  3 NA
 3 NA NA  3
 
 
 Regards,
 
 Marc Schwartz
 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide http://www.R-project.org/posting-
 guide.html
 and provide commented, minimal, self-contained, reproducible code.

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


Re: [R] Creating Multiple Repeating samples and Cross Correlating them.

2012-07-20 Thread John Kane
It would help if you supplied some code for us to see what you have tried. 
However something like this would presumably give you the sample data sets 
--note only n=100 in the example.


mymat  -  matrix(NA, ncol=  100, nrow= 100)
for (n in 1: 100) {
  mymat[,n] -  sample(c(0,1), 100, replace = TRUE)
}
mymat


John Kane
Kingston ON Canada


 -Original Message-
 From: baile...@ohsu.edu
 Sent: Thu, 19 Jul 2012 14:17:43 -0700 (PDT)
 To: r-help@r-project.org
 Subject: [R] Creating Multiple Repeating samples and Cross Correlating
 them.
 
 I'm having a lot of trouble getting this done and nothing I've written
 has
 been remotely successful. Basically I have about 64 binary data sets
 stored
 as vectors, with 72900 entries in each. I am not very familiar with cross
 correlations, but I was advised that I should create 1000 randomized date
 sets (used sample() for that) to use as a control, then compare the cross
 correlation of my real data sets to the controls. First of all I am
 having
 trouble creating 1000 different sample()'s of the same vector, and then I
 am
 not all too sure how to cross correlate them Nothing seems to produce
 what I
 want, all I get are copies of the same sample() of the vector repeated
 over
 and over. Please Help?
 
 --
 View this message in context:
 http://r.789695.n4.nabble.com/Creating-Multiple-Repeating-samples-and-Cross-Correlating-them-tp4637131.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.


Send any screenshot to your friends in seconds...
Works in all emails, instant messengers, blogs, forums and social networks.
TRY IM TOOLPACK at http://www.imtoolpack.com/default.aspx?rc=if2 for FREE

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

2012-07-20 Thread Peter Ehlers

On 2012-07-20 04:05, carla moreira wrote:


Hi,

I would like to evaluate a function, with 3 arguments, for instance,

myfunc-function(a,b,c) { sqrt(a)-exp(b)+4*c
 }

How to execute  myfunc(x,y,z), for all x, all y and all z, where x,y,z are
vectors?


Is this what you have in mind:

  myfunc - function(a, b, c){ sqrt(a)-exp(b)+4*c }
  myfunc2 - function(x){
a - x[1]
b - x[2]
c - x[3]
myfunc(a, b, c)
  }

  x - c(1, 4, 9)
  y - 1:2
  z - c(10, -10, 2, 20)
  d - expand.grid(x, y, z)
  d$value - apply(d, 1, myfunc2)

?

Peter Ehlers



Thank you very much in advance



--
View this message in context: 
http://r.789695.n4.nabble.com/Execute-a-function-tp4637182.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] Speeding up a loop

2012-07-20 Thread Jean V Adams
I've had to do something similar, so I wrote a small function to help.
This runs in about 1/4 the time of your code on my machine.
Others may have a more efficient approach.

all.combs - function(num, from=0, to=num) {
# create a matrix of all possible combinations of num items
# restrict each row to have between from and to items 
res - vector(list, to-from+1)
for(i in seq(from:to)) {
j - (from:to)[i]
if(j==0) res[[i]] - rep(FALSE, num)
comb - combn(num, j)
res[[i]] - t(apply(comb, 2, function(x) 
!is.na(match(1:num, x
}
do.call(rbind, res)
}

all.combs(20, 5, 13)

Jean


wwreith reith_will...@bah.com wrote on 07/20/2012 07:45:30 AM:

 General problem: I have 20 projects that can be invested in and I need 
to
 decide which combinations meet a certain set of standards. The total
 possible combinations comes out to 2^20. However I know for a fact that 
the
 number of projects must be greater than 5 and less than 13. So far the 
the
 code below is the best I can come up with for iteratively creating a set 
to
 check against my set of standards.
 
 Code
 x-matrix(0,nrow=1,ncol=20)
 for(i in 1:2^20)
 {
 x[1]-x[1]+1
   for(j in 1:20)
   {
 if(x[j]1)
 {
   x[j]=0
   if(j20)
   {
 x[j+1]=x[j+1]+1
   }
 }
   }
 if(sum(x)5  sum(x)13)
 {
 # insert criteria here.
 }
 }
 
 my code forces me to create all 2^20 x's and then use an if statement to
 decide if x is within my range of projects. Is there a faster way to
 increment x. Any ideas on how to kill the for loop so that it won't 
attempt
 to process an x where the sum is greater than 12 or less than 6?

[[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] Execute a function

2012-07-20 Thread Bert Gunter
Inline.

-- Bert

On Fri, Jul 20, 2012 at 6:59 AM, Peter Ehlers ehl...@ucalgary.ca wrote:
 On 2012-07-20 04:05, carla moreira wrote:


 Hi,

 I would like to evaluate a function, with 3 arguments, for instance,

 myfunc-function(a,b,c) { sqrt(a)-exp(b)+4*c
  }

 How to execute  myfunc(x,y,z), for all x, all y and all z, where x,y,z are
 vectors?


 Is this what you have in mind:

   myfunc - function(a, b, c){ sqrt(a)-exp(b)+4*c }
   myfunc2 - function(x){
 a - x[1]
 b - x[2]
 c - x[3]
 myfunc(a, b, c)
   }

   x - c(1, 4, 9)
   y - 1:2
   z - c(10, -10, 2, 20)
   d - expand.grid(x, y, z)

Peter, what's wrong with
with(d,myfunc(x,y,z))??

I realize this depends on the function be vectorizable, but isn't that
the point? It could be orders of magnitude faster than looping via
apply.

-- Bert

   d$value - apply(d, 1, myfunc2)

 ?

 Peter Ehlers


 Thank you very much in advance



 --
 View this message in context:
 http://r.789695.n4.nabble.com/Execute-a-function-tp4637182.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.



-- 

Bert Gunter
Genentech Nonclinical Biostatistics

Internal Contact Info:
Phone: 467-7374
Website:
http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm

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


Re: [R] function for inverse normal transformation

2012-07-20 Thread David L Carlson
Rui's example included z-score data (drawn from rnorm). You converted your
data to z-scores so you need to compare your results to the z-scores not
the original data.

Change these lines:

tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)*sign(scale(tmp))
# sign is of scale(tmp) not tmp
equal(scale(tmp), tmp.qnorm)
# compare to scale(tmp) not tmp

--
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-bounces@r-
 project.org] On Behalf Of carol white
 Sent: Friday, July 20, 2012 7:29 AM
 To: Rui Barradas
 Cc: r-help@r-project.org
 Subject: Re: [R] function for inverse normal transformation
 
 Thanks Rui.
 I changed my scripts to the followings and I think that it still is not
 correct. See also the attached file.
 
 Thanks for your help,
 
 
  tmp
  [1]  2.502519  1.828576  3.755778 17.415000  3.779296  2.956850
 2.379663
  [8]  1.103559  8.920316  2.744500  2.938480  7.522174 10.629200
 8.552259
 [15]  5.425938  4.388906  0.00  0.723887 11.337860  3.763786
 
 
  tmp.p =2*pnorm(abs(scale(tmp)),lower.tail=FALSE)
   tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)
   tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)*sign(tmp)
  equal(tmp, tmp.qnorm)
 [1] FALSE
  par(mfrow = c(1,3))
  hist(tmp)
  hist(tmp.p)
  hist(tmp.qnorm)
 
 
 
 
  From: Rui Barradas ruipbarra...@sapo.pt
 To: carol white wht_...@yahoo.com
 Cc: r-help r-help@r-project.org
 Sent: Friday, July 20, 2012 2:02 PM
 Subject: Re: [R] function for inverse normal transformation
 
 
 Hello,
 
 No it's not correct, you are computing a what seems to be a
 two-tailed probabiity, so the inverse should account for it. Look
 closely: you take the absolute value, then the upper tail
 probability, then multiply 2 into it. Reverse these steps to get
 the
 correct value.
 
 # Helper function
 equal - function(x, y, tol=.Machine$double.eps^0.5) all(abs(x -
 y)   tol)
 
 m - rnorm(5)
 p - 2*pnorm(abs(m), lower.tail=FALSE)
 m2 - qnorm(p/2, lower.tail=FALSE)*sign(m)
 
 equal(m, m2)
 
 (The helper function is just to test floating point values computed
 differently for equality.)
 
 Hope this helps,
 
 Rui Barradas
 
 
 Em 20-07-2012 12:36, carol white escreveu:
 
 Thanks for your reply. So to derive it from a given data set, is the
 following correct to do? my_data.p
 =2*pnorm(abs(my_data),lower.tail=FALSE) my_data.q = qnorm(my_data.p)
 Cheers,  From: Duncan Murdoch
 murdoch.dun...@gmail.com Cc: r-h...@stat.math.ethz.ch r-
 h...@stat.math.ethz.ch Sent: Friday, July 20, 2012 1:23 PM
 Subject: Re: [R] function for inverse normal transformation On 12-07-20
 6:21 AM, carol white wrote:
 Hi,
 What is the function for inverse normal transformation?
 qnorm Duncan Murdoch
 Thanks, Carol     [[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 for inverse normal transformation

2012-07-20 Thread Duncan Murdoch

On 20/07/2012 8:28 AM, carol white wrote:

Thanks Rui.
I changed my scripts to the followings and I think that it still is not correct.


You haven't told us clearly what you are trying to achieve, and you 
don't tell us what is wrong with what you have.  How do you expect 
anyone to help?


Duncan Murdoch



See also the attached file.

Thanks for your help,


  tmp
  [1]  2.502519  1.828576  3.755778 17.415000  3.779296  2.956850  2.379663
  [8]  1.103559  8.920316  2.744500  2.938480  7.522174 10.629200  8.552259
[15]  5.425938  4.388906  0.00  0.723887 11.337860  3.763786


  tmp.p =2*pnorm(abs(scale(tmp)),lower.tail=FALSE)
  tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)
  tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)*sign(tmp)
 equal(tmp, tmp.qnorm)
[1] FALSE
 par(mfrow = c(1,3))
 hist(tmp)
 hist(tmp.p)
 hist(tmp.qnorm)




  From: Rui Barradas ruipbarra...@sapo.pt
To: carol white wht_...@yahoo.com
Cc: r-help r-help@r-project.org
Sent: Friday, July 20, 2012 2:02 PM
Subject: Re: [R] function for inverse normal transformation
  


Hello,

No it's not correct, you are computing a what seems to be a
 two-tailed probabiity, so the inverse should account for it. Look
 closely: you take the absolute value, then the upper tail
 probability, then multiply 2 into it. Reverse these steps to get the
 correct value.

# Helper function
equal - function(x, y, tol=.Machine$double.eps^0.5) all(abs(x -
 y)   tol)

m - rnorm(5)
p - 2*pnorm(abs(m), lower.tail=FALSE)
m2 - qnorm(p/2, lower.tail=FALSE)*sign(m)

equal(m, m2)

(The helper function is just to test floating point values computed
 differently for equality.)

Hope this helps,

Rui Barradas


Em 20-07-2012 12:36, carol white escreveu:

Thanks for your reply. So to derive it from a given data set, is the following correct to do? 
my_data.p =2*pnorm(abs(my_data),lower.tail=FALSE) my_data.q = qnorm(my_data.p) Cheers, 
 From: Duncan Murdoch murdoch.dun...@gmail.com Cc: 
r-h...@stat.math.ethz.ch r-h...@stat.math.ethz.ch Sent: Friday, July 20, 2012 
1:23 PM
Subject: Re: [R] function for inverse normal transformation On 12-07-20 6:21 
AM, carol white wrote:
Hi,
What is the function for inverse normal transformation?
qnorm Duncan Murdoch
Thanks, Carol [[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.


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

2012-07-20 Thread Jean V Adams
Michael,

Use dput() to output your data (or perhaps a small subset).  Then paste 
the result of that call and your R code (just those lines of code that are 
needed to reproduce the problem) right in your message to R-help.  That 
makes it easier for the R-help list readers to help you troubleshoot.

Jean


Michael Eisenring eimic...@ethz.ch wrote on 07/20/2012 05:15:32 AM:

 Dear r-help members.
 I would like to compare species numbers of moths between eight 
 different forests (each sampled for six nights). I would like to do 
 a nested anova to compare species numbers between forests and nights.
 For more site specific details I wanted to do a Tukey test 
 (TukeyHSD). Unfortunately this test is not working (message:  nicht
 anwendbare Methode für 'TukeyHSD' auf Objekt der Klasse c('anova', 
 'data.frame ) /(message: method not appropriate for TukeyHSD for 
 objects of the class anova, data.frame).
 I tried it also with the aov instead of the anova function but 
 neither the aov nor the anova approach is working.
 Could anyone help me with this  problem?
 Attached you can find the txt.file and the r-script.
 I'm working with a MAC OSX snow leopard, with R-Studio 0.96 316.
 
 Thank you very much.
 Michael 

[[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] Speeding up a loop

2012-07-20 Thread Petr Savicky
On Fri, Jul 20, 2012 at 05:45:30AM -0700, wwreith wrote:
 General problem: I have 20 projects that can be invested in and I need to
 decide which combinations meet a certain set of standards. The total
 possible combinations comes out to 2^20. However I know for a fact that the
 number of projects must be greater than 5 and less than 13. So far the the
 code below is the best I can come up with for iteratively creating a set to
 check against my set of standards.
 
 Code
 x-matrix(0,nrow=1,ncol=20)
 for(i in 1:2^20)
 {
 x[1]-x[1]+1
   for(j in 1:20)
   {
 if(x[j]1)
 {
   x[j]=0
   if(j20)
   {
 x[j+1]=x[j+1]+1
   }
 }
   }
 if(sum(x)5  sum(x)13)
 {
 # insert criteria here.
 }
 }
 
 my code forces me to create all 2^20 x's and then use an if statement to
 decide if x is within my range of projects. Is there a faster way to
 increment x. Any ideas on how to kill the for loop so that it won't attempt
 to process an x where the sum is greater than 12 or less than 6?

Hi.

The restriction on the sum of the rows between 6 and 12 eliminates the
tails of the distribution, not the main part. So, the final number of
rows is not much smaller than 2^20. More exactly, it is

  sum(choose(20, 6:12))

which is about 0.8477173 * 2^20. On the other hand, all combinations
may be created using expand.grid() faster than using a for loop.

Try the following

  g - as.matrix(expand.grid(rep(list(0:1), times=20)))
  s - rowSums(g)
  x - g[s  5  s  13, ]
  nrow(x)

  [1] 96

Hope this helps.

Petr Savicky.

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


Re: [R] Creating a pdf with layers?

2012-07-20 Thread Christoph Scherber
Dear Prof Ripley,

Many thanks for your response. In fact, latticeExtra and ggplot2 both state 
they work with layers,
but essentially this is nothing more than sequentially adding graphical output 
to an existing device.

But given that PDF can include layers since version 1.5 it would be interesting 
to be able to divert
output to different layers inside a PDF structure.

Best regards
Christoph Scherber





Am 20.07.2012 13:48, schrieb Prof Brian Ripley:
 On 20/07/12 12:07, Christoph Scherber wrote:
 Dear all,

 Is it possible to create a pdf file with layers using the pdf() device in R?
 
 No.  Is it possible to specify layers in the R graphics language or any 
 device?  (From what I understand by 'layers', no.)
 
 The author of pdf()
 

 Many thanks for your help!
 Christoph


 (using R 2.15.1 on Windows 7 64-Bit)
 
 
 


-- 
PD Dr Christoph Scherber
Georg-August University Goettingen
Department of Crop Science
Agroecology
Grisebachstrasse 6
D-37077 Goettingen
Germany
phone 0049 (0)551 39 8807
fax 0049 (0)551 39 8806
http://www.gwdg.de/~cscherb1
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] convert date to a factor

2012-07-20 Thread Yolande Tra
Hello,

I would like to convert date as a factor to represent time in a repeated
measure situation in the following code. How would I do that?
  d - read.csv(file.path(dataDir,data.csv), as.is=T,stringsAsFactors =
FALSE)
 d[1:2,]
id date ab  c y
1   1 8/6/2008 Red15 B  22
2   1 8/6/2008   Green   15 B  22
Thank you,
Y

[[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] Speeding up a loop

2012-07-20 Thread Jean V Adams
Petr,

This is great.
MUCH faster than the code I provided.
And much more elegant code.
Thanks for posting!

Jean


Petr Savicky savi...@cs.cas.cz wrote on 07/20/2012 09:26:34 AM:

 On Fri, Jul 20, 2012 at 05:45:30AM -0700, wwreith wrote:
  General problem: I have 20 projects that can be invested in and I need 
to
  decide which combinations meet a certain set of standards. The total
  possible combinations comes out to 2^20. However I know for a fact 
that the
  number of projects must be greater than 5 and less than 13. So far the 
the
  code below is the best I can come up with for iteratively creating a 
set to
  check against my set of standards.
  
  Code
  x-matrix(0,nrow=1,ncol=20)
  for(i in 1:2^20)
  {
  x[1]-x[1]+1
for(j in 1:20)
{
  if(x[j]1)
  {
x[j]=0
if(j20)
{
  x[j+1]=x[j+1]+1
}
  }
}
  if(sum(x)5  sum(x)13)
  {
  # insert criteria here.
  }
  }
  
  my code forces me to create all 2^20 x's and then use an if statement 
to
  decide if x is within my range of projects. Is there a faster way to
  increment x. Any ideas on how to kill the for loop so that it won't 
attempt
  to process an x where the sum is greater than 12 or less than 6?
 
 Hi.
 
 The restriction on the sum of the rows between 6 and 12 eliminates the
 tails of the distribution, not the main part. So, the final number of
 rows is not much smaller than 2^20. More exactly, it is
 
   sum(choose(20, 6:12))
 
 which is about 0.8477173 * 2^20. On the other hand, all combinations
 may be created using expand.grid() faster than using a for loop.
 
 Try the following
 
   g - as.matrix(expand.grid(rep(list(0:1), times=20)))
   s - rowSums(g)
   x - g[s  5  s  13, ]
   nrow(x)
 
   [1] 96
 
 Hope this helps.
 
 Petr Savicky.

[[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] Speeding up a loop

2012-07-20 Thread Petr Savicky
On Fri, Jul 20, 2012 at 04:26:34PM +0200, Petr Savicky wrote:
 On Fri, Jul 20, 2012 at 05:45:30AM -0700, wwreith wrote:
  General problem: I have 20 projects that can be invested in and I need to
  decide which combinations meet a certain set of standards. The total
  possible combinations comes out to 2^20. However I know for a fact that the
  number of projects must be greater than 5 and less than 13. So far the the
  code below is the best I can come up with for iteratively creating a set to
  check against my set of standards.
  
  Code
  x-matrix(0,nrow=1,ncol=20)
  for(i in 1:2^20)
  {
  x[1]-x[1]+1
for(j in 1:20)
{
  if(x[j]1)
  {
x[j]=0
if(j20)
{
  x[j+1]=x[j+1]+1
}
  }
}
  if(sum(x)5  sum(x)13)
  {
  # insert criteria here.
  }
  }
  
  my code forces me to create all 2^20 x's and then use an if statement to
  decide if x is within my range of projects. Is there a faster way to
  increment x. Any ideas on how to kill the for loop so that it won't attempt
  to process an x where the sum is greater than 12 or less than 6?
 
 Hi.
 
 The restriction on the sum of the rows between 6 and 12 eliminates the
 tails of the distribution, not the main part. So, the final number of
 rows is not much smaller than 2^20. More exactly, it is
 
   sum(choose(20, 6:12))
 
 which is about 0.8477173 * 2^20. On the other hand, all combinations
 may be created using expand.grid() faster than using a for loop.
 
 Try the following
 
   g - as.matrix(expand.grid(rep(list(0:1), times=20)))
   s - rowSums(g)
   x - g[s  5  s  13, ]

Hi.

The above code creates a matrix, whose rows are vectors of 0,1, which
contain between 6 and 12 ones. Using this matrix, it is possible to
go through all these combinations using a for loop as follows.

  for (i in seq.int(length=nrow(x))) {
  here, x[i, ] is a row of the matrix
  }

Another option is to use ifelse() function, which allows to evaluate
a condition on the whole columns of the matrix. If this is possible,
then it is more efficient than a for loop.

Instead of using expand.grid() to create all 2^20 combinations, it is
possible to create only rows with a specified number of ones. The
rows of length n with exactly k ones can be created as follows.

  n - 5
  k - 2
  ind - combn(n, k)
  m - ncol(ind)
  x - matrix(0, nrow=m, ncol=n)
  x[cbind(rep(1:m, each=k), c(ind))] - 1
  x

   [1,]11000
   [2,]10100
   [3,]10010
   [4,]10001
   [5,]01100
   [6,]01010
   [7,]01001
   [8,]00110
   [9,]00101
  [10,]00011

Hope this helps.

Petr Savicky.

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

2012-07-20 Thread Peter Ehlers

Bert,

The only thing wrong is that I'm still 75% asleep! Yikes!!
Thanks for the heads-up.

Carla: See Bert's solution.

Peter Ehlers

On 2012-07-20 07:10, Bert Gunter wrote:

Inline.

-- Bert

On Fri, Jul 20, 2012 at 6:59 AM, Peter Ehlers ehl...@ucalgary.ca wrote:

On 2012-07-20 04:05, carla moreira wrote:



Hi,

I would like to evaluate a function, with 3 arguments, for instance,

myfunc-function(a,b,c) { sqrt(a)-exp(b)+4*c
  }

How to execute  myfunc(x,y,z), for all x, all y and all z, where x,y,z are
vectors?



Is this what you have in mind:

   myfunc - function(a, b, c){ sqrt(a)-exp(b)+4*c }
   myfunc2 - function(x){
 a - x[1]
 b - x[2]
 c - x[3]
 myfunc(a, b, c)
   }

   x - c(1, 4, 9)
   y - 1:2
   z - c(10, -10, 2, 20)
   d - expand.grid(x, y, z)


Peter, what's wrong with
with(d,myfunc(x,y,z))??

I realize this depends on the function be vectorizable, but isn't that
the point? It could be orders of magnitude faster than looping via
apply.

-- Bert


   d$value - apply(d, 1, myfunc2)

?

Peter Ehlers



Thank you very much in advance



--
View this message in context:
http://r.789695.n4.nabble.com/Execute-a-function-tp4637182.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] [RsR] How does rlm in R decide its w weights for each IRLS iteration?

2012-07-20 Thread S Ellison
 

 -Original Message-
 Subject: [RsR] How does rlm in R decide its w weights for 
 each IRLS iteration?
 I am also confused about the manual:
 
a.  The input arguments:
 
 wt.method are the weights case weights (giving the relative 
 importance of case, so a weight of 2 means there are two of 
 these) or the inverse of the variances, so a weight of two 
 means this error is half as variable?

When you give rlm weights (called 'weights', not 'w' on input, though you can 
abbreviate to 'w'), you need to tell it which of these two possibilities you 
used. 
If you gave it case numbers, say wt.method=case; if you gave it inverse 
variance weights, say wt.method=inv.var.
The default is inv.var.


 The input argument w is used for the initial values of the 
 rlm IRLS weighting and the output value w is the converged w.
There is no input argument 'w' for rlm (see above). 
The output w are a  calculated using the psi function, so between 0 and 1.
The effective weights for the final estimate would then be something like 
w*weights, using the full name of the input argument (and if I haven't 
forgotten a square root somewhere). At least, that would work for a simple 
location estimate (eg rlm(x~1)).

 If my understanding above is correct, how does rlm decide 
 its w for each IRLS iteration then?
It uses the given psi functions to calculate the iterative weights based on the 
scaled residuals.

 Any pointers/tutorials/notes to the calculation of these 
 w's in each IRLS iteration?
Read the cited references for a detailed guide. Or, of course, MASS - the 
package is, after all, intended to support the book, not replace it.



S Ellison

***
This email and any attachments are confidential. Any use...{{dropped:8}}

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


Re: [R] why does this simple example NOT work?

2012-07-20 Thread Berend Hasselman

Duncan Murdoch-2 wrote
 
 On 12-07-20 2:50 AM, Martin Ivanov wrote:
   Hello,

 I want to create and save objects in a loop, but this is precluded by the
 following obstacle:
 this part of the script fails to work:

 assign(x=paste(a, 1, sep=), value=1);
 save(paste(a, 1, sep=), file=paste(paste(a, 1, sep=), .RData,
 sep=))
 
 paste(a, 1, sep=) is equivalent to a1.  So you are saving the 
 string a1 to the file a1.RData.  You presumably want to save the 
 value of the variable a1, not that string.  (Though you didn't say that. 
   It's helpful to describe what you don't like, don't just say it fails 
 to work.)
 
 If you have the name of a variable in a string, you can use get() to get 
 the value.  So I think this code would do what you want:
 
 name - paste0(a, 1)
 assign(x=name, value=1)
 save(get(name), file=paste0(name, .RData))
 

I tried that in R2.15.1 and got the error message:

Error in save(get(name), file = paste0(name, .RData)) : 
  object ‘get(name)’ not found

Using the list argument worked:

save(list=name, file=paste0(name, .RData))

Berend




--
View this message in context: 
http://r.789695.n4.nabble.com/why-does-this-simple-example-NOT-work-tp4637158p4637228.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] TukeyHSD not working

2012-07-20 Thread Jean V Adams
Michael Eisenring eimic...@ethz.ch wrote on 07/20/2012 09:35:03 AM:

 Dear Jean,
 thanks for this email. I'm grateful for this instruction Just to 
 make sure that I understood you properly: something like this?:


Michael,
Yes, this is perfect.  Very helpful.


 My data:
 
 dput(geo_anova_nested_input):
 
 structure(list(Site = c(1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 
 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 7L, 7L, 
 8L, 8L, 8L, 8L, 8L, 8L), Night = 1:48, Species = c(21L, 22L, 
 13L, 12L, 26L, 12L, 28L, 22L, 25L, 27L, 27L, 10L, 21L, 31L, 33L, 
 36L, 32L, 17L, 31L, 26L, 27L, 38L, 32L, 14L, 27L, 36L, 37L, 40L, 
 38L, 5L, 20L, 31L, 29L, 36L, 30L, 30L, 18L, 19L, 16L, 43L, 31L, 
 27L, 24L, 27L, 28L, 49L, 34L, 30L)), .Names = c(Site, Night, 
 Species), class = data.frame, row.names = c(NA, -48L))
 
 RCode
 
 attach(geo_anova_nested_input)
 geo_anova_nested_input
 data-anova(lm(Species~Site/Night))
 TukeyHSD(data,Site, conf.level=0.90)


There are a number of issues with the line of code data-anova...
First of all, look at the output from
lm(Species~Site/Night)
You are fitting a linear regression model with Site and Site x Night 
interaction as the two independent continuous variables.  You have no 
categorical variables.  I assume that you want Site to be a categorical 
variable (it represents forests, right?).  I'm not sure how you want to 
deal with Night.  You may want to talk over your design with a 
statistician to be sure the model you are using is appropriate for the way 
your data were collected.

Second, look at the output from
anova(lm(Species~Site/Night))
This is simply an ANOVA table for the regression model you just fit.  But 
TukeyHSD(), in your next line of code, wants a ... fitted model object, 
usually an aov fit.  This is clearly stated in the help page for 
TukeyHSD().
?TukeyHSD

So, here's an example of some code that works (gives output without 
warning/error messages) using your data.  I have ignored Night, since I'm 
not sure how it should be dealt with (again, talk to a statistician). But, 
hopefully this will help you to understand how to make things work in R.
# create a new variable for Site, as a categorical variable (aka 
factor)
geo_anova_nested_input$Sitef - 
as.factor(geo_anova_nested_input$Site)
# fit an ANOVA model to the data
fit - aov(Species ~ Sitef, data=geo_anova_nested_input)
# compute Tukey honest significant differences for the Site factor
TukeyHSD(fit, Sitef, conf.level=0.90)

Jean


 Michael, 
 
 Use dput() to output your data (or perhaps a small subset).  Then 
 paste the result of that call and your R code (just those lines of 
 code that are needed to reproduce the problem) right in your message
 to R-help.  That makes it easier for the R-help list readers to help
 you troubleshoot. 
 
 Jean 
 
 
 Michael Eisenring eimic...@ethz.ch wrote on 07/20/2012 05:15:32 AM:
 
  Dear r-help members.
  I would like to compare species numbers of moths between eight 
  different forests (each sampled for six nights). I would like to do 
  a nested anova to compare species numbers between forests and nights.
  For more site specific details I wanted to do a Tukey test 
  (TukeyHSD). Unfortunately this test is not working (message:  nicht
  anwendbare Methode für 'TukeyHSD' auf Objekt der Klasse c('anova', 
  'data.frame ) /(message: method not appropriate for TukeyHSD for 
  objects of the class anova, data.frame).
  I tried it also with the aov instead of the anova function but 
  neither the aov nor the anova approach is working.
  Could anyone help me with this  problem?
  Attached you can find the txt.file and the r-script.
  I'm working with a MAC OSX snow leopard, with R-Studio 0.96 316.
  
  Thank you very much.
  Michael 
[[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] [External] Re: Speeding up a loop

2012-07-20 Thread Jean V Adams
Reith, William [USA] reith_will...@bah.com wrote on 07/20/2012 
09:52:02 AM:

 Would this matrix eat up memory making the rest of my program 
 slower? Each x needs to be multiplied by a matrix and the results 
 checked against a set of thresholds. Doing them one at a time takes 
 at least 24 hours right now. 
 
 Optimizing a program is not my thing. 
 
 Sent from my Verizon Wireless 4GLTE smartphone


It's not mine either.
Better to ask the group.
I'm ccing R-help on this message.

Jean


 - Reply message -
 From: Jean V Adams jvad...@usgs.gov
 To: Reith, William [USA] reith_will...@bah.com
 Cc: r-help@r-project.org r-help@r-project.org
 Subject: [External] Re: [R] Speeding up a loop
 Date: Fri, Jul 20, 2012 10:05 am

 
 
 I've had to do something similar, so I wrote a small function to help. 
 This runs in about 1/4 the time of your code on my machine. 
 Others may have a more efficient approach. 
 
 all.combs - function(num, from=0, to=num) { 
 # create a matrix of all possible combinations of num items 
 # restrict each row to have between from and to items 
 res - vector(list, to-from+1) 
 for(i in seq(from:to)) { 
 j - (from:to)[i] 
 if(j==0) res[[i]] - rep(FALSE, num) 
 comb - combn(num, j) 
 res[[i]] - t(apply(comb, 2, function(x) !is.na
 (match(1:num, x 
 } 
 do.call(rbind, res) 
 } 
 
 all.combs(20, 5, 13) 
 
 Jean 
 
 
 wwreith reith_will...@bah.com wrote on 07/20/2012 07:45:30 AM:
 
  General problem: I have 20 projects that can be invested in and I need 
to
  decide which combinations meet a certain set of standards. The total
  possible combinations comes out to 2^20. However I know for a fact 
that the
  number of projects must be greater than 5 and less than 13. So far the 
the
  code below is the best I can come up with for iteratively creating a 
set to
  check against my set of standards.
  
  Code
  x-matrix(0,nrow=1,ncol=20)
  for(i in 1:2^20)
  {
  x[1]-x[1]+1
for(j in 1:20)
{
  if(x[j]1)
  {
x[j]=0
if(j20)
{
  x[j+1]=x[j+1]+1
}
  }
}
  if(sum(x)5  sum(x)13)
  {
  # insert criteria here.
  }
  }
  
  my code forces me to create all 2^20 x's and then use an if statement 
to
  decide if x is within my range of projects. Is there a faster way to
  increment x. Any ideas on how to kill the for loop so that it won't 
attempt
  to process an x where the sum is greater than 12 or less than 6?
[[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] Speeding up a loop

2012-07-20 Thread Jean V Adams
Reith, William [USA] reith_will...@bah.com wrote on 07/20/2012 
09:52:02 AM:

 Would this matrix eat up memory making the rest of my program 
 slower? Each x needs to be multiplied by a matrix and the results 
 checked against a set of thresholds. Doing them one at a time takes 
 at least 24 hours right now. 
 
 Optimizing a program is not my thing. 
 
 Sent from my Verizon Wireless 4GLTE smartphone


It's not mine either.
Better to ask the group.
I'm ccing R-help on this message.

Jean


 Jean V Adams jvad...@usgs.gov wrote on 07/20/2012 10:05 AM:
 
 I've had to do something similar, so I wrote a small function to help. 
 This runs in about 1/4 the time of your code on my machine. 
 Others may have a more efficient approach. 
 
 all.combs - function(num, from=0, to=num) { 
 # create a matrix of all possible combinations of num items 
 # restrict each row to have between from and to items 
 res - vector(list, to-from+1) 
 for(i in seq(from:to)) { 
 j - (from:to)[i] 
 if(j==0) res[[i]] - rep(FALSE, num) 
 comb - combn(num, j) 
 res[[i]] - t(apply(comb, 2, function(x) !is.na
 (match(1:num, x 
 } 
 do.call(rbind, res) 
 } 
 
 all.combs(20, 5, 13) 
 
 Jean 
 
 
 wwreith reith_will...@bah.com wrote on 07/20/2012 07:45:30 AM:
 
  General problem: I have 20 projects that can be invested in and I need 
to
  decide which combinations meet a certain set of standards. The total
  possible combinations comes out to 2^20. However I know for a fact 
that the
  number of projects must be greater than 5 and less than 13. So far the 
the
  code below is the best I can come up with for iteratively creating a 
set to
  check against my set of standards.
  
  Code
  x-matrix(0,nrow=1,ncol=20)
  for(i in 1:2^20)
  {
  x[1]-x[1]+1
for(j in 1:20)
{
  if(x[j]1)
  {
x[j]=0
if(j20)
{
  x[j+1]=x[j+1]+1
}
  }
}
  if(sum(x)5  sum(x)13)
  {
  # insert criteria here.
  }
  }
  
  my code forces me to create all 2^20 x's and then use an if statement 
to
  decide if x is within my range of projects. Is there a faster way to
  increment x. Any ideas on how to kill the for loop so that it won't 
attempt
  to process an x where the sum is greater than 12 or less than 6?
[[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] Changing ungrouped cases to grouped cases

2012-07-20 Thread Christopher Desjardins
Thanks the aggregate() command is what I was looking for.
Chris

On Thu, Jul 19, 2012 at 9:03 PM, David L Carlson dcarl...@tamu.edu wrote:

  dtf - read.table(text=y A   B   C
 + 0 11   2
 + 0 12   1
 + 1 11   2
 + 0 11   2
 + 1 11   2
 + 1 12   1
 + 0 12   2,
 + header=TRUE)
  dtagroup - aggregate(y~A+B+C, dtf, sum)
 # Gets you the groups. If you need the column/row order:

  dtagroup - dtagroup[order(dtagroup$y, decreasing=TRUE),c(4, 1:3)]

 --
 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-bounces@r-
  project.org] On Behalf Of Christopher Desjardins
  Sent: Thursday, July 19, 2012 7:35 PM
  To: R help
  Subject: [R] Changing ungrouped cases to grouped cases
 
  Hi,
  I have my data the following way:
 
  y A   B   C
  0 11   2
  0 12   1
  1 11   2
  0 11   2
  1 11   2
  1 12   1
  0 12   2
  .
  .
  .
  And so on.  How can I make my data look like the following:
  y   A  B  C
  2   1   1  2
  1   1   2  1
  0   1   2   2
  .
  .
  .
 
  In other words how can I change my ungrouped cases into grouped cases?
  Thanks!
  Chris
 
[[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] convert date to a factor

2012-07-20 Thread R. Michael Weylandt
d$date - factor(d$date)

Best,
Michael

On Fri, Jul 20, 2012 at 9:44 AM, Yolande Tra yolande@gmail.com wrote:
 Hello,

 I would like to convert date as a factor to represent time in a repeated
 measure situation in the following code. How would I do that?
  d - read.csv(file.path(dataDir,data.csv), as.is=T,stringsAsFactors =
 FALSE)
 d[1:2,]
 id date ab  c y
 1   1 8/6/2008 Red15 B  22
 2   1 8/6/2008   Green   15 B  22
 Thank you,
 Y

 [[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] Extracting standard errors for adjusted fixed effect sizes in lmer

2012-07-20 Thread Matthew Ouellette
Dear R help list,

I have done a lot of searching but have not been able to find an answer to
my problem.  I apologize in advance if this has been asked before.

I am applying a mixed model to my data using lmer.  I will use sample data
to illustrate my question:

library(lme4)
library(arm)
data(HR, package = SASmixed)
 str(HR)
'data.frame': 120 obs. of  5 variables:
 $ Patient: Factor w/ 24 levels 201,202,203,..: 1 1 1 1 1 2 2 2 2 2
...
 $ Drug   : Factor w/ 3 levels a,b,p: 3 3 3 3 3 2 2 2 2 2 ...
 $ baseHR : num  92 92 92 92 92 54 54 54 54 54 ...
 $ HR : num  76 84 88 96 84 58 60 60 60 64 ...
 $ Time   : num  0.0167 0.0833 0.25 0.5 1 ...

 fm1 - lmer(HR ~ baseHR + Time + Drug + (1 | Patient), HR)

 fixef(fm1)  ##Extract estimates of fixed effects

(Intercept)  baseHRTime   Drugb   Drugp

 32.6037923   0.5881895  -7.0272873   4.6795262  -1.0027581

 se.fixef(fm1)  ##Extract standard error of estimates of fixed effects

(Intercept)  baseHRTime   Drugb   Drugp

  9.9034008   0.1184529   1.4181457   3.5651679   3.5843026

##Because the estimate of the fixed effects are displayed as differences
from the intercept (I think?), I can back calculate the actual effect sizes
easily enough.  However, how would I do a similar calculation for the
standard error for these effect sizes (since these error estimates are for
the difference in means of effects) if my design isn't balanced (which
confuses things tremendously when working with a data set as large as
mine)?  It may help to point out that I'm working with microarray data;
applying the same model for each gene (hundreds of genes total) across
multiple samples (hundreds of samples total), but as an R beginner I like
to start with small data samples and work my way up.

I appreciate the help,

MO

[[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] Dissolve polygon

2012-07-20 Thread Celine
Hi,

I am working with  a SpatialPolygonsDataFrame of many islands. There are a
lot of polygons (islands) composing my SpatialPolygonsDataFrame.
I want to extract the elevation of each island.
I need to separate the different polygons (like dissolve function in
arcgis), to have the elevation of each island.
Do you have any idea how can I do that ? 
I already read a lot of forum, and read the maptools package but I did not
find something useful for my problem.
Coordinates:
min   max
x  40.32140 63.500179
y -25.60895 -3.711519
Is projected: FALSE 
proj4string : [+proj=longlat +ellps=WGS84 +no_defs]
Data attributes:
  NAME TYPE  
 Madagascar and the Indian Ocean Islands:1   hotspot_area:1  
NAME_TYPE
 Madagascar and the Indian Ocean Islands_hotspot_area:1  
 LEGEND 
 Madagascar and the Indian Ocean Islands:1  
 class(hot)
[1] SpatialPolygonsDataFrame
attr(,package)
[1] sp

Thanks a lot for your help,

Céline




--
View this message in context: 
http://r.789695.n4.nabble.com/Dissolve-polygon-tp4637234.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] Subsetting problem data, 2

2012-07-20 Thread arun
Hi,

Just a doubt regarding the dataset.

Suppose, I include two more patients F and G with different missing values as 
in this new dataset and run the code.
dat1-read.table(text=
Patient  Cycle  V1  V2  V3  V4  V5
A  1  0.4  0.1  0.5  1.5  NA
A  2  0.3  0.2  0.5  1.6  NA
A  3  0.3  NA  0.6  1.7  NA
A  4  0.4  NA  0.4  1.8  NA
A  5  0.5  0.2  0.5  1.5  NA
B  1  0.4  NA  NA  NA  NA
B  2  0.4  NA  NA  NA  NA
C  1  0.9  0.9  0.9  NA  NA
C  3  0.3  0.5  0.6  NA  NA
C  4  NA  NA  NA  NA  NA
C  5  0.4  NA  NA  NA  NA
D  1  0.2  0.5  NA  NA  NA
D  2  0.5  0.7  NA  NA  NA
D  4  0.6  0.4  NA  NA  NA
D  5  0.5  0.5  NA  NA  NA
E  1  0.1  NA  NA  NA  NA
E  2  0.5  0.3  NA  NA  NA
E  3  0.4  0.3  NA  NA  NA
F  1  0.2  NA   0.2 0.5 0.1  
F  2  0.5  NA   0.4 NA   0.3
F  3  0.6  NA   NA  0.3  0.2
G  1  0.2   0.5  NA  0.5  0.2
G  3  0.4   0.3  0.4 NA  0.3
G  4  0.6   0.2  0.2  0.4 NA
,sep=,header=TRUE)


nms - names(dat1)[grep(^V[1-9]$, names(dat1))]
dd - split(dat1, dat1$Patient)
fun - function(x) any(is.na(x))  any(!is.na(x))
ix - sapply(dd, function(x) Reduce(`|`, lapply(x[, nms], fun)))

dd[ix]
do.call(rbind, dd[ix])
 Patient Cycle  V1  V2  V3  V4  V5
A.1    A 1 0.4 0.1 0.5 1.5  NA
A.2    A 2 0.3 0.2 0.5 1.6  NA
A.3    A 3 0.3  NA 0.6 1.7  NA
A.4    A 4 0.4  NA 0.4 1.8  NA
A.5    A 5 0.5 0.2 0.5 1.5  NA
C.8    C 1 0.9 0.9 0.9  NA  NA
C.9    C 3 0.3 0.5 0.6  NA  NA
C.10   C 4  NA  NA  NA  NA  NA
C.11   C 5 0.4  NA  NA  NA  NA
E.16   E 1 0.1  NA  NA  NA  NA
E.17   E 2 0.5 0.3  NA  NA  NA
E.18   E 3 0.4 0.3  NA  NA  NA
F.19   F 1 0.2  NA 0.2 0.5 0.1
F.20   F 2 0.5  NA 0.4  NA 0.3
F.21   F 3 0.6  NA  NA 0.3 0.2
G.22   G 1 0.2 0.5  NA 0.5 0.2
G.23   G 3 0.4 0.3 0.4  NA 0.3
G.24   G 4 0.6 0.2 0.2 0.4  NA



Then, patients F and G are included in the list.  But, according to your 
initial statement, V1 and V2 are the most important variables.  If B is not 
included in the list because B has missing values for both cycles of B, then do 
you know think F or G should be included in the list.  Only difference is that 
F and G have missing values in other variables which do not behave 
consistently.  Do you have situations like that?

A.K.








- Original Message -
From: Lib Gray libgray3...@gmail.com
To: Rui Barradas ruipbarra...@sapo.pt
Cc: r-help r-help@r-project.org
Sent: Thursday, July 19, 2012 8:17 PM
Subject: Re: [R] Subsetting problem data, 2

I'm still getting the message (if this is what you were suggesting I try).
The data set I'm using has many more columns other than these variables;
could that be a problem? I didn't think it would affect it.

pattern - L[1-8][12]
 nms-names(data)[grep(vars,names(data))]
Warning message:
In grep(vars, names(data)) :
  argument 'pattern' has length  1 and only the first element will be used


On Thu, Jul 19, 2012 at 6:55 PM, Rui Barradas ruipbarra...@sapo.pt wrote:

 Hello,

 Sorry, forgot about that. It's trickier to write code without a dataset to
 test it.

 Try

 pattern - L[1-8][12]

 and after the grep print nms to see if it's right.

 Rui Barradas

 Em 20-07-2012 00:33, Lib Gray escreveu:

 I'm getting this error message:

 nms-names(data)[grep(vars,**names(data))]
 Warning message:
 In grep(vars, names(data)) :
    argument 'pattern' has length  1 and only the first element will be
 used

 Is there a way around this?


 On Thu, Jul 19, 2012 at 6:17 PM, Rui Barradas ruipbarra...@sapo.pt
 wrote:

  Hello,

 I guess so, and I can save you some typing.

 vars - sort(apply(expand.grid(L, 1:8, 1:2), 1, paste, collapse=))


 Then use it and see the result.

 Rui Barradas

 Em 20-07-2012 00:00, Lib Gray escreveu:

  The variables are actually L11, L12, L21, L22, ... , L81, L82. Would
 just
 creating a vector c(L11,... ,L82) be fine? (I'm about to try it, but I
 wanted to check to see if that was going to be a big issue).

 On Thu, Jul 19, 2012 at 3:27 PM, Rui Barradas ruipbarra...@sapo.pt
 wrote:

   Hello,

 Try the following. The data is your example of Patient A through E, but
 from the output of dput().

 dat - structure(list(Patient = structure(c(1L, 1L, 1L, 1L, 1L, 2L,
 2L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 5L, 5L, 5L), .Label = c(A,
 B, C, D, E), class = factor), Cycle = c(1L, 2L, 3L,
 4L, 5L, 1L, 2L, 1L, 3L, 4L, 5L, 1L, 2L, 4L, 5L, 1L, 2L, 3L),
       V1 = c(0.4, 0.3, 0.3, 0.4, 0.5, 0.4, 0.4, 0.9, 0.3, NA, 0.4,
       0.2, 0.5, 0.6, 0.5, 0.1, 0.5, 0.4), V2 = c(0.1, 0.2, NA,
       NA, 0.2, NA, NA, 0.9, 0.5, NA, NA, 0.5, 0.7, 0.4, 0.5, NA,
       0.3, 0.3), V3 = c(0.5, 0.5, 0.6, 0.4, 0.5, NA, NA, 0.9, 0.6,
       NA, NA, NA, NA, NA, NA, NA, NA, NA), V4 = c(1.5, 1.6, 1.7,
       1.8, 1.5, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
       NA), V5 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
       NA, NA, NA, NA, NA, NA)), .Names = c(Patient, Cycle,
 V1, V2, V3, V4, V5), class = data.frame, row.names = c(NA,
 

Re: [R] TukeyHSD not working

2012-07-20 Thread arun
Hi,

Check this link.
http://stackoverflow.com/questions/5904513/tukeyhsd-after-within-factors-anova
A.K.



- Original Message -
From: Michael Eisenring eimic...@ethz.ch
To: r-help@r-project.org
Cc: 
Sent: Friday, July 20, 2012 6:15 AM
Subject: [R] TukeyHSD not working

Dear r-help members.
I would like to compare species numbers of moths between eight different 
forests (each sampled for six nights). I would like to do a nested anova to 
compare species numbers between forests and nights.
For more site specific details I wanted to do a Tukey test (TukeyHSD). 
Unfortunately this test is not working (message:  nicht anwendbare Methode für 
'TukeyHSD' auf Objekt der Klasse c('anova', 'data.frame ) /(message: method 
not appropriate for TukeyHSD for objects of the class anova, data.frame).
I tried it also with the aov instead of the anova function but neither the 
aov nor the anova approach is working.
Could anyone help me with this  problem?
Attached you can find the txt.file and the r-script.
I'm working with a MAC OSX snow leopard, with R-Studio 0.96 316.

Thank you very much.
Michael 
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] r-h...@stat.math.ethz.ch

2012-07-20 Thread Norman....
Contact me with your full names, phone numbers and residential address if you 
own this email r-h...@stat.math.ethz.ch. 

I have a proposal for you

 

 

[[alternative HTML version deleted]]

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


Re: [R] Crosstab with Average and Count

2012-07-20 Thread arun
Hi,

Try this:
dat1-data.frame(x,y,z)


 xtabs(z~x+y,aggregate(z~x+y,dat1,mean))
   y
x    10  20  30
  1 100   0   0
  2   0 200   0
  3   0   0 300

table(dat1$z,dat1$y)
 
  10 20 30
  100  2  0  0
  200  0  3  0
  300  0  0  3
or,
table(dat1$x,dat1$z)
   
    100 200 300
  1   2   0   0
  2   0   3   0
  3   0   0   3




A.K.



- Original Message -
From: vioravis viora...@gmail.com
To: r-help@r-project.org
Cc: 
Sent: Friday, July 20, 2012 6:30 AM
Subject: [R] Crosstab with Average and Count

I have the following data:

x - as.factor(c(1,1,1,2,2,2,3,3,3))
y - as.factor(c(10,10,10,20,20,20,30,30,30))
z - c(100,100,NA,200,200,200,300,300,300)

I could create the cross tab of x and y with Sum of z as its elements using
the xtabs function as follows:

# X Vs. Y with Sum Z

xtabs(z ~ x + y)

   y
x    10  20  30
  1 200   0   0
  2   0 600   0
  3   0   0 900

How do I replace the sum with average and count so that I can get the
following outputs??

# X Vs. Y with Average of Z
   y
x      10  20  30
  1    100 0   0
  2    0   200 0
  3    0   0   300

# X Vs. Y with Count Z
  y
x    10  20  30
1    2   0   0
2    0   3   0
3    0   0   3

Would appreciate any help on these? Thank you.

Ravi





--
View this message in context: 
http://r.789695.n4.nabble.com/Crosstab-with-Average-and-Count-tp4637180.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] TukeyHSD not working

2012-07-20 Thread arun
HI,

I don't have much experience with laercio package.  If it works for you, it is 
good.

A.K.



- Original Message -
From: Michael Eisenring eimic...@ethz.ch
To: arun smartpink...@yahoo.com
Cc: 
Sent: Friday, July 20, 2012 8:38 AM
Subject: Re: [R] TukeyHSD not working

Hi
Thanks for the quick response. I installed the laercio package. Now I have a 
slightly different output of the Tukey Test but better than nothing.
Thank  you!
Michi
Am 20.07.2012 um 13:58 schrieb arun:

 Hi,
 
 Check this link.
 http://stackoverflow.com/questions/5904513/tukeyhsd-after-within-factors-anova
 A.K.
 
 
 
 - Original Message -
 From: Michael Eisenring eimic...@ethz.ch
 To: r-help@r-project.org
 Cc: 
 Sent: Friday, July 20, 2012 6:15 AM
 Subject: [R] TukeyHSD not working
 
 Dear r-help members.
 I would like to compare species numbers of moths between eight different 
 forests (each sampled for six nights). I would like to do a nested anova to 
 compare species numbers between forests and nights.
 For more site specific details I wanted to do a Tukey test (TukeyHSD). 
 Unfortunately this test is not working (message:  nicht anwendbare Methode 
 für 'TukeyHSD' auf Objekt der Klasse c('anova', 'data.frame ) /(message: 
 method not appropriate for TukeyHSD for objects of the class anova, 
 data.frame).
 I tried it also with the aov instead of the anova function but neither 
 the aov nor the anova approach is working.
 Could anyone help me with this  problem?
 Attached you can find the txt.file and the r-script.
 I'm working with a MAC OSX snow leopard, with R-Studio 0.96 316.
 
 Thank you very much.
 Michael 
 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/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] Execute a function

2012-07-20 Thread carla moreira
Yes,

that's this.

Thank you very much.



2012/7/20 Peter Ehlers [via R] ml-node+s789695n463721...@n4.nabble.com

 On 2012-07-20 04:05, carla moreira wrote:

 
  Hi,
 
  I would like to evaluate a function, with 3 arguments, for instance,
 
  myfunc-function(a,b,c) { sqrt(a)-exp(b)+4*c
   }
 
  How to execute  myfunc(x,y,z), for all x, all y and all z, where x,y,z
 are
  vectors?

 Is this what you have in mind:

myfunc - function(a, b, c){ sqrt(a)-exp(b)+4*c }
myfunc2 - function(x){
  a - x[1]
  b - x[2]
  c - x[3]
  myfunc(a, b, c)
}

x - c(1, 4, 9)
y - 1:2
z - c(10, -10, 2, 20)
d - expand.grid(x, y, z)
d$value - apply(d, 1, myfunc2)

 ?

 Peter Ehlers

 
  Thank you very much in advance
 
 
 
  --
  View this message in context:
 http://r.789695.n4.nabble.com/Execute-a-function-tp4637182.html
  Sent from the R help mailing list archive at Nabble.com.
 
  __
  [hidden email] 
  http://user/SendEmail.jtp?type=nodenode=4637212i=0mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
  and provide commented, minimal, self-contained, reproducible code.
 

 __
 [hidden email] http://user/SendEmail.jtp?type=nodenode=4637212i=1mailing 
 list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide
 http://www.R-project.org/posting-guide.html
 and provide commented, minimal, self-contained, reproducible code.


 --
  If you reply to this email, your message will be added to the discussion
 below:
 http://r.789695.n4.nabble.com/Execute-a-function-tp4637182p4637212.html
  To unsubscribe from Execute a function, click 
 herehttp://r.789695.n4.nabble.com/template/NamlServlet.jtp?macro=unsubscribe_by_codenode=4637182code=Y2FybGFtZ21tQGdtYWlsLmNvbXw0NjM3MTgyfDE3NjAxODQ0MA==
 .
 NAMLhttp://r.789695.n4.nabble.com/template/NamlServlet.jtp?macro=macro_viewerid=instant_html%21nabble%3Aemail.namlbase=nabble.naml.namespaces.BasicNamespace-nabble.view.web.template.NabbleNamespace-nabble.view.web.template.NodeNamespacebreadcrumbs=notify_subscribers%21nabble%3Aemail.naml-instant_emails%21nabble%3Aemail.naml-send_instant_email%21nabble%3Aemail.naml




-- 
Carla Moreira




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

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


Re: [R] cenbox(): Changing Default x-axis Group Labels

2012-07-20 Thread MacQueen, Don
Hi Rich,

I don't have a cenbox() function that I can find, but your solution will
(probably? hopefully?) be along the lines of:


foo - boxplot( y ~ x, data=sdf, plot=FALSE)

foo$names - ifelse(foo$names, Label for TRUE, Label for FALSE)

bxp(foo)


where sdf is a dataframe containing your data and y and x are appropriate
variables in it.


-- 
Don MacQueen

Lawrence Livermore National Laboratory
7000 East Ave., L-627
Livermore, CA 94550
925-423-1062





On 7/18/12 2:30 PM, Rich Shepard rshep...@appl-ecosys.com wrote:

   I've looked at the lattice book and the 'R Graphics Cookbook' without
seeing how to change the labels along the x axis for groups in a box plot,
specifically cenbox().

   The attached example has a main and axes labels with default group
labels.
Please point me to a reference on how I can change 'FALSE' and 'TRUE' to
something more easily understood by viewers.

TIA,

Rich

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] convert date to a factor

2012-07-20 Thread Yolande Tra
Thanks,
Y

On Fri, Jul 20, 2012 at 11:47 AM, R. Michael Weylandt 
michael.weyla...@gmail.com wrote:

 d$date - factor(d$date)

 Best,
 Michael

 On Fri, Jul 20, 2012 at 9:44 AM, Yolande Tra yolande@gmail.com
 wrote:
  Hello,
 
  I would like to convert date as a factor to represent time in a repeated
  measure situation in the following code. How would I do that?
   d - read.csv(file.path(dataDir,data.csv), as.is=T,stringsAsFactors
 =
  FALSE)
  d[1:2,]
  id date ab  c y
  1   1 8/6/2008 Red15 B  22
  2   1 8/6/2008   Green   15 B  22
  Thank you,
  Y
 
  [[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.htmlhttp://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] cenbox(): Changing Default x-axis Group Labels

2012-07-20 Thread Rich Shepard

On Fri, 20 Jul 2012, MacQueen, Don wrote:


I don't have a cenbox() function that I can find, but your solution will
(probably? hopefully?) be along the lines of:


Don,

  Well, I must have mis-typed that although I'm sure I read about it in the
NADA.pdf or Dennis' book. I'll look again. I don't get an error when calling
this method. Will look further and provide more information.


foo - boxplot( y ~ x, data=sdf, plot=FALSE)
foo$names - ifelse(foo$names, Label for TRUE, Label for FALSE)
bxp(foo)
where sdf is a dataframe containing your data and y and x are appropriate
variables in it.


  I will definitely work more on this. Right now I need to get my client to
explain the many discrepancies in the data they sent before I go further
with analyses.

Thanks very much,

Rich

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] convert date to a factor

2012-07-20 Thread Yolande Tra
Thanks,
Y

On Fri, Jul 20, 2012 at 1:02 PM, arun smartpink...@yahoo.com wrote:

 Hi,

 Try this:
 dat1-read.table(text=
iddateab  cy
 1  1 8/6/2008Red15B  22
 2  1 8/6/2008  Green  15B  22
 ,sep=,header=TRUE)
  dat1$date-with(dat1,as.factor(date))
  dat1
   id date a  b c  y
 1  1 8/6/2008   Red 15 B 22
 2  1 8/6/2008 Green 15 B 22
  is.factor(dat1$date)
 [1] TRUE
 A.K.



 - Original Message -
 From: Yolande Tra yolande@gmail.com
 To: r-help@r-project.org
 Cc:
 Sent: Friday, July 20, 2012 10:44 AM
 Subject: [R] convert date to a factor

 Hello,

 I would like to convert date as a factor to represent time in a repeated
 measure situation in the following code. How would I do that?
   d - read.csv(file.path(dataDir,data.csv), as.is=T,stringsAsFactors =
 FALSE)
  d[1:2,]
 id date ab  c y
 1   1 8/6/2008 Red15 B  22
 2   1 8/6/2008   Green   15 B  22
 Thank you,
 Y

  [[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.htmlhttp://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] any package available to calculate systematic root mean square error (RMSEs)?

2012-07-20 Thread Peng Zhao
Hi. Is there any package available to calculate the systematic root mean
square error (RMSEs) and unsystematic RMSE (RMSEu)? I could not find such a
package. I am too lazy to calculate them step by step.

The calculation could be found in equation 3 and 4 in the Appendix A of the
following link:
http://www.cobblestoneconcepts.com/ucgis2summer/anderson/anderson.htm
http://www.cobblestoneconcepts.com/ucgis2summer/anderson/anderson.htm 



--
View this message in context: 
http://r.789695.n4.nabble.com/any-package-available-to-calculate-systematic-root-mean-square-error-RMSEs-tp4637247.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] convert date to a factor

2012-07-20 Thread arun
Hi,

Try this:
dat1-read.table(text=
   id    date    a    b  c    y
1  1 8/6/2008    Red    15    B  22
2  1 8/6/2008  Green  15    B  22
,sep=,header=TRUE)
 dat1$date-with(dat1,as.factor(date))
 dat1
  id date a  b c  y
1  1 8/6/2008   Red 15 B 22
2  1 8/6/2008 Green 15 B 22
 is.factor(dat1$date)
[1] TRUE
A.K.



- Original Message -
From: Yolande Tra yolande@gmail.com
To: r-help@r-project.org
Cc: 
Sent: Friday, July 20, 2012 10:44 AM
Subject: [R] convert date to a factor

Hello,

I would like to convert date as a factor to represent time in a repeated
measure situation in the following code. How would I do that?
  d - read.csv(file.path(dataDir,data.csv), as.is=T,stringsAsFactors =
FALSE)
 d[1:2,]
        id     date         a        b      c         y
1       1 8/6/2008     Red    15     B          22
2       1 8/6/2008   Green   15     B          22
Thank you,
Y

    [[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] findAssocs()

2012-07-20 Thread rtw30606
Hi

Here is some code to illustrate how the correlations are calculated.

 data -  c(word1, word1 word2,word1 word2 word3,word1 word2 word3
 word4,word1 word2 word3 word4 word5)
 frame -  data.frame(data)
 frame
   data
1 word1
2   word1 word2
3 word1 word2 word3
4   word1 word2 word3 word4
5 word1 word2 word3 word4 word5
 test -  Corpus(DataframeSource(frame, encoding = UTF-8))
 dtm -  DocumentTermMatrix(test)
 as.matrix(dtm)
Terms
Docs word1 word2 word3 word4 word5
   1 1 0 0 0 0
   2 1 1 0 0 0
   3 1 1 1 0 0
   4 1 1 1 1 0
   5 1 1 1 1 1
 
 findAssocs(dtm, word2, 0.1)
word2 word3 word4 word5 
 1.00  0.61  0.41  0.25 
 # Correlation word2 with word3
 cor(c(0,1,1,1,1),c(0,0,1,1,1))
[1] 0.6123724
 # Correlation word2 with word4
 cor(c(0,1,1,1,1),c(0,0,0,1,1))
[1] 0.4082483
 # Correlation word2 with word5
 cor(c(0,1,1,1,1),c(0,0,0,0,1))
[1] 0.25

Cheers

Rick




--
View this message in context: 
http://r.789695.n4.nabble.com/findAssocs-tp3845751p4637248.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] R packages installation error in Ubuntu

2012-07-20 Thread uday
hi ,

last time I did not mention the specific version of linux and it was problem
with only syntax but with correction I got new error. 
I found that some people has  faced this problem before and some of them
have mentioned that it may be because of dependencies issue and  
install.packages(Rcmdr)  could fix this problem. 
But Still I have not fixed the problem. 

could you please tell me how to resolve this issue, this is the first time I
am using R in ubuntu. 
Thank you 





--
View this message in context: 
http://r.789695.n4.nabble.com/R-packages-installation-error-in-Ubuntu-tp4637105p4637250.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] cenbox(): Changing Default x-axis Group Labels

2012-07-20 Thread Peter Ehlers

On 2012-07-20 09:41, Rich Shepard wrote:

On Fri, 20 Jul 2012, MacQueen, Don wrote:


I don't have a cenbox() function that I can find, but your solution will
(probably? hopefully?) be along the lines of:


Don,

Well, I must have mis-typed that although I'm sure I read about it in the
NADA.pdf or Dennis' book. I'll look again. I don't get an error when calling
this method. Will look further and provide more information.


foo - boxplot( y ~ x, data=sdf, plot=FALSE)
foo$names - ifelse(foo$names, Label for TRUE, Label for FALSE)
bxp(foo)
where sdf is a dataframe containing your data and y and x are appropriate
variables in it.


I will definitely work more on this. Right now I need to get my client to
explain the many discrepancies in the data they sent before I go further
with analyses.


Well, You didn't say (in your original request) that you were using
the NADA package. The function is cenboxplot() and it's just a
wrapper for boxplot() and hence passes arguments to boxplot().
Thus the solution to your problem is just to add the 'names='
argument, as in (using the example in ?cenboxplot):

  with(Golden, cenboxplot(Blood, BloodCen, DosageGroup,
   names = c(yabba, doo)))

Peter Ehlers
p.s. I can sympathize with your client data troubles. 'Tis ever thus.



Thanks very much,

Rich

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

2012-07-20 Thread Barry Rowlingson
On Fri, Jul 20, 2012 at 4:51 PM, Celine bellard.cel...@gmail.com wrote:
 Hi,

 I am working with  a SpatialPolygonsDataFrame of many islands. There are a
 lot of polygons (islands) composing my SpatialPolygonsDataFrame.
 I want to extract the elevation of each island.

 Is each island a row in your SPDF? Or is it just one row, with lots
of island polygons? What does dim(hot) (assuming 'hot' is your SPDF)
say?

 I need to separate the different polygons (like dissolve function in
 arcgis), to have the elevation of each island.

 If each island is a row, then its easy, if each island is a separate
loop in a single row, a bit trickier.

 How are you going to get the elevation? Do you have a separate data
source, such as a grid of elevations?

 Also, try hopping over to the R-sig-geo mailing list. For everything spatial.

 Barry

blog: http://geospaced.blogspot.com/
web: http://www.maths.lancs.ac.uk/~rowlings
web: http://www.rowlingson.com/
twitter: http://twitter.com/geospacedman
pics: http://www.flickr.com/photos/spacedman

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

2012-07-20 Thread Rui Barradas

Hello,

Yes, Carol is comparing what can't be compared. Your code should do it, 
I hope.


Rui Barradas

Em 20-07-2012 15:12, David L Carlson escreveu:

Rui's example included z-score data (drawn from rnorm). You converted your
data to z-scores so you need to compare your results to the z-scores not
the original data.

Change these lines:

tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)*sign(scale(tmp))
# sign is of scale(tmp) not tmp
equal(scale(tmp), tmp.qnorm)
# compare to scale(tmp) not tmp

--
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-bounces@r-
project.org] On Behalf Of carol white
Sent: Friday, July 20, 2012 7:29 AM
To: Rui Barradas
Cc: r-help@r-project.org
Subject: Re: [R] function for inverse normal transformation

Thanks Rui.
I changed my scripts to the followings and I think that it still is not
correct. See also the attached file.

Thanks for your help,


  tmp
  [1]  2.502519  1.828576  3.755778 17.415000  3.779296  2.956850
2.379663
  [8]  1.103559  8.920316  2.744500  2.938480  7.522174 10.629200
8.552259
[15]  5.425938  4.388906  0.00  0.723887 11.337860  3.763786


  tmp.p =2*pnorm(abs(scale(tmp)),lower.tail=FALSE)

   tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)
   tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)*sign(tmp)
equal(tmp, tmp.qnorm)

[1] FALSE

par(mfrow = c(1,3))
hist(tmp)
hist(tmp.p)
hist(tmp.qnorm)




  From: Rui Barradas ruipbarra...@sapo.pt
To: carol white wht_...@yahoo.com
Cc: r-help r-help@r-project.org
Sent: Friday, July 20, 2012 2:02 PM
Subject: Re: [R] function for inverse normal transformation


Hello,

No it's not correct, you are computing a what seems to be a
 two-tailed probabiity, so the inverse should account for it. Look
 closely: you take the absolute value, then the upper tail
 probability, then multiply 2 into it. Reverse these steps to get
the
 correct value.

# Helper function
equal - function(x, y, tol=.Machine$double.eps^0.5) all(abs(x -
 y)   tol)

m - rnorm(5)
p - 2*pnorm(abs(m), lower.tail=FALSE)
m2 - qnorm(p/2, lower.tail=FALSE)*sign(m)

equal(m, m2)

(The helper function is just to test floating point values computed
 differently for equality.)

Hope this helps,

Rui Barradas


Em 20-07-2012 12:36, carol white escreveu:

Thanks for your reply. So to derive it from a given data set, is the
following correct to do? my_data.p
=2*pnorm(abs(my_data),lower.tail=FALSE) my_data.q = qnorm(my_data.p)
Cheers,  From: Duncan Murdoch
murdoch.dun...@gmail.com Cc: r-h...@stat.math.ethz.ch r-
h...@stat.math.ethz.ch Sent: Friday, July 20, 2012 1:23 PM
Subject: Re: [R] function for inverse normal transformation On 12-07-20
6:21 AM, carol white wrote:

Hi,

What is the function for inverse normal transformation?

qnorm Duncan Murdoch
Thanks, Carol [[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] Speeding up a loop

2012-07-20 Thread wwreith
That is faster than what I was doing and reducing 15% of my iterations it
still very helpful.

Next question.

I need to multiply each row x[i,] of the matrix x by another matrix A.
Specifically

for(i in 1:n)
{
If (x[i,]%*%A[,1].5 || x[i,]%*%A[,2]42 || x[i,]%*%A[,3]150) 
{
x-x[-i,] 
n-n-1
}. #In other words remove row i from x if it does not meet criteria (=.5,
=42, =150). When multiplied to A
}
Is there a better way than using a for loop for this or x-x[-i,] for that
matter? I assume building a new matrix would be worse. 

Ideally I want to also exclude some x[,i] as well example if x[1,] is better
than x[2,] in all three categories i.e. bigger, bigger, and smaller than
x[2,] when multiplied to A then I want to exclude x[2,] as well. Any
suggestions on whether it is better to do this all at once or in stages?

Thanks for helping!



--
View this message in context: 
http://r.789695.n4.nabble.com/Speeding-up-a-loop-tp4637201p4637255.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] Speeding up a loop

2012-07-20 Thread Richard M. Heiberger
This works to multiply the ith row of a by the ith value of b.
It might be what you can use

a - matrix(1:30, 6, 5)
b - 1:6

a
a*b



To simplify your code, I think you can do this---one multiplication

xA - x %*% A

Now you can do the tests on xA and not have any matrix multiplications
inside the loop.
The next line is not tested but I hope gives you the idea

remove.set - (xA[,1]  .5) || (xA[,2]  52) || (xA[,3]  150)

new.x - x[remove.set,]

You now have new.x with no explicit loops at all, because this takes
advantage of vector operations.

Rich

On Fri, Jul 20, 2012 at 2:34 PM, wwreith reith_will...@bah.com wrote:

 That is faster than what I was doing and reducing 15% of my iterations it
 still very helpful.

 Next question.

 I need to multiply each row x[i,] of the matrix x by another matrix A.
 Specifically

 for(i in 1:n)
 {
 If (x[i,]%*%A[,1].5 || x[i,]%*%A[,2]42 || x[i,]%*%A[,3]150)
 {
 x-x[-i,]
 n-n-1
 }. #In other words remove row i from x if it does not meet criteria (=.5,
 =42, =150). When multiplied to A
 }
 Is there a better way than using a for loop for this or x-x[-i,] for that
 matter? I assume building a new matrix would be worse.

 Ideally I want to also exclude some x[,i] as well example if x[1,] is
 better
 than x[2,] in all three categories i.e. bigger, bigger, and smaller than
 x[2,] when multiplied to A then I want to exclude x[2,] as well. Any
 suggestions on whether it is better to do this all at once or in stages?

 Thanks for helping!



 --
 View this message in context:
 http://r.789695.n4.nabble.com/Speeding-up-a-loop-tp4637201p4637255.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.


Re: [R] Changing ungrouped cases to grouped cases

2012-07-20 Thread Christopher Desjardins
As a follow up is there any way to also get the count for each combination?
For example

 y A   B   C
 0 11   2
 0 12   1
 1 11   2
 0 11   2
 1 11   2
 1 12   1
 0 12   2

Should become:

 y   A  B  C  count
 2   1   1  24
 1   1   2  12
 0   1   2   2   1
.
.
.
So I would know there were 2 successes out of 4.
Thanks!
Chris



On Fri, Jul 20, 2012 at 10:41 AM, Christopher Desjardins 
cddesjard...@gmail.com wrote:

 Thanks the aggregate() command is what I was looking for.
 Chris


 On Thu, Jul 19, 2012 at 9:03 PM, David L Carlson dcarl...@tamu.eduwrote:

  dtf - read.table(text=y A   B   C
 + 0 11   2
 + 0 12   1
 + 1 11   2
 + 0 11   2
 + 1 11   2
 + 1 12   1
 + 0 12   2,
 + header=TRUE)
  dtagroup - aggregate(y~A+B+C, dtf, sum)
 # Gets you the groups. If you need the column/row order:

  dtagroup - dtagroup[order(dtagroup$y, decreasing=TRUE),c(4, 1:3)]

 --
 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-bounces@r-
  project.org] On Behalf Of Christopher Desjardins
  Sent: Thursday, July 19, 2012 7:35 PM
  To: R help
  Subject: [R] Changing ungrouped cases to grouped cases
 
  Hi,
  I have my data the following way:
 
  y A   B   C
  0 11   2
  0 12   1
  1 11   2
  0 11   2
  1 11   2
  1 12   1
  0 12   2
  .
  .
  .
  And so on.  How can I make my data look like the following:
  y   A  B  C
  2   1   1  2
  1   1   2  1
  0   1   2   2
  .
  .
  .
 
  In other words how can I change my ungrouped cases into grouped cases?
  Thanks!
  Chris
 
[[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] [RsR] How does rlm in R decide its w weights for each IRLS iteration?

2012-07-20 Thread Valentin Todorov
Hi Michael, S Ellison,

I do not actually understand what you want to achieve with the M
estimates of rlm in MASS, but why you do not give a try of lmrob in
'robustbase'. Please have a llok in the references (?lmrob) about the
advantages of MM estimators over the M estimators.

Best regards,
Valentin




On Fri, Jul 20, 2012 at 5:11 PM, S Ellison s.elli...@lgcgroup.com wrote:


 -Original Message-
 Subject: [RsR] How does rlm in R decide its w weights for
 each IRLS iteration?
 I am also confused about the manual:

a.  The input arguments:

 wt.method are the weights case weights (giving the relative
 importance of case, so a weight of 2 means there are two of
 these) or the inverse of the variances, so a weight of two
 means this error is half as variable?

 When you give rlm weights (called 'weights', not 'w' on input, though you can 
 abbreviate to 'w'), you need to tell it which of these two possibilities you 
 used.
 If you gave it case numbers, say wt.method=case; if you gave it inverse 
 variance weights, say wt.method=inv.var.
 The default is inv.var.


 The input argument w is used for the initial values of the
 rlm IRLS weighting and the output value w is the converged w.
 There is no input argument 'w' for rlm (see above).
 The output w are a  calculated using the psi function, so between 0 and 1.
 The effective weights for the final estimate would then be something like 
 w*weights, using the full name of the input argument (and if I haven't 
 forgotten a square root somewhere). At least, that would work for a simple 
 location estimate (eg rlm(x~1)).

 If my understanding above is correct, how does rlm decide
 its w for each IRLS iteration then?
 It uses the given psi functions to calculate the iterative weights based on 
 the scaled residuals.

 Any pointers/tutorials/notes to the calculation of these
 w's in each IRLS iteration?
 Read the cited references for a detailed guide. Or, of course, MASS - the 
 package is, after all, intended to support the book, not replace it.



 S Ellison

 ***
 This email and any attachments are confidential. Any u...{{dropped:6}}

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


Re: [R] Changing ungrouped cases to grouped cases

2012-07-20 Thread Rui Barradas

Hello,

Still with aggregate, use length() not sum().


dtagroup - aggregate(y ~ A + B + C, data=dtf, sum)
dtalength - aggregate(y ~ A + B + C, data=dtf, length)

# Now merge the two
names(dtalength)[4] - count
mm - merge(dtagroup, dtalength)

# And make it pretty
mm - mm[, c(y, A, B, C, count)]
mm

Hope this helps,

Rui Barradas

Em 20-07-2012 19:52, Christopher Desjardins escreveu:

As a follow up is there any way to also get the count for each combination?
For example

  y A   B   C
  0 11   2
  0 12   1
  1 11   2
  0 11   2
  1 11   2
  1 12   1
  0 12   2

Should become:

  y   A  B  C  count
  2   1   1  24
  1   1   2  12
  0   1   2   2   1
.
.
.
So I would know there were 2 successes out of 4.
Thanks!
Chris



On Fri, Jul 20, 2012 at 10:41 AM, Christopher Desjardins 
cddesjard...@gmail.com wrote:


Thanks the aggregate() command is what I was looking for.
Chris


On Thu, Jul 19, 2012 at 9:03 PM, David L Carlson dcarl...@tamu.eduwrote:


dtf - read.table(text=y A   B   C

+ 0 11   2
+ 0 12   1
+ 1 11   2
+ 0 11   2
+ 1 11   2
+ 1 12   1
+ 0 12   2,
+ header=TRUE)

dtagroup - aggregate(y~A+B+C, dtf, sum)

# Gets you the groups. If you need the column/row order:


dtagroup - dtagroup[order(dtagroup$y, decreasing=TRUE),c(4, 1:3)]

--
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-bounces@r-
project.org] On Behalf Of Christopher Desjardins
Sent: Thursday, July 19, 2012 7:35 PM
To: R help
Subject: [R] Changing ungrouped cases to grouped cases

Hi,
I have my data the following way:

y A   B   C
0 11   2
0 12   1
1 11   2
0 11   2
1 11   2
1 12   1
0 12   2
.
.
.
And so on.  How can I make my data look like the following:
y   A  B  C
2   1   1  2
1   1   2  1
0   1   2   2
.
.
.

In other words how can I change my ungrouped cases into grouped cases?
Thanks!
Chris

   [[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] Changing ungrouped cases to grouped cases

2012-07-20 Thread Peter Ehlers

On 2012-07-20 12:09, Rui Barradas wrote:

Hello,

Still with aggregate, use length() not sum().


dtagroup - aggregate(y ~ A + B + C, data=dtf, sum)
dtalength - aggregate(y ~ A + B + C, data=dtf, length)

# Now merge the two
names(dtalength)[4] - count
mm - merge(dtagroup, dtalength)

# And make it pretty
mm - mm[, c(y, A, B, C, count)]
mm

Hope this helps,

Rui Barradas


Package plyr's ddply() with summarize does this sort
of thing nicely:

 library(plyr)
 ddply(dtf, .(A,B,C), summarize,
  y = sum(y),
  count = length(y))

Peter Ehlers



Em 20-07-2012 19:52, Christopher Desjardins escreveu:

As a follow up is there any way to also get the count for each combination?
For example

   y A   B   C
   0 11   2
   0 12   1
   1 11   2
   0 11   2
   1 11   2
   1 12   1
   0 12   2

Should become:

   y   A  B  C  count
   2   1   1  24
   1   1   2  12
   0   1   2   2   1
.
.
.
So I would know there were 2 successes out of 4.
Thanks!
Chris



On Fri, Jul 20, 2012 at 10:41 AM, Christopher Desjardins 
cddesjard...@gmail.com wrote:


Thanks the aggregate() command is what I was looking for.
Chris


On Thu, Jul 19, 2012 at 9:03 PM, David L Carlson dcarl...@tamu.eduwrote:


dtf - read.table(text=y A   B   C

+ 0 11   2
+ 0 12   1
+ 1 11   2
+ 0 11   2
+ 1 11   2
+ 1 12   1
+ 0 12   2,
+ header=TRUE)

dtagroup - aggregate(y~A+B+C, dtf, sum)

# Gets you the groups. If you need the column/row order:


dtagroup - dtagroup[order(dtagroup$y, decreasing=TRUE),c(4, 1:3)]

--
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-bounces@r-
project.org] On Behalf Of Christopher Desjardins
Sent: Thursday, July 19, 2012 7:35 PM
To: R help
Subject: [R] Changing ungrouped cases to grouped cases

Hi,
I have my data the following way:

y A   B   C
0 11   2
0 12   1
1 11   2
0 11   2
1 11   2
1 12   1
0 12   2
.
.
.
And so on.  How can I make my data look like the following:
y   A  B  C
2   1   1  2
1   1   2  1
0   1   2   2
.
.
.

In other words how can I change my ungrouped cases into grouped cases?
Thanks!
Chris

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



__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] cenbox(): Changing Default x-axis Group Labels

2012-07-20 Thread Rich Shepard

On Fri, 20 Jul 2012, Peter Ehlers wrote:


Well, You didn't say (in your original request) that you were using
the NADA package. The function is cenboxplot() and it's just a
wrapper for boxplot() and hence passes arguments to boxplot().
Thus the solution to your problem is just to add the 'names='
argument, as in (using the example in ?cenboxplot):

 with(Golden, cenboxplot(Blood, BloodCen, DosageGroup,
  names = c(yabba, doo)))


Peter,

  Thank you. I wrote off the top of my head and obviously missed the correct
function name.


p.s. I can sympathize with your client data troubles. 'Tis ever thus.


  Yep. They all use spreadsheets rather than databases and there's no
integrity checks. Sigh. Well, the pay's the same regardless.

Carpe weekend,

Rich

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

2012-07-20 Thread John Fox
Dear Carol,

 -Original Message-
 From: carol white [mailto:wht_...@yahoo.com]
 Sent: July-20-12 2:45 PM
 To: John Fox
 Subject: Re: inverse normal transformation
 
 Thanks John for your quick reply.
 
 The purpose of applying inverse normal transformation is to reduce the
 impact of outliers and deviations from normality on statistical
 analysis.

In other words, you're forcing the variable to follow a normal distribution
and making the units of measurement uninterpretable. I'll assume that this
somehow makes sense.

 
 Indeed, it includes the steps that you went through. However, I don't
 know why you calculated  (rank - 0.5)/20 to get the p-value. Then, how
 could we convert the quantiles (Q) into normal deviates?

They *are* quantiles on the standard normal scale -- that's what qnorm()
provides (with the default mean of 0 and standard deviation of 1). The
cumulative probabilities (not p-values) are calculated from the order
statistics of your data, where subtracting 0.5 avoids cumulative
probabilities of 0 or 1. This (or something close to it) is standard for
computing comparison quantiles.

I'm copying this message to r-help (with the original subject line) since
the discussion there continues.

Best,
 John

 
 Many thanks,
 
 Carol
 
 
 
 
 From: John Fox j...@mcmaster.ca
 To: 'carol white' wht_...@yahoo.com
 Sent: Friday, July 20, 2012 4:43 PM
 Subject: RE: inverse normal transformation
 
 
 Dear Carol,
 
 Like the people on r-help list who tried to help you, I have no idea
 why you
 want to do this. If you're trying to get the corresponding standard
 normal
 quantiles for your data, as for a QQ plot (and why else you might want
 them
 isn't clear to me), you can simply compute
 
 rank - rank(tmp)
 P - (rank - 0.5)/20
 Q - qnorm(P)
 
 Then, the QQ plot is
 
 plot(Q, tmp)
 
 Best,
 John
 
 
 John Fox
 Senator William McMaster
   Professor of Social Statistics
 Department of Sociology
 McMaster University
 Hamilton, Ontario, Canada
 http://socserv.mcmaster.ca/jfox
 
 
 
 
  -Original Message-
  From: carol white [mailto:wht_...@yahoo.com]
  Sent: July-20-12 9:08 AM
  To: j...@mcmaster.ca
  Subject: inverse normal transformation
 
  Dear John,
 
 
  Are the following scripts correct to get the inverse normal
  transformation of a data set?
 
 
  Thanks for your help,
 
 
  Carol
  ---
 
   tmp
   [1]  2.502519  1.828576  3.755778 17.415000  3.779296  2.956850
  2.379663  [8]  1.103559  8.920316  2.744500  2.938480  7.522174
  10.629200  8.552259 [15]  5.425938  4.388906  0.00  0.723887
  11.337860  3.763786
 
 
   tmp.p =2*pnorm(abs(scale(tmp)),lower.tail=FALSE)
tmp.qnorm = qnorm(tmp.p/2,lower.tail=FALSE)  tmp.qnorm =
   qnorm(tmp.p/2,lower.tail=FALSE)
 
   par(mfrow = c(1,3))
   hist(tmp)
   hist(tmp.p)
   hist(tmp.qnorm)
 
 
 
 


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