Re: [R] plm: plm.data vs pdata.frame

2009-05-10 Thread Achim Zeileis

Stephen:


I am trying to use the plm package for panel econometrics. I am just
trying to get started and load my data.


Look at the Journal of Statistical Software paper that introduces the 
"plm" package:

  http://www.jstatsoft.org/v27/i02/


It seems from most of the
sample documentation that I need to use the pdata.frame function to
get my data loaded. However, even after installing the "plm" package,
my R installation cannot find the function.


pdata.frame() was removed from the package. plm.data() took over some of 
its purposes. See the paper mentioned above for more details.



I am trying to follow the
example in plmEN.pdf ( cran.mirroring.de/doc/vignettes/plm/plmEN.pdf )


That version is two years old, plm changed (and improved) a lot since 
them. You shouldn't use that mirror in general (unless you want to use 
tools from 2007...).


Yves & Giovanni: The latest version of plm has lost its "inst" directory. 
Hence, there is no vignette and no pointer to the JSS paper via the 
CITATION. It seems that the inst directory wasn't moved to R-Forge...


hth,
Z

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


[R] Traceback NA's to their first appearance

2009-05-10 Thread Benjamin Leutner

Hello,

I've written a relatively complex interconnected population model in R. 
When changing a certain parameter, the outputs end up containing NA's.
I would like to find out, in which step the model (in form of a loop) 
starts to produce NA's. Does anyone know how to achieve this?

Since it does not give an error I can't use the traceback() function.

Thank you very much in advance
Ben

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Predictive Analytics Seminar: May 27-28, New York City

2009-05-10 Thread Elise Johnson

Hi all, 

I wanted to let you know about our training seminar on predictive analytics
- coming May, Oct, and Nov in NYC, Stockholm, DC and other cities.  This is
intensive training for marketers, managers and business professionals to
make actionable sense of customer data by predicting buying behavior, churn,
etc.  Past attendees provided rave reviews.

Here's more info:
--

Training Program: Predictive Analytics for Business, Marketing and Web

A two-day intensive seminar brought to you by Prediction Impact, Inc. 

Dates: May 27-28, Oct 14-15, Oct 18-19, and Nov 11-12, 2009
Locations: NYC (May), Stockholm (Oct), DC (Oct), San Francisco (Nov) 

93% rate this program Excellent or Very Good. 
**The official training program of Predictive Analytics World**
**Offered in conjunction with eMetrics events** 

Also see our Online Training: Predictive Analytics Applied - immediate
access at any time:
www.predictionimpact.com/predictive-analytics-online-training.html


ABOUT THIS SEMINAR:

Business metrics do a great job summarizing the past. But if you want to
predict how customers will respond in the future, there is one place to
turn--predictive analytics. By learning from your abundant historical data,
predictive analytics provides the marketer something beyond standard
business reports and sales forecasts: actionable predictions for each
customer. These predictions encompass all channels, both online and off,
foreseeing which customers will buy, click, respond, convert or cancel. If
you predict it, you own it. 

The customer predictions generated by predictive analytics deliver more
relevant content to each customer, improving response rates, click rates,
buying behavior, retention and overall profit. For online applications such
as e-marketing and customer care recommendations, predictive analytics acts
in real-time, dynamically selecting the ad, web content or cross-sell
product each visitor is most likely to click on or respond to, according to
that visitor's profile. This is AB selection, rather than just AB testing. 

Predictive Analytics for Business, Marketing and Web is a concentrated
training program that includes interactive breakout sessions and a brief
hands-on exercise. In two days we cover: 

- The techniques, tips and pointers you need in order to run a successful
predictive analytics and data mining initiative

- How to strategically position and tactically deploy predictive analytics
and data mining at your company

- How to bridge the prevalent gap between technical understanding and
practical use

- How a predictive model works, how it's created and how much revenue it
generates

- Several detailed case studies that demonstrate predictive analytics in
action and make the concepts concrete 

- NEW TOPIC: Five Ways to Lower Costs with Predictive Analytics 


No background in statistics or modeling is required. The only specific
knowledge assumed for this training program is moderate experience with
Microsoft Excel or equivalent. 

For more information, visit
www.predictionimpact.com/predictive-analytics-training.html, or e-mail us at
train...@predictionimpact.com.  You may also call (415) 683-1146. 

Cross-Registration Special: Attendees earn $250 off the Predictive Analytics
World Conference 

SNEAK PREVIEW VIDEO:
www.predictionimpact.com/predictive-analytics-times.html 

$100 off early registration, 3 weeks ahead

-- 
View this message in context: 
http://www.nabble.com/Predictive-Analytics-Seminar%3A-May-27-28%2C-New-York-City-tp23467596p23467596.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] Sweave \Sexpr{} advice please

2009-05-10 Thread Dieter Menne
Kyle Matoba  ucdavis.edu> writes:

> A question in my work I use \Sexpr{} statements scalar values and the xtable
> package for all manner of tables.  What I'd like to do is to use a vector
> inline, rather than a whole separate table.   Something like:
> 
> % Sweave block:
> <<>>=
> covmat <- cov(matrix(runif(100),ncol=3))
> @
> 
> % back to Latex, typing up a report, my homework, etc.
> The first column of the covariance matrix is $(\Sexpr{covmat[1,1]},
> \Sexpr{covmat[2,1]}, \Sexpr{covmat[3,1]})^T$
> 

To unclutter the text in Sexpr, I normally prepare the results in the 
Sweave block in advance

Dieter

covmat <- cov(matrix(runif(30),ncol=3))
ft = paste(round(covmat[,3]),collapse=", ")
ft
named = c(Val=3,StdDev=4,p=0.4)
ft = paste(names(named),named,sep="= ",collapse=", ")
ft

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


Re: [R] What does it mean by skip=2 and skip=7?

2009-05-10 Thread Dieter Menne
Tung86  gmail.com> writes:

> 
> Can anyone tell me what is skip=2, skip =7 and %in% mean here?
> 
> fromcsv=read.csv
> ('2_2005_top200_postdoc.csv',header=FALSE,skip=7,stringsAsFactors=FALSE)

Did you check the docs?

skipinteger: the number of lines of the data file to skip before beginning 
to
read data.
 
> fromreadxls = fromreadxls[fromreadxls$V7 %in% c('Public','Private'),]

See the documentation of "match"

Dieter

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


Re: [R] Sweave \Sexpr{} advice please

2009-05-10 Thread Tobias Verbeke

Hi Kyle,


First off, my deepest gratitude to the Sweave developers: this tool has
improved my quality greatly.

A question in my work I use \Sexpr{} statements scalar values and the xtable
package for all manner of tables.  What I'd like to do is to use a vector
inline, rather than a whole separate table.   Something like:

 begin code
% Latex junk

% Sweave block:
<<>>=
covmat <- cov(matrix(runif(100),ncol=3))
@

% back to Latex, typing up a report, my homework, etc.
The first column of the covariance matrix is $(\Sexpr{covmat[1,1]},
\Sexpr{covmat[2,1]}, \Sexpr{covmat[3,1]})^T$

% end code

but, of course, this is poor way of going about it.  Any suggestions?


<>=
  require(xtable)
@

<<>>=
covmat <- cov(matrix(runif(99), ncol=3)) # 99 not 100
@

The first row of the\dots is
<>=
xtable(covmat[1,,drop=FALSE])
@

This should get you started. Be sure to
explore the facilities of the xtable package
(you can use captions, labels etc.):

?xtable
?print.xtable

The results=tex option to the chunk is
important as well, of course.

HTH,
Tobias

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


Re: [R] Spatstat

2009-05-10 Thread Peter Dalgaard

Richard Chirgwin wrote:

Hi all,

I am trying to install Spatstat on OpenSUSE 11.1.
install.packages("spatstat", dependencies = TRUE)
fails on the basis of various compiler packages (full message below).

I have gcc version 4.3.2, which should include gfortran and g++ - so I'm not
sure what to do!


First check your assumptions...

If gfortran doesn't work on the command line, then you have an 
installation issue, nothing to do with R, much less spatstat.


Hint (on SUSE 11.0):

viggo:~/>rpm -qf `which gfortran`
gcc-fortran-4.3-39.1

-pd




Richard

* Installing *source* package ‘deldir’ ...
** libs
gfortran   -fpic  -O2 -c acchk.f -o acchk.o
make: gfortran: Command not found
make: *** [acchk.o] Error 127
ERROR: compilation failed for package ‘deldir’
* Removing ‘/home/richard/R/i686-pc-linux-gnu-library/2.9/deldir’
* Installing *source* package ‘spatstat’ ...
** libs
gcc -std=gnu99 -I/usr/lib/R/include  -I/usr/local/include-fpic  -O2 -c
Kborder.c -o Kborder.o
gcc -std=gnu99 -I/usr/lib/R/include  -I/usr/local/include-fpic  -O2 -c
Kwborder.c -o Kwborder.o
g++ -I/usr/lib/R/include  -I/usr/local/include-fpic  -O2 -c
PerfectStrauss.cc -o PerfectStrauss.o
make: g++: Command not found
make: *** [PerfectStrauss.o] Error 127
ERROR: compilation failed for package ‘spatstat’
* Removing ‘/home/richard/R/i686-pc-linux-gnu-library/2.9/spatstat’

The downloaded packages are in
‘/tmp/RtmpdcNYyo/downloaded_packages’
Warning messages:
1: In install.packages("spatstat") :
  installation of package 'deldir' had non-zero exit status
2: In install.packages("spatstat") :
  installation of package 'spatstat' had non-zero exit status

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



--
   O__   Peter Dalgaard Øster Farimagsgade 5, Entr.B
  c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
 (*) \(*) -- University of Copenhagen   Denmark  Ph:  (+45) 35327918
~~ - (p.dalga...@biostat.ku.dk)  FAX: (+45) 35327907

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


Re: [R] sqlSave()

2009-05-10 Thread Dieter Menne
Felipe Carrillo  yahoo.com> writes:

> I have created a MS Access table named 'PredictedValues' through the statement
below:
> myDB <- odbcConnectAccess("C:/Documents and Settings/Owner/Desktop/Rpond
> Farming.mdb",uid="admin",pwd="")
> sqlSave(myDB,PredictedValues,rownames=FALSE)
>   close(myDB) 
> 
> But if I run the code again with new values I get the message below:
> Error in sqlSave(myDB, PredictedValues, rownames = FALSE) : 
>   table ‘PredictedValues’ already exists
> and my new records don't get updated.
> 
> I was under the impression that 'sqlSave' would copy new data on top of the
existing one or if the table didn't
> exist it would create one with the new values. I tried 'sqlUpdate' but my
existing 'PredictedValues'
> didn't update. What am I doing wrong.

Either try safer = FALSE (great white shark) or append=TRUE (depending on what
you want).

sqlSave(safer = FALSE) uses the 'great white shark' method of testing tables
(bite it and see). The logic will unceremoniously DROP the table and create it
anew with its own choice of column types in its attempt to find a writable
solution. test = TRUE will not necessarily predict this behaviour. Attempting to
write indexed columns or writing to pseudo-columns are less obvious causes of
failed writes followed by a DROP. If your table structure is precious to you
back it up.

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


Re: [R] Traceback NA's to their first appearance

2009-05-10 Thread Duncan Murdoch

On 10/05/2009 1:01 AM, Benjamin Leutner wrote:

Hello,

I've written a relatively complex interconnected population model in R. 
When changing a certain parameter, the outputs end up containing NA's.
I would like to find out, in which step the model (in form of a loop) 
starts to produce NA's. Does anyone know how to achieve this?

Since it does not give an error I can't use the traceback() function.



There's no automatic way that I know of, but modifications to your code 
are relatively easy.  any(is.na(x)) will give TRUE if there are any NAs 
in x.  So adding the line


if (any(is.na(x))) browser()

to your code will drop to the browser if this is ever detected; you can 
then debug your code to work out what happened.


Duncan Murdoch

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


Re: [R] "Overloading" some non-dispatched S3 methods for new classes

2009-05-10 Thread Gavin Simpson
On Sat, 2009-05-09 at 20:50 +0200, Carlos J. Gil Bellosta wrote:
> Hello,
> 
> I am building a package that creates a new kind of object not unlike a
> dataframe. However, it is not an extension of a dataframe, as the data
> themselves reside elsewhere. It only contains "metadata".
> 
> I would like to be able to retrieve data from my objects such as the
> number of rows, the number of columns, the colnames, etc.
> 
> I --quite naively-- thought that ncol, nrow, colnames, etc. would be
> dispatched, so I would only need to create a, say, ncol.myclassname
> function so as to be able to invoke "ncol" directly and transparently.
> 
> However, it is not the case. The only alternative I can think about is
> to create decorated versions of ncol, nrow, etc. to avoid naming
> conflicts. But I would still prefer my package users to be able to use
> the undecorated function names.
> 
> Do I have a chance?

Yes, if I understand you correctly. nrow, ncol are not S3 generics. You
can make them generic in your package and assign to the default method
the function definition in the R base package. For example:

ncol <- function(object, ...) UseMethod("ncol")
ncol.default <- base::ncol

mat <- matrix(0, ncol = 10, nrow = 5)
class(mat) <- "myclass"

ncol.myclass <- function(object, ...) {
## as example, perversely, swap rows and cols
dim(object)[1]
## but your code would go in here
}
ncol(mat)

This is covered in the Writing R Extensions manual, section 7.1

http://cran.r-project.org/doc/manuals/R-exts.html#Adding-new-generics

Also, this is probably not the right list for such questions. R-Devel
would have been more appropriate.

HTH

G

> 
> Best regards,
> 
> Carlos J. Gil Bellosta
> http://www.datanalytics.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.
-- 
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
 Dr. Gavin Simpson [t] +44 (0)20 7679 0522
 ECRC, UCL Geography,  [f] +44 (0)20 7679 0565
 Pearson Building, [e] gavin.simpsonATNOSPAMucl.ac.uk
 Gower Street, London  [w] http://www.ucl.ac.uk/~ucfagls/
 UK. WC1E 6BT. [w] http://www.freshwaters.org.uk
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%

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


Re: [R] how to get design matrix?

2009-05-10 Thread Gavin Simpson
On Sat, 2009-05-09 at 19:29 -0700, linakpl wrote:
> If I was doing an ANOVA analysis how can I get the design matrix R used?

You can do ANOVA several different ways in R, and that is just the few
functions to do this that I am aware of. Showing us what function/code
you used would be really helpful and is requested by the posting guide.

If you fitted the model using lm, then use the model.matrix method:

dat <- data.frame(Y = rnorm(100), X = gl(4,25))
mod <- lm(Y ~ X, data = dat)
mod
summary(mod)
model.matrix(mod)

This also works if you used aov()

If you used something else, then you'll have to say what...

HTH

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

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


Re: [R] Citing R/Packages Question

2009-05-10 Thread Jim Lemon
Most common styles (e.g. APA, Harvard) include the date of access for an 
electronic reference. While this may be an artifact of history, both 
reviewers and editors are justified in asking authors to adhere to the 
style used in a particular journal. That said, I don't see why Nature or 
any other journal would refuse to include a reference that was properly 
formatted.


Jim

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


[R] ade4 and inertia of axis 2

2009-05-10 Thread John Poulsen


Hello,

Could someone please tell me how to find the estimate of inertia for the 
second axis of an RLQ analysis using ade4? Using the example from the 
ade4 package, I *suspect* that the inertia for the 2nd axis of the R 
table would be 4.139332 (see below results summary from rlq1).  However, 
the row is labeled "12" not "2" which suggests this is not correct.  In 
addition, this would mean that axis 2 of the RLQ would explain more than 
100% of the inertia of the separate ordination for dudimil 
((100*4.139)/dudimil$eig[2]=176%) which I believe is not supposed to occur.


Could anyone set me straight by either letting me know how to get the 
inertia for axis 2 or that my thinking is perhaps wrong?


Thanks,
John

data(aviurba)
 coa1 <- dudi.coa(aviurba$fau, scannf = FALSE, nf = 2)
 dudimil <- dudi.hillsmith(aviurba$mil, scannf = FALSE, nf = 2, row.w = coa1$lw)
 duditrait <- dudi.hillsmith(aviurba$traits, scannf = FALSE, nf = 2, row.w = 
coa1$cw)
 rlq1 <- rlq(dudimil, coa1, duditrait, scannf = FALSE, nf = 2)
 plot(rlq1)
 summary(rlq1)
 randtest.rlq(rlq1)


   summary(rlq1)


Eigenvalues decomposition:
  eig covar  sdR  sdQ  corr
1 0.4782826 0.6915798 1.558312 1.158357 0.3831293
2 0.1418508 0.3766308 1.308050 1.219367 0.2361331

Inertia & coinertia R:
  inertia  max ratio
1  2.428337 2.996911 0.8102800
12 4.139332 5.345110 0.7744148

Inertia & coinertia Q:
  inertia  max ratio
1  1.341791 2.603139 0.5154512
12 2.828648 4.202981 0.6730098

Correlation L:
 corr   max ratio
1 0.3831293 0.6435487 0.5953384
2 0.2361331 0.5220054 0.4523576

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


[R] ggplot2: recommended workaround for broken legend.position="top"

2009-05-10 Thread Zeljko Vrba
Searching the mail archives I found that using legend.position as in
p.ring.3 + opts(legend.position="top")

is a known bug.  I tried doing
p.ring.3 + opts(legend.position=c(0.8, 0.2))

which works, but the legend background is transparent, i.e. I see the
plot background through the legend.  Adding additional option

opts(legend.background=theme_rect(fill=TRUE,colour="white"))

fills the whole rectangle black(!), making text invisible, but leaves the shape
symbols visible.

So, how can I obtain a graph with legend positioned within the plot boundaries
(that's OK, I don't even mind manually positioning the legend), but on a white
background, i.e., so that the plot underneath is not visible?

Thanks.

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


Re: [R] ggplot2: recommended workaround for broken legend.position="top"

2009-05-10 Thread hadley wickham
On Sun, May 10, 2009 at 10:32 AM, Zeljko Vrba  wrote:
> Searching the mail archives I found that using legend.position as in
> p.ring.3 + opts(legend.position="top")
>
> is a known bug.  I tried doing
> p.ring.3 + opts(legend.position=c(0.8, 0.2))
>
> which works, but the legend background is transparent, i.e. I see the
> plot background through the legend.  Adding additional option
>
> opts(legend.background=theme_rect(fill=TRUE,colour="white"))
>
> fills the whole rectangle black(!), making text invisible, but leaves the 
> shape
> symbols visible.
>
> So, how can I obtain a graph with legend positioned within the plot boundaries
> (that's OK, I don't even mind manually positioning the legend), but on a white
> background, i.e., so that the plot underneath is not visible?

opts(legend.background=theme_rect(fill="white"))

Hadley

-- 
http://had.co.nz/

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


Re: [R] need help with chisq

2009-05-10 Thread Stephanie Kovalchik

JC,

If each row are the counts for a 2 x 2 contingency table - so for the  
ith contingency table you have counts for row 1 c(Y08[i],Z08[i]) and  
row 2 (Y09[i],Z09[i]) then you could use apply:


X <- cbind(vdata$Y08,vdata$X08-vdata$Y08,vdata$Y09,vdata$X09-vdata$Y98)

f.chisq <- function(x){
m <- matrix(x,2,2)
chisq.test(m)$p.value
}

apply(X,1,f.chisq)



Quoting JC :



I am very new to R. I have some data from a CVS stored in vdata with 4
columns labeled:
X08, Y08, X09, Y09.

I have created two new "columns" like so:

Z08 <- (vdata$X08-vdata$Y08)

Z09 <- (vdata$X09-vdata$Y09)

I would like to use chisq.test for each "row" and output the p-value
for each in a stored variable. I don't know how to do it. Can you
help?

so far I have done it for one row (but I want it done automatically
for all my data):

chidata=rbind(c(vdata$Y08[1],Z08[1]),c(vdata$Y09[1],Z09[1]))
results <- chisq.test(chidata)
results$p.value

I tried removing the [1] and the c() but that didn't work...  Any
ideas?

THANKS!

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



__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
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 general way to select a subset of matrix rows?

2009-05-10 Thread Stephanie Kovalchik

One simple adjustment is the following:

apply(matrix(x[rows,],nr=length(rows)),2,mean)


Quoting Peter Kharchenko :


Dear fellow R users,
I can't figure out how to do a simple thing properly: apply an
operation to matrix columns on a selected subset of rows. Things go
wrong when only one row is being selected. I am sure there's a way to
do this properly.

 Here's an example:
# define a 3-by-4 matrix x

x <- matrix(runif(12),ncol=4)
str(x)

num [1:3, 1:4] 0.568 0.217 0.309 0.859 0.651 ...

# calculate column means for selected rows

rows <- c(1,2)
apply(x[rows,],2,mean)

[1] 0.3923531 0.7552746 0.3661532 0.1069531
# now the same thing, but the "rows" vector is actually just one row

rows <- c(2)
apply(x[rows,],2,mean)

Error in apply(x[rows, ], 2, mean) : dim(X) must have a positive length

The problem is that while x[rows,] in the first case returned a matrix,
in the second case, when only one row was selected, it returned a
vector (and the apply obviously failed).  Is there a general way to
subset a matrix so it still returns a matrix even if it's one row?
Unfortunately doing as.matrix(x[rows,]) doesn't work either, as it
returns a transposed matrix in the case of a single row.

Is there a way to do this properly without writing out hideous if
statements accounting for single row exception?

thanks,
-peter.

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

2009-05-10 Thread Martin Batholdy

hi,


I am searching for a nice function which computes correlations out of  
a data.frame and adds asterix-signs after each correlation when they  
are significant...


...or a function which just show correlations greater than x in the  
output.





thanks!

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


[R] Function recommendation for this study...

2009-05-10 Thread Paul Heinrich Dietrich

Hi,
I'm not used to thinking along these lines, and wanted to ask your advice:

Suppose you have a sample of around 100, consisting of patients according to
doctors, in which patients and doctors are given a questionnaire with
categorical responses.  Each patient somehow has roughly 3 doctors, or 3
rows of data.  The goal is to assess by category of each question or DV the
agreement between the patient and 3 doctors.  For example, a question may be
asked about how well the treatment is understood by the patient, and the
patient answers with their perception, while the 3 doctors each answer with
their perception.

The person currently working on this has used a Wilcoxon Sign Rank test, and
asked what I thought.  Personally, I shy away from nonparametrics and prefer
parametric Bayesian methods, but of course am up for whatever is most
appropriate.  I was concerned about using multiple Wilcoxon tests, one for
each question, and wondering if there is a parametric method in R for
something like this, and a method which is multivariate?  Thanks for any
suggestions.
-- 
View this message in context: 
http://www.nabble.com/Function-recommendation-for-this-study...-tp23469646p23469646.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] ade4 and inertia of axis 2

2009-05-10 Thread John Poulsen

Hello,

Could someone please tell me how to find the estimate of inertia for the 
second axis of an RLQ analysis using ade4? Using the example from the 
ade4 package, I *suspect* that the inertia for the 2nd axis of the R 
table would be 4.139332 (see below results summary from rlq1).  However, 
the row is labeled "12" not "2" which suggests this is not correct.  In 
addition, this would mean that axis 2 of the RLQ would explain more than 
100% of the inertia of the separate ordination for dudimil 
((100*4.139)/dudimil$eig[2]=176%) which I believe is not supposed to occur.


Could anyone set me straight by either letting me know how to get the 
inertia for axis 2 or that my thinking is perhaps wrong?


Thanks,
John

data(aviurba)
  coa1 <- dudi.coa(aviurba$fau, scannf = FALSE, nf = 2)
  dudimil <- dudi.hillsmith(aviurba$mil, scannf = FALSE, nf = 2, row.w = 
coa1$lw)
  duditrait <- dudi.hillsmith(aviurba$traits, scannf = FALSE, nf = 2, row.w = 
coa1$cw)
  rlq1 <- rlq(dudimil, coa1, duditrait, scannf = FALSE, nf = 2)
  plot(rlq1)
  summary(rlq1)
  randtest.rlq(rlq1)


   summary(rlq1)


Eigenvalues decomposition:
   eig covar  sdR  sdQ  corr
1 0.4782826 0.6915798 1.558312 1.158357 0.3831293
2 0.1418508 0.3766308 1.308050 1.219367 0.2361331

Inertia & coinertia R:
   inertia  max ratio
1  2.428337 2.996911 0.8102800
12 4.139332 5.345110 0.7744148

Inertia & coinertia Q:
   inertia  max ratio
1  1.341791 2.603139 0.5154512
12 2.828648 4.202981 0.6730098

Correlation L:
  corr   max ratio
1 0.3831293 0.6435487 0.5953384
2 0.2361331 0.5220054 0.4523576

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


[R] Help with kalman-filterd betas using the dlm package

2009-05-10 Thread tom81

Hi all R gurus out there, 
Im a kind of newbie to kalman-filters after some research I have found that
the dlm package is the easiest to start with. So be patient if some of my
questions are too basic.

I would like to set up a beta estimation between an asset and a market index
using a kalman-filter. Much littarture says it gives superior estimates
compared to OLS estimates. So I would like to learn and to use the filter.

I would like to run two types of kalman-filters, one with using a
random-walk model (RW) and one with a stationary model, in other worlds the
transition equition either follow a RW or AR(1) model.

This is how I think it would be set up;

I will have my time-series Y,X, where Y is the response variable

this setup should give me a RW process if I have understood the example
correctly
mydlmModel = dlmModReg(X)  + dlmModPoly(order=1)

and then run on the dlm model
dlmFilter(Y,mydlmModel )

but setting up a AR(1) process is unclear, should I use dlmModPoly or the
dlmModARMA to set up the model.

And at last but not the least, how do I set up a proper build function to
use with dlmMLE to optimize the starting values.

Regards Tom
-- 
View this message in context: 
http://www.nabble.com/Help-with-kalman-filterd-betas-using-the-dlm-package-tp23473796p23473796.html
Sent from the R help mailing list archive at Nabble.com.

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


[R] Partial Derivatives in R

2009-05-10 Thread Paul Heinrich Dietrich

Quick question:

Which function do you use to calculate partial derivatives from a model
equation?

I've looked at deriv(), but think it gives derivatives, not partial
derivatives.  Of course my equation isn't this simple, but as an example,
I'm looking for something that let's you control whether it's a partial or
not, such as:

somefunction(y~a+bx, with respect to x, partial=TRUE)

Is there anything like this in R?
-- 
View this message in context: 
http://www.nabble.com/Partial-Derivatives-in-R-tp23470413p23470413.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] 2 dimension convolve

2009-05-10 Thread John Ros

Did you try the imgFFTConvolve function in the biOps package?



GreenBrower wrote:
> 
> Anybody know how to do a fast 2 dimension convolve in R? I used a
> traditional method, but it's very slow, even can not stand it!
> 

-- 
View this message in context: 
http://www.nabble.com/2-dimension-convolve-tp22440837p23470879.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] function for nice correlation output with significance symbols

2009-05-10 Thread Chuck Cleland
On 5/10/2009 3:26 PM, Martin Batholdy wrote:
> hi,
> 
> I am searching for a nice function which computes correlations out of a
> data.frame and adds asterix-signs after each correlation when they are
> significant...
> 
> ...or a function which just show correlations greater than x in the output.
> 
> thanks!

  This might be a starting point:

corstars <- function(x){
require(Hmisc)
x <- as.matrix(x)
R <- rcorr(x)$r
p <- rcorr(x)$P
mystars <- ifelse(p < .01, "**|", ifelse(p < .05, "* |", "  |"))
R <- format(round(cbind(rep(-1.111, ncol(x)), R), 3))[,-1]
Rnew <- matrix(paste(R, mystars, sep=""), ncol=ncol(x))
diag(Rnew) <- paste(diag(R), "  |", sep="")
rownames(Rnew) <- colnames(x)
colnames(Rnew) <- paste(colnames(x), "|", sep="")
Rnew <- as.data.frame(Rnew)
return(Rnew)
}

corstars(swiss[,1:4])

Fertility| Agriculture| Examination| Education|
Fertility 1.000  | 0.353* |-0.646**|  -0.664**|
Agriculture   0.353* | 1.000  |-0.687**|  -0.640**|
Examination  -0.646**|-0.687**| 1.000  |   0.698**|
Education-0.664**|-0.640**| 0.698**|   1.000  |

  At the very least, you could add a note that indicates what different
numbers of stars mean.

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

-- 
Chuck Cleland, Ph.D.
NDRI, Inc. (www.ndri.org)
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894

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


Re: [R] I'm offering $300 for someone who know R-programming to do the assignments for me.

2009-05-10 Thread Bob Doherty
On Sat, 09 May 2009 12:17:57 -0400, Carl Witthoft 
wrote:

>Sorry, but your professor offered me $500 NOT to do your assignments.
Actually all he wanted was your name so you could be expelled for
plagiarism!
-- 
Bob Doherty

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


Re: [R] Spatstat

2009-05-10 Thread Richard Chirgwin
Peter Dalgaard wrote:
> Richard Chirgwin wrote:
>> Hi all,
>>
>> I am trying to install Spatstat on OpenSUSE 11.1.
>> install.packages("spatstat", dependencies = TRUE)
>> fails on the basis of various compiler packages (full message below).
>>
>> I have gcc version 4.3.2, which should include gfortran and g++ -
>> so I'm not
>> sure what to do!
> 
> First check your assumptions...
> 
> If gfortran doesn't work on the command line, then you have an
> installation issue, nothing to do with R, much less spatstat.
> 
> Hint (on SUSE 11.0):
> 
> viggo:~/>rpm -qf `which gfortran`
> gcc-fortran-4.3-39.1
> 
> -pd

Thanks!

Richard

> 
> 
>>
>> Richard
>>
>> * Installing *source* package ‘deldir’ ...
>> ** libs
>> gfortran   -fpic  -O2 -c acchk.f -o acchk.o
>> make: gfortran: Command not found
>> make: *** [acchk.o] Error 127
>> ERROR: compilation failed for package ‘deldir’
>> * Removing ‘/home/richard/R/i686-pc-linux-gnu-library/2.9/deldir’
>> * Installing *source* package ‘spatstat’ ...
>> ** libs
>> gcc -std=gnu99 -I/usr/lib/R/include  -I/usr/local/include   
>> -fpic  -O2 -c
>> Kborder.c -o Kborder.o
>> gcc -std=gnu99 -I/usr/lib/R/include  -I/usr/local/include   
>> -fpic  -O2 -c
>> Kwborder.c -o Kwborder.o
>> g++ -I/usr/lib/R/include  -I/usr/local/include-fpic  -O2 -c
>> PerfectStrauss.cc -o PerfectStrauss.o
>> make: g++: Command not found
>> make: *** [PerfectStrauss.o] Error 127
>> ERROR: compilation failed for package ‘spatstat’
>> * Removing ‘/home/richard/R/i686-pc-linux-gnu-library/2.9/spatstat’
>>
>> The downloaded packages are in
>> ‘/tmp/RtmpdcNYyo/downloaded_packages’
>> Warning messages:
>> 1: In install.packages("spatstat") :
>>   installation of package 'deldir' had non-zero exit status
>> 2: In install.packages("spatstat") :
>>   installation of package 'spatstat' had non-zero exit status
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/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] Partial Derivatives in R

2009-05-10 Thread Ravi Varadhan
If you want `numerical' partial derivatives, check out:

  require(numDeriv)
  ?grad

Ravi.



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

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


- Original Message -
From: Paul Heinrich Dietrich 
Date: Sunday, May 10, 2009 4:45 pm
Subject: [R]  Partial Derivatives in R
To: r-help@r-project.org


>  Quick question:
>  
>  Which function do you use to calculate partial derivatives from a model
>  equation?
>  
>  I've looked at deriv(), but think it gives derivatives, not partial
>  derivatives.  Of course my equation isn't this simple, but as an example,
>  I'm looking for something that let's you control whether it's a 
> partial or
>  not, such as:
>  
>  somefunction(y~a+bx, with respect to x, partial=TRUE)
>  
>  Is there anything like this in R?
>  -- 
>  View this message in context: 
>  Sent from the R help mailing list archive at Nabble.com.
>  
>  __
>  R-help@r-project.org mailing list
>  
>  PLEASE do read the posting guide 
>  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] Partial Derivatives in R

2009-05-10 Thread David Winsemius


On May 10, 2009, at 10:05 AM, Paul Heinrich Dietrich wrote:



Quick question:

Which function do you use to calculate partial derivatives from a  
model

equation?

I've looked at deriv(), but think it gives derivatives, not partial
derivatives.


Your reading of the help page and the examples differs from mine. It  
says:
" It returns a call for computing the expr and its (partial)  
derivatives, simultaneously."


> dx2x <- deriv(~ x^2, "x") ; dx2x  # the example on the help page
expression({
.value <- x^2
.grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
.grad[, "x"] <- 2 * x
attr(.value, "gradient") <- .grad
.value
})
> dyx2.x <- deriv(~ y*x^2, "x") ; dyx2.x
expression({
.value <- y * x^2
.grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
.grad[, "x"] <- y * (2 * x)
attr(.value, "gradient") <- .grad
.value
})
> dy2x2.x <- deriv(~ y^2*x^2, "x") ; dy2x2.x
expression({
.expr1 <- y^2
.value <- .expr1 * x^2
.grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
.grad[, "x"] <- .expr1 * (2 * x)
attr(.value, "gradient") <- .grad
.value
})
> dy2x2.xy <- deriv(~ y^2*x^2, c("x","y")) ; dy2x2.xy
expression({
.expr1 <- y^2
.expr2 <- x^2
.value <- .expr1 * .expr2
.grad <- array(0, c(length(.value), 2L), list(NULL, c("x",
"y")))
.grad[, "x"] <- .expr1 * (2 * x)
.grad[, "y"] <- 2 * y * .expr2
attr(.value, "gradient") <- .grad
.value
})



Of course my equation isn't this simple, but as an example,
I'm looking for something that let's you control whether it's a  
partial or

not, such as:

somefunction(y~a+bx, with respect to x, partial=TRUE)




That appears to be precisely what your are offered with deriv,   
just not needing the partial=TRUE


deriv(
David Winsemius, MD
Heritage Laboratories
West Hartford, CT

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


Re: [R] Partial Derivatives in R

2009-05-10 Thread spencerg
 Have you considered genD{numDeriv}? 

 If this does not answer your question, I suggest you try the 
"RSiteSearch" package.  The following will open a list of options in a 
web browser, sorted by package most often found with your search term: 



library(RSiteSearch)
pd <- RSiteSearch.function('partial derivative')
pds <- RSiteSearch.function('partial derivatives')
attr(pd, 'hits') # 58
attr(pds, 'hits')# 52
summary(pd)
HTML(pd)
HTML(pds)

  
 The development version available via 
'install.packages("RSiteSearch", repos="http://R-Forge.R-project.org";)' 
also supports the following: 



pd. <- unionRSiteSearch(pd, pds)
attr(pd., 'hits')# 94
HTML(pd.)


 Hope this helps. 
 Spencer Graves


Paul Heinrich Dietrich wrote:

Quick question:

Which function do you use to calculate partial derivatives from a model
equation?

I've looked at deriv(), but think it gives derivatives, not partial
derivatives.  Of course my equation isn't this simple, but as an example,
I'm looking for something that let's you control whether it's a partial or
not, such as:

somefunction(y~a+bx, with respect to x, partial=TRUE)

Is there anything like this in R?



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

2009-05-10 Thread Andreas Christoffersen
Maybe pairs.panels(df,scale=T) from the psych library - se more here:
http://www.personality-project.org/r/html/00.psych-package.html

setting scale=T scales the cor coefficient according to their value. I
have seen an implementation with added asterix' but couldn't find it
right now.

On Sun, May 10, 2009 at 10:36 PM, Chuck Cleland  wrote:
> On 5/10/2009 3:26 PM, Martin Batholdy wrote:
>> hi,
>>
>> I am searching for a nice function which computes correlations out of a
>> data.frame and adds asterix-signs after each correlation when they are
>> significant...
>>
>> ...or a function which just show correlations greater than x in the output.
>>
>> thanks!
>

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


[R] Unintended loading of package:datasets

2009-05-10 Thread David A Vavra
The dataset package is being loaded apparently by one of the packages that I
am using. The loading of the datasets takes a long time and I would like to
eliminate it. I thought the datasets were effectively examples so don't
understand why they would be required at all.

1) How can I determine what is causing the datasets to be loaded?
2) How can I stop them from doing so?

I am using the following:

Rpart, grDevices, graphics, stats, utils, methods, base
There is also an environment named 'Autoloads'

TIA

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


[R] plot(survfit(fitCox)) graph shows one line - should show two

2009-05-10 Thread John Sorkin
R 2.8.1
Windows XP

I am trying to plot the results of a coxph using plot(survfit()). The plot 
should, I believe, show two lines one for survival in each of two treatment 
(Drug) groups, however my plot shows only one line. What am I doing wrong?

My code is reproduced below, my figure is attached to this EMail message.
John


> #Create simple survival object
> GVHDdata<-list(Time=GVHD$Time,Time30=(GVHD$Time)/30,
+ Age=GVHD$Age,Drug=GVHD$Drug,Died=GVHD$Died,
+ AgeGrp=cut(GVHD$Age,breaks=c(0,15,25,45)))
> 
> summary(GVHD$Drug)
MTX MXT+CSP 
 32  32 
> 
> 
> 
> fit0<-coxph(Surv(Time30,Died)~Drug,data=GVHDdata)
> summary(fit0)
Call:
coxph(formula = Surv(Time30, Died) ~ Drug, data = GVHDdata)

  n= 64 
 coef exp(coef) se(coef) z p
Drug[T.MXT+CSP] -1.15 0.3160.518 -2.23 0.026

exp(coef) exp(-coef) lower .95 upper .95
Drug[T.MXT+CSP] 0.316   3.16 0.115 0.871

Rsquare= 0.086   (max possible= 0.915 )
Likelihood ratio test= 5.75  on 1 df,   p=0.0165
Wald test= 4.96  on 1 df,   p=0.026
Score (logrank) test = 5.52  on 1 df,   p=0.0188

> 
> 
> plot(survfit(fit0))
 

John David Sorkin M.D., Ph.D.
Chief, Biostatistics and Informatics
University of Maryland School of Medicine Division of Gerontology
Baltimore VA Medical Center
10 North Greene Street
GRECC (BT/18/GR)
Baltimore, MD 21201-1524
(Phone) 410-605-7119
(Fax) 410-605-7913 (Please call phone number above prior to faxing)

Confidentiality Statement:
This email message, including any attachments, is for the sole use of the 
intended recipient(s) and may contain confidential and privileged information.  
Any unauthorized use, disclosure or distribution is prohibited.  If you are not 
the intended recipient, please contact the sender by reply email and destroy 
all copies of the original message. 
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Select the rows in a dataframe that matches a criteria in another dataframe

2009-05-10 Thread Cecilia Carmo
Hi everyone! Thank you for the help you have been given to 
me, and here I'm with another problem with my dataframes:
 
I have two dataframes (with much more observations), like 
these:

Dataframe1
Firm Year  cash
5004002002007 100
5004002002006 200
5004002002005 400
5004003002007 300
5004003002006 240
5004003002005 120
5004004002007 340
5004004002006 890
5004004002005 250

Dataframe 2
FirmAudited consolidate
500400200 yes   no
500400300 yes  yes
500400400 nono

I want to make another dataframe equal to the dataframe1, 
but just with the firms «audited», or with the firms 
«audited» and «consolidate». For example, with the audited 
and consolidated, the output would be just firm 500400300, 
like this:

Firm Year  cash
5004003002007 300
5004003002006 240
5004003002005 120

I’ve tried intersect () but it gives me just the number of 
the firm, and it is not what I want. What I want is a 
dataframe with all the information, but just the firms 
that match my criteria.


Could anyone help me with another function?

Thank you in advance,

Cecília Carmo (Portugal)

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


[R] How does one plot the -log(log(survival))

2009-05-10 Thread John Sorkin
R 2.8.1
Windows XP

How does one plot the -log(log(survival)) from a coxph? Survfit does not seem 
to be up to the task.
Thanks,
John

John David Sorkin M.D., Ph.D.
Chief, Biostatistics and Informatics
University of Maryland School of Medicine Division of Gerontology
Baltimore VA Medical Center
10 North Greene Street
GRECC (BT/18/GR)
Baltimore, MD 21201-1524
(Phone) 410-605-7119
(Fax) 410-605-7913 (Please call phone number above prior to faxing)

Confidentiality Statement:
This email message, including any attachments, is for th...{{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] plot(survfit(fitCox)) graph shows one line - should show two

2009-05-10 Thread Finak Greg
   I did not see an attached figure to compare against, perhaps it was removed 
by the server.  I suspect that:
plot(survfit(fit0,newdata=GVHDdata)
will produce the plot you're looking for. It will generate survival curves for 
the data in GVHDdata using the estimates from the fit0 coxph model. At least, 
that's what I gather from the documentation..

Cheers,
On 10/05/09 5:43 PM, "John Sorkin"  wrote:

R 2.8.1
Windows XP

I am trying to plot the results of a coxph using plot(survfit()). The plot 
should, I believe, show two lines one for survival in each of two treatment 
(Drug) groups, however my plot shows only one line. What am I doing wrong?

My code is reproduced below, my figure is attached to this EMail message.
John


> #Create simple survival object
> GVHDdata<-list(Time=GVHD$Time,Time30=(GVHD$Time)/30,
+ Age=GVHD$Age,Drug=GVHD$Drug,Died=GVHD$Died,
+ AgeGrp=cut(GVHD$Age,breaks=c(0,15,25,45)))
>
> summary(GVHD$Drug)
MTX MXT+CSP
 32  32
>
>
>
> fit0<-coxph(Surv(Time30,Died)~Drug,data=GVHDdata)
> summary(fit0)
Call:
coxph(formula = Surv(Time30, Died) ~ Drug, data = GVHDdata)

  n= 64
 coef exp(coef) se(coef) z p
Drug[T.MXT+CSP] -1.15 0.3160.518 -2.23 0.026

exp(coef) exp(-coef) lower .95 upper .95
Drug[T.MXT+CSP] 0.316   3.16 0.115 0.871

Rsquare= 0.086   (max possible= 0.915 )
Likelihood ratio test= 5.75  on 1 df,   p=0.0165
Wald test= 4.96  on 1 df,   p=0.026
Score (logrank) test = 5.52  on 1 df,   p=0.0188

>
>
> plot(survfit(fit0))


John David Sorkin M.D., Ph.D.
Chief, Biostatistics and Informatics
University of Maryland School of Medicine Division of Gerontology
Baltimore VA Medical Center
10 North Greene Street
GRECC (BT/18/GR)
Baltimore, MD 21201-1524
(Phone) 410-605-7119
(Fax) 410-605-7913 (Please call phone number above prior to faxing)

Confidentiality Statement:
This email message, including any attachments, is for th...{{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] plot(survfit(fitCox)) graph shows one line - should show two

2009-05-10 Thread Thomas Lumley

On Sun, 10 May 2009, John Sorkin wrote:


R 2.8.1
Windows XP

I am trying to plot the results of a coxph using plot(survfit()). The plot 
should, I believe, show two lines one for survival in each of two treatment 
(Drug) groups, however my plot shows only one line. What am I doing wrong?


Expecting two lines.  To get survival curves from a Cox model you need to 
specify the covariate values.  If you don't, you get the baseline survival (ie, 
at the mean covariate values).

Look at the example on ?survfit.coxph
  fit <- coxph(Surv(futime, fustat) ~ age, data = ovarian)
  plot(survfit(fit, newdata=data.frame(age=60)),
 xscale=365.25, xlab = "Years", ylab="Survival")

This plots the fitted curve for age 60. If you wanted multiple ages you can 
specify them, eg,
  plot(survfit(fit, newdata=data.frame(age=c(40,50,60))),
 xscale=365.25, xlab = "Years", ylab="Survival") 
gives three curves.


  -thomas



My code is reproduced below, my figure is attached to this EMail message.
John



#Create simple survival object
GVHDdata<-list(Time=GVHD$Time,Time30=(GVHD$Time)/30,

+ Age=GVHD$Age,Drug=GVHD$Drug,Died=GVHD$Died,
+ AgeGrp=cut(GVHD$Age,breaks=c(0,15,25,45)))


summary(GVHD$Drug)

   MTX MXT+CSP
32  32




fit0<-coxph(Surv(Time30,Died)~Drug,data=GVHDdata)
summary(fit0)

Call:
coxph(formula = Surv(Time30, Died) ~ Drug, data = GVHDdata)

 n= 64
coef exp(coef) se(coef) z p
Drug[T.MXT+CSP] -1.15 0.3160.518 -2.23 0.026

   exp(coef) exp(-coef) lower .95 upper .95
Drug[T.MXT+CSP] 0.316   3.16 0.115 0.871

Rsquare= 0.086   (max possible= 0.915 )
Likelihood ratio test= 5.75  on 1 df,   p=0.0165
Wald test= 4.96  on 1 df,   p=0.026
Score (logrank) test = 5.52  on 1 df,   p=0.0188




plot(survfit(fit0))



John David Sorkin M.D., Ph.D.
Chief, Biostatistics and Informatics
University of Maryland School of Medicine Division of Gerontology
Baltimore VA Medical Center
10 North Greene Street
GRECC (BT/18/GR)
Baltimore, MD 21201-1524
(Phone) 410-605-7119
(Fax) 410-605-7913 (Please call phone number above prior to faxing)

Confidentiality Statement:
This email message, including any attachments, is for ...{{dropped:7}}


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


Re: [R] How does one plot the -log(log(survival))

2009-05-10 Thread Thomas Lumley

On Sun, 10 May 2009, John Sorkin wrote:


R 2.8.1
Windows XP

How does one plot the -log(log(survival)) from a coxph? Survfit does not seem 
to be up to the task.


It depends on what you mean.

You can convert the plots to a -log(log(survival)) scale with the option 
fun="cloglog". In the example from my previous reply
   plot(survfit(fit, newdata=data.frame(age=c(40,50,60))),
  xscale=365.25, xlab = "Years", ylab="Survival",fun="cloglog")

This is a bit pointless: since these curves are from a Cox model, they are 
exactly parallel on a -log(log(survival)) scale.

People sometimes want to fit a model where they stratify on one covariate and 
plot the two baseline hazards. You can do this like:
  sfit<-coxph(Surv(futime, fustat) ~ age+strata(rx), data = ovarian)
  plot(survfit(sfit),fun="cloglog")

One motivation given for doing this is examining the proportional hazards 
assumption, but this is not a very good way of doing it. The eye is not very 
good at judging whether there is a constant vertical separation between curves, 
even before you take into account that the uncertainty varies dramatically 
along the curve.  I would use cox.zph() and plot.cox.zph(), which give valid 
formal tests and reasonable estimates of how the log hazard ratio varies with 
time.  Still, if you want to do log-log plots, the survival package is up to 
the task.

-thomas

Thomas Lumley   Assoc. Professor, Biostatistics
tlum...@u.washington.eduUniversity of Washington, Seattle

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


Re: [R] Unintended loading of package:datasets

2009-05-10 Thread Rolf Turner


On 11/05/2009, at 9:17 AM, David A Vavra wrote:

The dataset package is being loaded apparently by one of the  
packages that I
am using. The loading of the datasets takes a long time and I would  
like to
eliminate it. I thought the datasets were effectively examples so  
don't

understand why they would be required at all.

1) How can I determine what is causing the datasets to be loaded?
2) How can I stop them from doing so?

I am using the following:

Rpart, grDevices, graphics, stats, utils, methods, base
There is also an environment named 'Autoloads'

TIA


The datasets (note the ``s'') is a required R package which is  
*always* loaded

automatically --- and in my experience instantaneously.

I don't know about a dataset (singular) package.  There does not  
appear to be

one on CRAN.

There is some confusion in what you are doing.

cheers,

Rolf Turner

##
Attention:\ This e-mail message is privileged and confid...{{dropped:9}}

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


Re: [R] Select the rows in a dataframe that matches a criteria in another dataframe

2009-05-10 Thread David Winsemius


?subset
?"%in%"

(I have gotten tired of converting dataframes that are presented in a  
non-executable form, such as is supported by the dput function.   
So, ... you should read those help pages and take the obvious path to  
success.)


Something along the lines of:

subset(df1, Firm %in% df2[df2$Audited==yes, "Firm"] )  #untested

--
David Winsemius

On May 10, 2009, at 5:52 PM, Cecilia Carmo wrote:

Hi everyone! Thank you for the help you have been given to me, and  
here I'm with another problem with my dataframes:

I have two dataframes (with much more observations), like these:
Dataframe1
Firm Year  cash
5004002002007 100
5004002002006 200
5004002002005 400
5004003002007 300
5004003002006 240
5004003002005 120
5004004002007 340
5004004002006 890
5004004002005 250

Dataframe 2
FirmAudited consolidate
500400200 yes   no
500400300 yes  yes
500400400 nono

I want to make another dataframe equal to the dataframe1, but just  
with the firms «audited», or with the firms «audited» and  
«consolidate». For example, with the audited and consolidated, the  
output would be just firm 500400300, like this:

Firm Year  cash
5004003002007 300
5004003002006 240
5004003002005 120

I’ve tried intersect () but it gives me just the number of the firm,  
and it is not what I want. What I want is a dataframe with all the  
information, but just the firms that match my criteria.


David Winsemius, MD
Heritage Laboratories
West Hartford, CT

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


Re: [R] Unintended loading of package:datasets

2009-05-10 Thread David Winsemius


On May 10, 2009, at 7:33 PM, Rolf Turner wrote:



On 11/05/2009, at 9:17 AM, David A Vavra wrote:

The dataset package is being loaded apparently by one of the  
packages that I
am using. The loading of the datasets takes a long time and I would  
like to
eliminate it. I thought the datasets were effectively examples so  
don't

understand why they would be required at all.

1) How can I determine what is causing the datasets to be loaded?
2) How can I stop them from doing so?

I am using the following:

Rpart, grDevices, graphics, stats, utils, methods, base
There is also an environment named 'Autoloads'

TIA


The datasets (note the ``s'') is a required R package which is  
*always* loaded

automatically --- and in my experience instantaneously.

I don't know about a dataset (singular) package.  There does not  
appear to be

one on CRAN.

There is some confusion in what you are doing.



My guess is that she has created a large object in .Rdata that she  
does not remember. The way forward would be to use one of the several  
implementations of ls() that report size. Searching on"

 ls object size

...in the r-search pages may offer ideas.



Here is one that Jim Holtman offered earlier this year:

my.ls <-
 function (pos = 1, sorted = F)
 {
 .result <- sapply(ls(pos = pos, all.names = TRUE), function(..x)
object.size(eval(as.symbol(..x
 if (sorted) {
 .result <- rev(sort(.result))
 }
 .ls <- as.data.frame(rbind(as.matrix(.result), `**Total` =  
sum(.result)))

 names(.ls) <- "Size"
 .ls$Size <- formatC(.ls$Size, big.mark = ",", digits = 0,
 format = "f")
 .ls$Mode <- c(unlist(lapply(rownames(.ls)[-nrow(.ls)],
function(x) mode(eval(as.symbol(x),
 "---")
 .ls
 }

--

David Winsemius, MD
Heritage Laboratories
West Hartford, CT

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


Re: [R] Function recommendation for this study...

2009-05-10 Thread Murray Cooper

Paul,

I suggest looking up "observer agreement". The description of your study 
sounds like a classical

categorical observer agreement problem. I can't give
a reference off the top of my head, but if you get
stuck, e-mail me and I'll try and find a ref to get you started.

Murray M Cooper, Ph.D.
Richland Statistics
9800 N 24th St
Richland, MI, USA 49083
Mail: richs...@earthlink.net

- Original Message - 
From: "Paul Heinrich Dietrich" 

To: 
Sent: Sunday, May 10, 2009 8:25 AM
Subject: [R] Function recommendation for this study...




Hi,
I'm not used to thinking along these lines, and wanted to ask your advice:

Suppose you have a sample of around 100, consisting of patients according 
to

doctors, in which patients and doctors are given a questionnaire with
categorical responses.  Each patient somehow has roughly 3 doctors, or 3
rows of data.  The goal is to assess by category of each question or DV 
the
agreement between the patient and 3 doctors.  For example, a question may 
be

asked about how well the treatment is understood by the patient, and the
patient answers with their perception, while the 3 doctors each answer 
with

their perception.

The person currently working on this has used a Wilcoxon Sign Rank test, 
and
asked what I thought.  Personally, I shy away from nonparametrics and 
prefer

parametric Bayesian methods, but of course am up for whatever is most
appropriate.  I was concerned about using multiple Wilcoxon tests, one for
each question, and wondering if there is a parametric method in R for
something like this, and a method which is multivariate?  Thanks for any
suggestions.
--
View this message in context: 
http://www.nabble.com/Function-recommendation-for-this-study...-tp23469646p23469646.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] Unintended loading of package:datasets

2009-05-10 Thread David A Vavra
Thanks. Perhaps something else is going on. There is a large time period
(about 20 sec.) after the message about loading the package. More
investigation, I suppose.

Thanks again,
DAV


-Original Message-
From: Rolf Turner [mailto:r.tur...@auckland.ac.nz] 
Sent: Sunday, May 10, 2009 7:34 PM
To: David A Vavra
Cc: r-help@r-project.org
Subject: Re: [R] Unintended loading of package:datasets


On 11/05/2009, at 9:17 AM, David A Vavra wrote:

> The dataset package is being loaded apparently by one of the  
> packages that I
> am using. The loading of the datasets takes a long time and I would  
> like to
> eliminate it. I thought the datasets were effectively examples so  
> don't
> understand why they would be required at all.
>
> 1) How can I determine what is causing the datasets to be loaded?
> 2) How can I stop them from doing so?
>
> I am using the following:
>
> Rpart, grDevices, graphics, stats, utils, methods, base
> There is also an environment named 'Autoloads'
>
> TIA

The datasets (note the ``s'') is a required R package which is  
*always* loaded
automatically --- and in my experience instantaneously.

I don't know about a dataset (singular) package.  There does not  
appear to be
one on CRAN.

There is some confusion in what you are doing.

cheers,

Rolf Turner

##
Attention:\ This e-mail message is privileged and confid...{{dropped:13}}

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


Re: [R] Reading large files quickly

2009-05-10 Thread Rob Steele
At the moment I'm just reading the large file to see how fast it goes.
Eventually, if I can get the read time down, I'll write out a processed
version.  Thanks for suggesting scan(); I'll try it.

Rob

jim holtman wrote:
> Since you are reading it in chunks, I assume that you are writing out each
> segment as you read it in.  How are you writing it out to save it?  Is the
> time you are quoting both the reading and the writing?  If so, can you break
> down the differences in what these operations are taking?
> 
> How do you plan to use the data?  Is it all numeric?  Are you keeping it in
> a dataframe?  Have you considered using 'scan' to read in the data and to
> specify what the columns are?  If you would like some more help, the answer
> to these questions will help.
> 
> On Sat, May 9, 2009 at 10:09 PM, Rob Steele 
> wrote:
> 
>> Thanks guys, good suggestions.  To clarify, I'm running on a fast
>> multi-core server with 16 GB RAM under 64 bit CentOS 5 and R 2.8.1.
>> Paging shouldn't be an issue since I'm reading in chunks and not trying
>> to store the whole file in memory at once.  Thanks again.
>>
>> Rob Steele wrote:
>>> I'm finding that readLines() and read.fwf() take nearly two hours to
>>> work through a 3.5 GB file, even when reading in large (100 MB) chunks.
>>>  The unix command wc by contrast processes the same file in three
>>> minutes.  Is there a faster way to read files in R?
>>>
>>> Thanks!
>>  >
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
> 
> 
>

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


Re: [R] Unintended loading of package:datasets

2009-05-10 Thread Berwin A Turlach
G'day David,

On Sun, 10 May 2009 17:17:38 -0400
"David A Vavra"  wrote:

> The dataset package is being loaded apparently by one of the packages
> that I am using. The loading of the datasets takes a long time and I
> would like to eliminate it. I thought the datasets were effectively
> examples so don't understand why they would be required at all.
> 
> 1) How can I determine what is causing the datasets to be loaded?

help(options)
and then read the part on "defaultPackages"

> 2) How can I stop them from doing so?

Create an appropriate .Rprofile.  R-DownUnder had long time ago a
discussion on what people put into their .Rprofile,  and Bill
Venables' .Rprofile seem to contain the following snippet:

### puts four more packages on to the default
### search path.  I use them all the time
local({
  old <- getOption("defaultPackages")
  options(defaultPackages = c(old, "splines",
"lattice", "mgcv", "MASS"))
}) 
 
> I am using the following:
> 
> Rpart, grDevices, graphics, stats, utils, methods, base
> There is also an environment named 'Autoloads'

Rpart??  Or rpart??  I know of a package with the latter name but not
the former, and R is case sensitive.

So you may try:

local({
  options(defaultPackages=c("rpart", "grDevices", "graphics", "stats",
   "utils", "methods")
})

FWIW, note that for command XXX, help(XXX) will give you a help page
that has usually some example code at the end, that code can be run by
example(XXX).  Typically, such example code is using data sets from
the datasets package; so if you do not load it, the example() command
might not work anymore for some functions.  

Don't complain if your R installation doesn't work anymore if you mess
around with defaultPackages, in particular if you remove packages that
are usually in the default of defaultPackages. :)

And I agree with Rolf (Turner), it is hard to believe that the
datasets package would produce a noticeable delay on start-up; for me
it also loads instantaneously.  I guess that your problem is more
along David's (Winsemuis) guess. 

Cheers,

Berwin

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


Re: [R] Unintended loading of package:datasets

2009-05-10 Thread Berwin A Turlach
G'day David,

On Sun, 10 May 2009 21:35:30 -0400
"David A Vavra"  wrote:

> Thanks. Perhaps something else is going on. There is a large time
> period (about 20 sec.) after the message about loading the package.
> More investigation, I suppose.

What R version are you using?  I do not remember ever getting a message
that the package datasets is loaded, nor a message for any of the other
packages that are loaded by default.  (There once, long ago, may have
been such a message, but if so, I do not remember this anymore.)

Cheers,

Berwin

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


[R] Building US maps in R

2009-05-10 Thread Dirty D
Hi,

I'm trying to build some maps of the US by county that will have the 
following characteristics:

Feature/Map
Map 1
Map2
Both
Broken out by county
Yes
Yes
Yes
Heatmaps of US Census Data for income by county
Yes
No
Yes
Heatmaps of US Census Data for race by county (recoded as white and 
&non-white, with each county color coded based on the majority)
No
Yes
No
Polygon markers showing the coordinates of a list of restaurants that 
have geographical coordinates attached
Yes
Yes
Yes


I started learning R a year ago, but had to drop it for work reasons. Is 
there someone on the list who can help guide me through the process or 
refer me to /very, very/ detailed instructions and tutorials?

Thanks,

DD

[[alternative HTML version deleted]]

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


Re: [R] Partial Derivatives in R

2009-05-10 Thread Paul Heinrich Dietrich

Thank you for suggesting other functions, I will look into them.

When I read the deriv() function, it did mention partial, but I (being a
newbie) wasn't able to get partials for a simple MNL equation.  I'm sure I
did something wrong then, but here's what I tried the following and got
different answers, concluding prematurely maybe that deriv actually just
gives regular derivatives.  Thanks for looking:

### Variables for an observation

x01 <- rnorm(1,0,1)

x02 <- rnorm(1,0,1)

### Parameters for an observation

b00.1 <- rnorm(1,0,1)

b00.2 <- rnorm(1,0,1)

b00.3 <- 0

b01.1 <- rnorm(1,0,1)

b01.2 <- rnorm(1,0,1)

b01.3 <- 0

b02.1 <- rnorm(1,0,1)

b02.2 <- rnorm(1,0,1)

b02.3 <- 0

### Predicted Probabilities for an observation

phat1 <- 0.6

phat2 <- 0.3

phat3 <- 0.1



### Correct way to calculate a partial derivative
 for MNL
partial.b01.1 <- phat1 * (b01.1 - (b01.1*phat1+b01.2*phat2+b01.3*phat3))

partial.b01.2 <- phat2 * (b01.2 - (b01.1*phat1+b01.2*phat2+b01.3*phat3))

partial.b01.3 <- phat3 * (b01.3 - (b01.1*phat1+b01.2*phat2+b01.3*phat3))

partial.b01.1; partial.b01.2; partial.b01.3



partial.b02.1 <- phat1 * (b02.1 - (b01.1*phat1+b01.2*phat2+b01.3*phat3))

partial.b02.2 <- phat2 * (b02.2 - (b01.1*phat1+b01.2*phat2+b01.3*phat3))

partial.b02.3 <- phat3 * (b02.3 - (b01.1*phat1+b01.2*phat2+b01.3*phat3))

partial.b02.1; partial.b02.2; partial.b02.3



### Derivatives for MNL
 according to (my interpretation of) the deriv() function
dp1.dx <- deriv(phat1 ~ exp(b00.1+b01.1*x01+b02.1*x02) / 

(exp(b00.1+b01.1*x01+b02.1*x02)+exp(b00.2+b01.2*x01+b02.2*x02)+

exp(b00.3+b01.3*x01+b02.3*x02)), c("x01","x02"))

dp2.dx <- deriv(phat2 ~ exp(b00.2+b01.2*x01+b02.2*x02) / 

(exp(b00.1+b01.1*x01+b02.1*x02)+exp(b00.2+b01.2*x01+b02.2*x02)+

exp(b00.3+b01.3*x01+b02.3*x02)), c("x01","x02"))

dp3.dx <- deriv(phat3 ~ exp(b00.3+b01.3*x01+b02.3*x02) / 

(exp(b00.1+b01.1*x01+b02.1*x02)+exp(b00.2+b01.2*x01+b02.2*x02)+

exp(b00.3+b01.3*x01+b02.3*x02)), c("x01","x02"))

attr(eval(dp1.dx), "gradient")

attr(eval(dp2.dx), "gradient")

attr(eval(dp3.dx), "gradient")
-- 
View this message in context: 
http://www.nabble.com/Partial-Derivatives-in-R-tp23470413p23475411.html
Sent from the R help mailing list archive at Nabble.com.

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


Re: [R] Function recommendation for this study...

2009-05-10 Thread spencerg
 To look up "observer agreement", you might consider the 
"RSiteSearch" package.  The "RSiteSearch.function" looks only for 
matches in help pages of contributed packages. 



library(RSiteSearch)
oa <- RSiteSearch.function("observer agreement")
attr(oa, "hits") # 4 functions matching this term
HTML(oa) # to open the results as a table in a web browser. 



 An alternative search term for this is "interrater reliability". 



ir.r <- RSiteSearch.function("inter-rater reliability")
attr(ir.r, "hits") # 1
irr <- RSiteSearch.function("interrater reliability")
attr(irr, "hits") # 19
ir..r <- RSiteSearch.function("inter rater reliability")
attr(ir..r, "hits") # 1


 In particular, this identified a package "irr" for interrater 
reliability. 



 The development version of the "RSiteSearch" package on R-Forge 
includes a "unionRSiteSearch" function that allows one to combine all 
these searches.  You can get this version using 
'install.packages("RSiteSearch", 
repos="http://R-Forge.R-project.org";)';  if this gives you version 
1.0-0, please wait 24 hours and try again, because this version contains 
a bug that I just fixed.  The new version 1.0-1 should be available in 
24 hours.  With it, the following just worked for me: 



IRR <- unionRSiteSearch(oa, unionRSiteSearch(ir.r,
unionRSiteSearch(irr, ir..r) ) )
attr(IRR, "hits")
HTML(IRR)


 Hope this helps. 
 Spencer Graves


Murray Cooper wrote:

Paul,

I suggest looking up "observer agreement". The description of your 
study sounds like a classical

categorical observer agreement problem. I can't give
a reference off the top of my head, but if you get
stuck, e-mail me and I'll try and find a ref to get you started.

Murray M Cooper, Ph.D.
Richland Statistics
9800 N 24th St
Richland, MI, USA 49083
Mail: richs...@earthlink.net

- Original Message - From: "Paul Heinrich Dietrich" 


To: 
Sent: Sunday, May 10, 2009 8:25 AM
Subject: [R] Function recommendation for this study...




Hi,
I'm not used to thinking along these lines, and wanted to ask your 
advice:


Suppose you have a sample of around 100, consisting of patients 
according to

doctors, in which patients and doctors are given a questionnaire with
categorical responses.  Each patient somehow has roughly 3 doctors, or 3
rows of data.  The goal is to assess by category of each question or 
DV the
agreement between the patient and 3 doctors.  For example, a question 
may be

asked about how well the treatment is understood by the patient, and the
patient answers with their perception, while the 3 doctors each 
answer with

their perception.

The person currently working on this has used a Wilcoxon Sign Rank 
test, and
asked what I thought.  Personally, I shy away from nonparametrics and 
prefer

parametric Bayesian methods, but of course am up for whatever is most
appropriate.  I was concerned about using multiple Wilcoxon tests, 
one for

each question, and wondering if there is a parametric method in R for
something like this, and a method which is multivariate?  Thanks for any
suggestions.
--
View this message in context: 
http://www.nabble.com/Function-recommendation-for-this-study...-tp23469646p23469646.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.