[R] hashmap in R

[Qiang Li (Jonathan)]

Hi friends on R list,

Have people tried to implement a hashmap in R? What is the generic way to 
implement a lookup table in R?

Best,
Jonathan

 
===
Jonathan Qiang Li

===


  


[[elided Yahoo spam]]
[[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] Finding dependencies and clusters in live survey data with a mix of independent variable types

[[EMAIL PROTECTED]]

I have a set of live data about customer satisfaction and desires of a
live ecommerce site. There are only 311 survey responses. There were
approximately 154 questions. A large fraction of these questions were
questions with numerical answers (e.g. on a scale of 1 to 10 how
satisfied are you with our service,  how many months have you been a
customer for, how old are you, how many computers do you own). A
second large fraction of the questions had binary answers (e.g. do you
own an ipod,  do you think blogging will be more or less popular in 5
years time than it is now, do you use online video sites). The
remaining data were multinomial answers (e.g. from which of these
sources did you first find out about this site,  which of these most
closely describes the industry you are in).

I am mostly interested in finding subsets of customers for whom some
subset of survey answers best correlate with their answer to the
question On a scale of 1 to 10, how would you rate our overall
service?  I am also interested in identifying market segments of
like-minded individuals with similar interests and views and find out
what they, as a group most want from the service in the future.

I am aware of how to perform multiple linear regression using R but I
am not sure how to
1. handle the binary variables and multinomial variables as
independent variables
2. find a set of canonical independent variables which most closely
correlate in combination to the overall service rating data
3. find market segments among the data by looking for clusters of like
interests and views

Are any of the above suitable for analysis by R? If so, do there exist
example programs available which achieve similar things that I can
study as guides?

Thanks in advance for your contemplation.

Charlie

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

Re: [R] hashmap in R

[Dieter Menne]

Qiang Li (Jonathan jonathanli_1999 at yahoo.com writes:

 Have people tried to implement a hashmap in R? What is the generic way to
implement a lookup table in R?

By default, this is implicit, for example, when using named vectors.

http://finzi.psych.upenn.edu/R/Rhelp02a/archive/51507.html

and ?environment

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] SVD for a 500, 000* 500, 000 matrix (singular value decomposition)

[Uwe Ligges]

kayj wrote:

Hi,


I tried to run SVD on a  500,000* 500,000 matrix and i get a message that it
can  not allocate a vector of length 270 mb 



Well, you will obviously need  1Tera(!)bytes of RAM just in order to 
store the matrix (or is it some sparse one?). I wonder how you managed 
that issue.


Uwe Ligges



doe snayone know how to solve this problem? any ideas on other softwares
where I can do this? 
I appreciate your help

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] RMySQL and R 2.7.0 - Error in field.types$row.names : $ operator is invalid for atomic vectors

[Rainer M Krug]

Hi

when executin the following code, I get an $ operator is invalid for
atomic vectors. I understand the meaning of the error (and have seen
the warnings in earlier R versions). My question is, is there an
updated version of RMySQL which deals with it, or is my on;ly option
to switch to RODBC (which I would not like to do as I am using Linux)?
If I have to use RODBC, how can I easily adapt the follolwing code?

Thanks

Rainer

 m - dbDriver(MySQL)
 conn - dbConnect(m, group = test)
tblname - dummy
 dat - data.frame(n1=runif(10), n2=runif(10))
 dbWriteTable(
   conn,
   tblName,
   dat,
   overwrite=FALSE,
   append=TRUE
   )

Error in field.types$row.names : $ operator is invalid for atomic vectors


-- 
Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation
Biology, UCT), Dipl. Phys. (Germany)

Plant Conservation Unit
Department of Botany
University of Cape Town
Rondebosch 7701
South Africa

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

Re: [R] help with the unique function

[ravi]

Hi,
Thanks for the different suggestions - Jorge, Erik, Marc, Ted and Phil. Thanks 
to you, I have learnt quite a few, for me new, tricks and methods.
Thanks,
Ravi
- Original Message 
From: Phil Spector [EMAIL PROTECTED]
To: ravi [EMAIL PROTECTED]
Cc: r-help@r-project.org
Sent: Wednesday, 7 May, 2008 7:14:28 PM
Subject: Re: [R] help with the unique function

Another way to produce the data frame is

subset(as.data.frame(table(x)),Freq0)

                                        - Phil Spector
                    Statistical Computing Facility
                    Department of Statistics
                    UC Berkeley
                    [EMAIL PROTECTED]


On Wed, 7 May 2008, ravi wrote:

 Hi,
 The unique function is easy to understand and use. Beyond that, I want to get 
 also the frequency of repetition of each individual row in a data frame
 Let me explain with an example :
 x-data.frame(a=c(1,2,3,1,2),b=c(2,3,4,2,3),c=c(10,20,30,10,20))
 xu-unique(x)
 We have,
 x
   a b  c
 1 1 2 10
 2 2 3 20
 3 3 4 30
 4 1 2 10
 5 2 3 20
 xu
   a b  c
 1 1 2 10
 2 2 3 20
 3 3 4 30

 I want to get the following data frame :
   a b  c    Freq
 1 1 2 10    2
 2 2 3 20    2
 3 3 4 30    1
 That is, in addition to the unique rows, I want to get the frequency of 
 repetion of each individual row.
 I will appreciate all the help that I can get.
 Thank You,
 Ravi

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

[Megh Dal]

I got following error while I was using solve.QP() in my problem:
   
   Dmat = matrix(c(0.0001741, 0.0001280, 0.0001280, 0.0002570), nrow=2)
 dvec = t(c(0,0))
 Amat = matrix(c(-1,1,0,-1,0, 1,0,1,0,-1), nrow=5)
 bvec = c(-2, 1, 1, -5, -5)
 solve.QP(Dmat,dvec,Amat,bvec=bvec)
Error in solve.QP(Dmat, dvec, Amat, bvec = bvec) : 
  Amat and dvec are incompatible!
 

   
  Can anyone tell me where is the error in my definition?
   
  Regards,

   
-
[[elided Yahoo spam]]
[[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] cpower and censoring

[Daniel Brewer]

I would like to do some power estimations for a log-rank two sample test
and cpower seems to fit the bill.  I am getting confused though by the
man page and what the arguments actually mean. I am also not sure
whether cpower takes into account censoring or not.

Could anyone provide a simple example of how I would get the power for a
set control/non-control clinical trial where censoring occurs at an
estimated rate with an estimated drop out rate.

Quite confused about this.

-- 
**
Daniel Brewer, Ph.D.

Institute of Cancer Research
Molecular Carcinogenesis
Email: [EMAIL PROTECTED]
**

The Institute of Cancer Research: Royal Cancer Hospital, a charitable Company 
Limited by Guarantee, Registered in England under Company No. 534147 with its 
Registered Office at 123 Old Brompton Road, London SW7 3RP.

This e-mail message is confidential and for use by the a...{{dropped:2}}

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

[Daniele Amberti]

Here is my problem:
Autoregressive models are very interesting in forecasting consumptions (eg 
water, gas etc).

Generally time series of this type have a long history with relatively simple 
patterns and can be useful to add external regressors for calendar events 
(holydays, vacations etc).

arima() is a very powerful function but kalman filter is very slow (and I foun 
difficulties of estimation) while ar() is too simple but fast (but do not have 
a method for forecasting I think)

Is there something like arima() but entirely implemented in C and efficient 
like ar() ???
Is there something like step() for ARIMAX? It would be very useful for external 
regressors.

Try the code below (imagine daily data for some years):

x - rep(c(15,20,20,20,20,12,10), 5*52)
set.seed(1234)
x - x + rnorm(length(x))

#plot(as.ts(x[1:21]))

#slow
arima(x, c(1,0,1), list(order = c(2,0,0), period = 7))
arima(x, c(2,0,0), list(order = c(3,0,0), period = 7))
#slower
arima(x, c(2,0,1), list(order = c(3,0,0), period = 7))
# do not converge
arima(x, c(2,0,0), list(order = c(3,0,1), period = 7))

#fast but not enough sophisticated
ar(x)

Thanks in advance
Daniele




ORS Srl

Via Agostino Morando 1/3 12060 Roddi (Cn) - Italy
Tel. +39 0173 620211
Fax. +39 0173 620299 / +39 0173 433111
Web Site www.ors.it


Qualsiasi utilizzo non autorizzato del presente messaggio e dei suoi allegati 
è vietato e potrebbe costituire reato.
Se lei avesse ricevuto erroneamente questo messaggio, Le saremmo grati se 
provvedesse alla distruzione dello stesso
e degli eventuali allegati.
Opinioni, conclusioni o altre informazioni riportate nella e-mail, che non 
siano relative alle attività e/o
alla missione aziendale di O.R.S. Srl si intendono non attribuibili alla 
società stessa, né la impegnano in alcun modo.

[[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] Bar Charts

[hoogeebear]

Hi,

I created some bar charts. My first one is concerned with males, and my
second concerned with females.
Is there a way I can put the charts into one chart? There is 2 different
columns in each file. Here is my new file containing males and females:

gender,familar
Female,Yes
Female,Yes
Female,Yes
Female,Yes
Female,Yes
Female,No
Female,Yes
Female,Yes
Female,Yes
Female,Yes
Female,Yes
Female,Yes
Female,No
Female,Yes
Female,No
Female,Yes
Female,Yes
Female,Yes
Female,Yes
Female,Yes
Male,Yes
Male,Yes
Male,Yes
Male,No
Male,No
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,Yes
Male,No
Male,Yes
Male,Yes

Here is the code I use for creating a female chart:
library(plotrix)
library(prettyR)
female_familar -read.table(C://females.csv, sep=,, header=TRUE)
barp(rbind(rep(length(female_familar$gender),2),freq(female_familar$familar)[[1]]),
ylab=20 Females participated in the survey,
col=4:5,names.arg=c(FemalesNo(3),Females   Yes(17)))
legend(topright,c(Females,Familarity),fill=4:5)

Does the above need to change much to include males and females in the one
bar chart?

Hope to hear from someone soon.

Regards,

Jack.
-- 
View this message in context: 
http://www.nabble.com/Bar-Charts-tp17124678p17124678.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] (no subject)

[patricia garcía gonzález]

Hi everyone, 

Is there any function to standardize a matrix. For sure it must, but i can't 
find it. For standardize, i just mean, to make the mean as zero and standard 
deviation as one.It is also call z-score.
Thanks in advance



_
MSN Video. 

[[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] Aling elmentos into Windows with TK

[eduardo san miguel]

Hi,

There are many ways to do that.

An example:


require(tcltk)

tt - tktoplevel()
te - tkentry(tt)
tl - tklabel(tt)
tb - tkbutton(tt)

tkconfigure(tl, text = 'Enter text')
tkconfigure(tb, text = 'Show', command = function()
{cat(as.character(tkget(te)))})

tkgrid(tl, row = 0, column = 0, sticky = 'e')
tkgrid(te, row = 0, column = 1, sticky = 'e')
tkgrid(tb, row = 1, column = 0, columnspan = 2)

Sys.sleep(2)

tkgrid.remove(tl, te, tb)

tkgrid(tl, row = 0, column = 0, sticky = 'w')
tkgrid(te, row = 1, column = 0, sticky = 'w')
tkgrid(tb, row = 0, column = 1, rowspan = 2)


Anyway you should have a look \Tcl\doc  which is included in the standard R
dist.

Greetings

[[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] LaTeX in system()

[Christoph Heibl]

Dear list,

I want to run latex from an R script:

system(latex mysource.tex)

or:

	texi2dvi(mysource.tex, pdf = TRUE, clean = FALSE, quiet = 		TRUE,  
texi2dvi = latex)


but latex does not seem to be on the search path:

/bin/sh: line 1: latex: command not found.

Although 'printenv PATH' tells me that the usr/texbin is looked for  
executables:


/Library/Frameworks/Python.framework/Versions/Current/bin:/sw/bin:/sw/ 
sbin:/bin:/sbin:/usr/bin:/usr/sbin:/usr/local/bin:/usr/texbin:/usr/ 
X11R6/bin


Or am I wrong here?

In any case this is strange because if I call latex from the Terminal  
shell it runs without problems.  On the other hand texi2dvi from the  
tools package also does not work (for the same reason?)


I use R version 2.6.2 on Intel Mac OS X 10.4.11

Why is there a difference between the way the call to latex behaves  
directly in the shell or via the R system () command?


Thanks in advance for any advice!
Christoph












Christoph Heibl

Systematic Botany
Ludwig-Maximilians-Universität München
Menzinger Str. 67
D-80638 München
GERMANY

phone: +49-(0)89-17861-251
e-mail:[EMAIL PROTECTED]

http://www.botanik.biologie.uni-muenchen.de/botsyst/heibl/ch-home.html

SAVE PAPER - THINK BEFORE YOU PRINT

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

[G. Draisma]

Dear R-users,
I have output files having a variable number of tables
in the following format:
-
1
Pietje
I1 I2 Value
1  1  0.11
1  2  0.12
2  1  0.21

2
Jantje
I1 I2 I3 Value
1  1  1  0.111
3  3  3  0.333
...
-

Would there be an easy way
of turning this into (a list of) data.frames
with names Pietje, Jantje
and variables I1,I2,...Value?
(I1,I2 are string or categorical,
Value is double).

I used an sed script to extract the tables
from the file into separte file,
but would rather be able to process
the output files directly.

Thanks,
Gerrit.

--
Gerrit Draisma
Department of Public Health
Erasmus MC, University Medical Center Rotterdam
Room AE-103
P.O. Box 2040 3000 CA  Rotterdam The Netherlands
Phone: +31 10 7043124 Fax: +31 10 010-7038474
http://mgzlx4.erasmusmc.nl/pwp/?gdraisma

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

[R] How to compare a value in a matrix with its neighbours

[Simon Parker]

Hello.

I'm trying to compare a value in a 2-D matrix with its 4 neighbours.

I think I could use row() and col() and use x[row(x)] and x[row(x)-1] etc..., 
but I can't see how it would work.


Also, any way of having a matrix fill itself with its own coordinates would 
also work.

Any ideas please.


Simon Parker
Imperial College




   
-

A Smarter Email.
[[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] orders of the arima model

[Shubha Vishwanath Karanth]

Hi R,

 

Was just checking the ?ar , the autoregressive model. It's great that
R can give the order of the autoregressive model. Suppose if I had 2000
observations to fit an AR model. Then, if I am correct R builds 33+1
autoregressive models (10*log10(2000)=33)  and select the order at which
the aic's of the models were the least.

 

My question is whether the procedure ?arima  can throw the orders (p
and q) of the arima model?

 

Many Thanks,

Shubha

 

This e-mail may contain confidential and/or privileged i...{{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.

[R] acf function

[Irene Mantzouni]

Dear all, 

 

I have an annual time-series of population numbers and I would like to
estimate the auto-correlation. Can I use acf() function and judge
whether auto-correlation is significant by the plots? The acf array, eg:

 

Autocorrelations of series 'x$log.s.r', by lag

 

 0  1  2  3  4  5  6  7  8  9
10 11 12 

 1.000  0.031 -0.171 -0.451  0.021  0.070  0.238 -0.079 -0.046  0.006
0.188 -0.134 -0.016 

13 14 15 

-0.153  0.146  0.096

 

 gives the auto-correlation at lags 1, 2 Is that right? 

 

Thank you!


[[alternative HTML version deleted]]

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

[R] Re turning variable names with variables (in a function)

[Stropharia]

Dear R Users,

I have written a function that returns 4 variables. I would like to have the
variables returned with their variable names, is this possible?

- R Code -
mc.error - function(T, p=0.05){
se - sqrt((p)*(1-(p))/T) # standard error
upper - p+(1.96*se) # upper CI
lower - p-(1.96*se) # lower CI
cv - se/p # coefficient of variation
results - c(se, upper, lower, cv)
return(results) 
}
- R Code -

This returns (or something like this, depending on the value of T):

[1] 0.004998685 0.059797422 0.040202578 0.099973695

I would like:

[1] se=0.004998685 upper=0.059797422 lower=0.040202578 cv=0.099973695

Or even better:

se upper lowercv
0.004998685   0.059797422   0.040202578   0.099973695

Any help is much appreciated, thanks.

Steve

~~
Steven Worthington
Ph.D. Candidate
New York Consortium in
Evolutionary Primatology 
Department of Anthropology
New York University
25 Waverly Place
New York, NY 10003
U.S.A.
~~
-- 
View this message in context: 
http://www.nabble.com/Returning-variable-names-with-variables-%28in-a-function%29-tp17124874p17124874.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] help with updating to R2.7

[ravi]

Hi,
There is a lot of information at the link shown below. But it is not easy to 
know where to look. Also, not easy to interpret correctly the directions. 
I followed the following two methods, both of which failed for me. 
(1) I first installed R2.7. I then followed the directions in, What's the best 
way to upgrade ? in R for Windows FAQ. I copied the library folder from R2.6 
into the R2.7 folder.
I then tried to run the following command :
update.packages(checkBuilt=TRUE, ask=FALSE)
I don't remember the exact error message now, but it did not work. In the 
instructions, it is stated that one should transfer the installed packages to 
the library folder. Does this mean that I can just retain the R.css file in the 
R2.7 library folder and then fill it up with the old package folders from R2.6? 
Any specific and detailed suugestions are welcome here.
(2) I then tried the following method :
#---run in previous version
packages - installed.packages()[,Package]
save(packages, file=Rpackages)
#---run in new version
load(Rpackages)
for (p in setdiff(packages, installed.packages()[,Package]))
install.packages(p)
I copied the Rpackages file from R2.6 into R2.7 folder. But this did not work 
either. I know that it would be best if I reproduced the exact error messages, 
but I have tried so many different things now and have lost track of the exact 
error messages at each stage.
Looking for help,
Thanking you,
Ravi

- Original Message 
From: Gabor Grothendieck [EMAIL PROTECTED]
To: ravi [EMAIL PROTECTED]
Cc: r-help@r-project.org
Sent: Wednesday, 7 May, 2008 8:46:07 PM
Subject: Re: [R] help with updating to R2.7

There are a few different methods discussed in the README to the batchfiles
distribution.  Home page at:

http://batchfiles.googlecode.com


On Wed, May 7, 2008 at 12:51 PM, ravi [EMAIL PROTECTED] wrote:
 Hi,
 From R 2.6, I would like to update to R2.7. I would like to have some tips 
 on the recommended method of installing the latest versions of an entire 
 list of packages in R2.7 - i.e. all the packages that I have presently 
 installed in R2.6.
 I am hoping that there is an easier method than fetching the packages 
 individually as I did, to begin with, for R2.6.
 Additionally, I would like to install R2.7 on another pc which is not 
 connected to the internet. I have previously installed R2.6 on it by fetching 
 the zip files from another pc and installed them from the local zip files. I 
 would like to know if there is a method of transferring an entire 
 installation (base+selected packages) from one pc to another. Will there be a 
 special problem if the first is a win vista pc and the 2nd has win 2000?
 One problem that I have in continuing to use R2.6 is in getting compatible 
 versions of a new package that happens to attract my interest. For example, I 
 was interested in getting the arm package for R2.6.1. I could not find the 
 package in the R Binaries list on CRAN webpage. In the packages list, I did 
 find old versions in the arm archive. However, the files were in 
 .tar.gz format. How does one install from this format? It does not 
 appear to be as straightforward as in installing from the zip format 
 files.Recently, I had the same problem in installing a compatible version of 
 the norm package also. I am interested to know how I can tackle these 
 problems.
 I will appreciate all the help that I can get.
 Thanks,
 Ravi

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

[Prof Brian Ripley]

On Thu, 8 May 2008, Daniele Amberti wrote:


Here is my problem:
Autoregressive models are very interesting in forecasting consumptions (eg 
water, gas etc).

Generally time series of this type have a long history with relatively simple 
patterns and can be useful to add external regressors for calendar events 
(holydays, vacations etc).

arima() is a very powerful function but kalman filter is very slow (and I foun 
difficulties of estimation) while ar() is too simple but fast (but do not have 
a method for forecasting I think)

Is there something like arima() but entirely implemented in C and efficient 
like ar() ???


You mean, like arima0()?

I am not sure arima() is inefficient, rather that you are asking for the 
solution to a computationally difficult problem (which in your example is 
looking to estimate structure that is not there!).



Is there something like step() for ARIMAX? It would be very useful for external 
regressors.

Try the code below (imagine daily data for some years):

x - rep(c(15,20,20,20,20,12,10), 5*52)
set.seed(1234)
x - x + rnorm(length(x))

#plot(as.ts(x[1:21]))

#slow
arima(x, c(1,0,1), list(order = c(2,0,0), period = 7))
arima(x, c(2,0,0), list(order = c(3,0,0), period = 7))
#slower
arima(x, c(2,0,1), list(order = c(3,0,0), period = 7))
# do not converge
arima(x, c(2,0,0), list(order = c(3,0,1), period = 7))

#fast but not enough sophisticated
ar(x)

Thanks in advance
Daniele




ORS Srl

Via Agostino Morando 1/3 12060 Roddi (Cn) - Italy
Tel. +39 0173 620211
Fax. +39 0173 620299 / +39 0173 433111
Web Site www.ors.it


Qualsiasi utilizzo non autorizzato del presente messaggio e dei suoi allegati 
?? vietato e potrebbe costituire reato.
Se lei avesse ricevuto erroneamente questo messaggio, Le saremmo grati se 
provvedesse alla distruzione dello stesso
e degli eventuali allegati.
Opinioni, conclusioni o altre informazioni riportate nella e-mail, che non 
siano relative alle attivit?? e/o
alla missione aziendale di O.R.S. Srl si intendono non attribuibili alla 
societ?? stessa, n?? la impegnano in alcun modo.

[[alternative HTML version deleted]]




--
Brian D. Ripley,  [EMAIL PROTECTED]
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] Re turning variable names with variables (in a function)

[Shubha Vishwanath Karanth]

You can just go for a small change in the function as,

results - cbind(se, upper, lower, cv)
OR
results - data.frame(se, upper, lower, cv)

-S-
-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED]
On Behalf Of Stropharia
Sent: Thursday, May 08, 2008 5:51 PM
To: r-help@r-project.org
Subject: [R] Re turning variable names with variables (in a function)


Dear R Users,

I have written a function that returns 4 variables. I would like to have
the
variables returned with their variable names, is this possible?

- R Code -
mc.error - function(T, p=0.05){
se - sqrt((p)*(1-(p))/T) # standard error
upper - p+(1.96*se) # upper CI
lower - p-(1.96*se) # lower CI
cv - se/p # coefficient of variation
results - c(se, upper, lower, cv)
return(results) 
}
- R Code -

This returns (or something like this, depending on the value of T):

[1] 0.004998685 0.059797422 0.040202578 0.099973695

I would like:

[1] se=0.004998685 upper=0.059797422 lower=0.040202578 cv=0.099973695

Or even better:

se upper lowercv
0.004998685   0.059797422   0.040202578   0.099973695

Any help is much appreciated, thanks.

Steve

~~
Steven Worthington
Ph.D. Candidate
New York Consortium in
Evolutionary Primatology 
Department of Anthropology
New York University
25 Waverly Place
New York, NY 10003
U.S.A.
~~
-- 
View this message in context:
http://www.nabble.com/Returning-variable-names-with-variables-%28in-a-fu
nction%29-tp17124874p17124874.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.
This e-mail may contain confidential and/or privileged i...{{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.

Re: [R] Re turning variable names with variables (in a function)

[Ian Fiske]

You want to set the names attribute of your results vector.  You can do
this with the names() function (see ?names).  Specifically, you might use
something like this:

results - c(se, upper, lower, cv)
names(results) - c(se, upper, lower, cv)

Good luck,
Ian


Stropharia wrote:
 
 Dear R Users,
 
 I have written a function that returns 4 variables. I would like to have
 the variables returned with their variable names, is this possible?
 
 - R Code -
 mc.error - function(T, p=0.05){
   se - sqrt((p)*(1-(p))/T) # standard error
   upper - p+(1.96*se) # upper CI
   lower - p-(1.96*se) # lower CI
   cv - se/p # coefficient of variation
   results - c(se, upper, lower, cv)
   return(results) 
 }
 - R Code -
 
 This returns (or something like this, depending on the value of T):
 
 [1] 0.004998685 0.059797422 0.040202578 0.099973695
 
 I would like:
 
 [1] se=0.004998685 upper=0.059797422 lower=0.040202578 cv=0.099973695
 
 Or even better:
 
 se upper lowercv
 0.004998685   0.059797422   0.040202578   0.099973695
 
 Any help is much appreciated, thanks.
 
 Steve
 
 ~~
 Steven Worthington
 Ph.D. Candidate
 New York Consortium in
 Evolutionary Primatology 
 Department of Anthropology
 New York University
 25 Waverly Place
 New York, NY 10003
 U.S.A.
 ~~
 

-- 
View this message in context: 
http://www.nabble.com/Returning-variable-names-with-variables-%28in-a-function%29-tp17124874p17125080.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] help with updating to R2.7

[Duncan Murdoch]

On 5/8/2008 8:30 AM, ravi wrote:
I know that it would be best if I reproduced the exact error messages, but I have tried so many different things now 

 and have lost track of the exact error messages at each stage.

This is not a reasonable request.  Rather than trying those two methods 
one more time, and recording what goes wrong, you expect someone on the 
list to do it for you?  If the advice was in the FAQ, presumably it 
works for others, so the reason it didn't work for you is likely that 
you misinterpreted part of it, or your system has a unique problem, etc. 
 The *only* way someone could diagnose that would be to see the error 
messages.


You need to use some common sense when you're asking for help.

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] cpower and censoring

[Frank E Harrell Jr]

Daniel Brewer wrote:

I would like to do some power estimations for a log-rank two sample test
and cpower seems to fit the bill.  I am getting confused though by the
man page and what the arguments actually mean. I am also not sure
whether cpower takes into account censoring or not.

Could anyone provide a simple example of how I would get the power for a
set control/non-control clinical trial where censoring occurs at an
estimated rate with an estimated drop out rate.

Quite confused about this.



cpower handles censoring by specification of an accrual and follow-up 
periods and the event probability by reference time tref.  I believe the 
censoring distribution is assumed to be uniform, i.e., that accrual 
occurs at a uniform rate.


To have more control over the situation, the simulation-based spower 
function is one to look into.


--
Frank E Harrell Jr   Professor and Chair   School of Medicine
 Department of Biostatistics   Vanderbilt University

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

Re: [R] coxph - weights- robust SE

[Terry Therneau]

--- begin included message -
I am using coxph with weights to represent sampling fraction of subjects. 
Our simulation results show that the robust SE of beta systematically 
under-estimate the empirical SD of beta. 
Does anyone know how the robust SE are estimated in coxph using weights? 
Is there any analytical formula for the “weighted” robust SE?

--- end include ---

The formal reference is Binder, D.A. (1992).
Fitting Cox's proportional hazards models from survey data.
 Biometrika {\bf 79}, 139-47.
I am not aware of a systematic study of how well this estimator does, but there 
could well be papers in the literature.

Let D = matrix of dfbeta residuals = residual(fit, type='dfbeta', weighted=F), 
one row per observation and one column per variable.  Then for a simple 
weighted 
model

  var = D' W^2 D
  
where W = diagonal matrix of weights.  Note that the dfbetas are the 
derivatives 
of beta with respect to the weights; if all the weights are doubled the dfbetas 
go down by half.

Terry Therneau


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

Re: [R] ARIMA, AR, STEP - [ ] Message is from an unknown sender

[Daniele Amberti]

I agree with You that I'm trying to estimate structure that is not there, it's 
just an example, anyway I mean like function ar {stats} not arima0(). Function 
ar() is very simple with an automatic selection criterion (AIC etc). but miss 
moving average, integration, seasonality and predict method.

Library forecasting seems to do something this way but documentation is not so 
exhaustive.

Thanks

-Original Message-
From: Prof Brian Ripley [mailto:[EMAIL PROTECTED]
Sent: giovedì 8 maggio 2008 14.48
To: Daniele Amberti
Cc: r-help@r-project.org
Subject: Re: [R] ARIMA, AR, STEP - [ ] Message is from an unknown sender

On Thu, 8 May 2008, Daniele Amberti wrote:

 Here is my problem:
 Autoregressive models are very interesting in forecasting consumptions (eg 
 water, gas etc).

 Generally time series of this type have a long history with relatively simple 
 patterns and can be useful to add external regressors for calendar events 
 (holydays, vacations etc).

 arima() is a very powerful function but kalman filter is very slow (and I 
 foun difficulties of estimation) while ar() is too simple but fast (but do 
 not have a method for forecasting I think)

 Is there something like arima() but entirely implemented in C and efficient 
 like ar() ???

You mean, like arima0()?

I am not sure arima() is inefficient, rather that you are asking for the
solution to a computationally difficult problem (which in your example is
looking to estimate structure that is not there!).

 Is there something like step() for ARIMAX? It would be very useful for 
 external regressors.

 Try the code below (imagine daily data for some years):

 x - rep(c(15,20,20,20,20,12,10), 5*52)
 set.seed(1234)
 x - x + rnorm(length(x))

 #plot(as.ts(x[1:21]))

 #slow
 arima(x, c(1,0,1), list(order = c(2,0,0), period = 7))
 arima(x, c(2,0,0), list(order = c(3,0,0), period = 7))
 #slower
 arima(x, c(2,0,1), list(order = c(3,0,0), period = 7))
 # do not converge
 arima(x, c(2,0,0), list(order = c(3,0,1), period = 7))

 #fast but not enough sophisticated
 ar(x)

 Thanks in advance
 Daniele



 
 ORS Srl

 Via Agostino Morando 1/3 12060 Roddi (Cn) - Italy
 Tel. +39 0173 620211
 Fax. +39 0173 620299 / +39 0173 433111
 Web Site www.ors.it

 
 Qualsiasi utilizzo non autorizzato del presente messaggio e dei suoi allegati 
 ?? vietato e potrebbe costituire reato.
 Se lei avesse ricevuto erroneamente questo messaggio, Le saremmo grati se 
 provvedesse alla distruzione dello stesso
 e degli eventuali allegati.
 Opinioni, conclusioni o altre informazioni riportate nella e-mail, che non 
 siano relative alle attivit?? e/o
 alla missione aziendale di O.R.S. Srl si intendono non attribuibili alla 
 societ?? stessa, n?? la impegnano in alcun modo.

   [[alternative HTML version deleted]]



--
Brian D. Ripley,  [EMAIL PROTECTED]
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

ORS Srl

Via Agostino Morando 1/3 12060 Roddi (Cn) - Italy
Tel. +39 0173 620211
Fax. +39 0173 620299 / +39 0173 433111
Web Site www.ors.it


Qualsiasi utilizzo non autorizzato del presente messaggio e dei suoi allegati è 
vietato e potrebbe costituire reato.
Se lei avesse ricevuto erroneamente questo messaggio, Le saremmo grati se 
provvedesse alla distruzione dello stesso
e degli eventuali allegati.
Opinioni, conclusioni o altre informazioni riportate nella e-mail, che non 
siano relative alle attività e/o
alla missione aziendale di O.R.S. Srl si intendono non  attribuibili alla 
società stessa, né la impegnano in alcun modo.

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

Re: [R] (no subject)

[Jorge Ivan Velez]

Hi Paticia,

Perhaps ?scale is what you are looking for:

X=matrix(1:9,ncol=3)
[,1] [,2] [,3]
[1,]147
[2,]258
[3,]369

scale(X)
[,1] [,2] [,3]
[1,]   -1   -1   -1
[2,]000
[3,]111
attr(,scaled:center)
[1] 2 5 8
attr(,scaled:scale)
[1] 1 1 1


HTH,

Jorge


On Thu, May 8, 2008 at 6:48 AM, patricia garcía gonzález 
[EMAIL PROTECTED] wrote:


 Hi everyone,

 Is there any function to standardize a matrix. For sure it must, but i
 can't find it. For standardize, i just mean, to make the mean as zero and
 standard deviation as one.It is also call z-score.
 Thanks in advance



 _
 MSN Video.

[[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] acf function

[Irene Mantzouni]

Dear all, 
 
 
I have an annual time-series of population numbers and I would like to
estimate the auto-correlation. Can I use acf() function and judge
whether auto-correlation is significant by the plots? The acf array
produced by this functions gives the auto-correlation at lags 1, 2
Is that right? 
 
 
 
Thank you!

 


[[alternative HTML version deleted]]

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

Re: [R] SVD for a 500, 000* 500, 000 matrix (singular value decomposition)

[Aimin Yan]

Use blzpack, it could work it out.

Aimin

At 02:44 AM 5/8/2008, Uwe Ligges wrote:


kayj wrote:
Hi,

I tried to run SVD on a  500,000* 500,000 matrix and i get a message that it
can  not allocate a vector of length 270 mb


Well, you will obviously need  1Tera(!)bytes of RAM just in order to 
store the matrix (or is it some sparse one?). I wonder how you managed 
that issue.

Uwe Ligges


doe snayone know how to solve this problem? any ideas on other softwares
where I can do this? I appreciate your help
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.

Aimin Yan
Baker Center for Bioinformatics and Biological Statistics
109 Office  Lab Bldg
Iowa State University,Ames, IA 50010
E-mail: mailto:[EMAIL PROTECTED][EMAIL PROTECTED]
Phone: 515-520-0821
My 
server:http://omega.psi.iastate.edu/httphttp://omega.psi.iastate.edu/://omega.psi.iastate.edu/

[[alternative HTML version deleted]]

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

Re: [R] LaTeX in system()

[Peter Dalgaard]

Christoph Heibl wrote:
 Dear list,

 I want to run latex from an R script:

 system(latex mysource.tex)

 or:

 texi2dvi(mysource.tex, pdf = TRUE, clean = FALSE, quiet
 = TRUE, texi2dvi = latex)

 but latex does not seem to be on the search path:

 /bin/sh: line 1: latex: command not found.

 Although 'printenv PATH' tells me that the usr/texbin is looked for
 executables:

 /Library/Frameworks/Python.framework/Versions/Current/bin:/sw/bin:/sw/sbin:/bin:/sbin:/usr/bin:/usr/sbin:/usr/local/bin:/usr/texbin:/usr/X11R6/bin


 Or am I wrong here?

 In any case this is strange because if I call latex from the Terminal
 shell it runs without problems.  On the other hand texi2dvi from the
 tools package also does not work (for the same reason?)

 I use R version 2.6.2 on Intel Mac OS X 10.4.11

 Why is there a difference between the way the call to latex behaves
 directly in the shell or via the R system () command?


There can be several reasons, all depending crucially on your particular
setup. You might be better off on the R-sig-Mac list, but I'd start by
looking at system(echo $PATH) or system(printenv PATH). Presumably
that does not include /usr/texbin (if that is where latex resides). One
potential reason that it could be different from the command line lies
in the shell startup files -- which ones get executed depends on whether
or not it is a login shell. Another possibility is that R itself is
setting the PATH: try Sys.getenv(PATH).

 Thanks in advance for any advice!
 Christoph










 

 Christoph Heibl

 Systematic Botany
 Ludwig-Maximilians-Universität München
 Menzinger Str. 67
 D-80638 München
 GERMANY

 phone: +49-(0)89-17861-251
 e-mail:[EMAIL PROTECTED]

 http://www.botanik.biologie.uni-muenchen.de/botsyst/heibl/ch-home.html

 SAVE PAPER - THINK BEFORE YOU PRINT

 __
 R-help@r-project.org mailing list
 https://stat.ethz.ch/mailman/listinfo/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
~~ - ([EMAIL PROTECTED])  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] poisson regression with robust error variance ('eyestudy

[Ted Harding]

The below is an old thread:

On 02-Jun-04 10:52:29, Lutz Ph. Breitling wrote:
 Dear all,
 
 i am trying to redo the 'eyestudy' analysis presented on the site
 http://www.ats.ucla.edu/stat/stata/faq/relative_risk.htm
 with R (1.9.0), with special interest in the section on relative
 risk estimation by poisson regression with robust error variance.
 
 so i guess rlm is the function to use. but what is its equivalent
 to the glm's argument family to indicate 'poisson'? or am i
 somehow totally wrong and this is not applicable here?
 
 thx a lot-
 lutz
 =
 Lutz Ph. Breitling, CMd
 Unité des Recherches Médicale
 Hôpital Albert Schweitzer
 B.P. 118 Lambaréné (GABON)

It seems it may have led to a solution. However, I have tried
to trace through the thread in the R-help archives, and have
failed to find anything which lays out how a solution can be
formulated.[*]

I'm interested in the same question. Basically, if I fit a GLM
to Y=0/1 response data, to obtain relative risks, as in

  GLM - glm(Y ~ A + B + X + Z, family=poisson(link=log))

I can get the estimated RRs from

  RRs - exp(summary(GLM)$coef[,1])

but do not see how to implement confidence intervals based
on robust error variances using the output in GLM.

If there was a solution arrived at, can someone spell it
out for me, please?

With thanks,
Ted.

[*] PS: I did the R Site Help and Archive Search on eyestudy,
and got 18 hits, of which 6 were postings on that thread.

However, despite the fact that I requested the results to be
sorted by date, earliest on top, they seem to come out
in somewhat random order:

Sat Jun 5 14:47:26 2004 (Frank E. Harrell)
Wed Jun 9 13:06:04 2004 (Lutz Ph. Breitling)
Sat Jun 5 13:51:45 2004 (Lutz Ph. Breitling)
Wed Jun 2 18:24:54 2004 (Frank E. Harrell)
Wed Jun 2 12:44:45 2004 (Lutz Ph. Breitling) [the initial posting]
Wed Jun 2 16:36:43 2004 (Thomas Lumley)

This made it difficult to follow the thread, and indeed I wonder
if all the postings have been found. Does Namazu have problems
in this respect? (And, if I click on How to search, I get the
response You don't have permission to access /namazu.html on this
server.)

And, if I click on Contemporary messages sorted: ... by thread
I get the response:
The requested URL /R/Rhelp02a/archive/index.html was not found
 on this server.


E-Mail: (Ted Harding) [EMAIL PROTECTED]
Fax-to-email: +44 (0)870 094 0861
Date: 08-May-08   Time: 14:38:36
-- XFMail --

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

[R] Microseconds for a zoo object?

[Creighton, Sean]

Hello

I have a string which contains microseconds, can anyone help on
constructing this in to a time object, with the microseconds, that I can
take to a ZOO file?

Thanks
Sean



 UK[1,3] 
[1] 17:09:53.824
 UK[1,1]
[1] 2007-12-11 00:00:00
 mydates - paste( substr(UK[,1], 1, 10), UK[,3])  
 mydates[1]
[1] 2007-12-11 17:09:53.824
 is(mydates)
[1] character vector   
 dt1 - as.POSIXct(strptime(as.character(mydates),%Y-%m-%d
%H:%M:%S,tz=GMT))
 dt1[1]
[1] 2007-12-11 17:09:53 GMT
 is(dt1)
[1] POSIXt   oldClass POSIXct  POSIXlt 
 
 

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

[stephen sefick]

I have dissolved oxygen traces that are continuous  (fifteen minutes) for
two years (save for a couple of days, weeks, or minutes there depending on
the perogative of the river).  These traces are spaced out by river mile.  I
have figured out how to prepare data as to the sunspot example, but I can
not figure out how to get multiple traces into the prepSEIS function and
this is the warning that I get

Warning messages:
1: In notes[j] = paste(sep =  , GG[[ima]]$sta, GG[[ima]]$comp) :
  number of items to replace is not a multiple of replacement length
2: In stns[j] = GG[[ima]]$sta :
  number of items to replace is not a multiple of replacement length

when I plot this I get what I think is a signal composed of all of the
traces, and if I look at the results of prepSEIS it has compressed (I think)
all of the signals to one.

x - read.csv(testDO.csv, header=T)
x.ts - ts(x, deltat=4/60)
a -x.ts[,1]
b - x.ts[,2]
c - x.ts[,3]
d - x.ts[,4]
GG - prep1wig(wig = c(a, b, c, d), sta=c(RM215, RM202, RM198,
RM190, comp=DO, units=mg/L)) # I have also tried this without breaking
it up into individual data frames
xx - prepREIS(GG)

plot(x.ts) is what I want to go into PICK.GEN

Any help would be greatly appreciated- I am sure I am missing something.
thank you very much

Stephen Sefick



-- 
Let's not spend our time and resources thinking about things that are so
little or so large that all they really do for us is puff us up and make us
feel like gods. We are mammals, and have not exhausted the annoying little
problems of being mammals.

-K. Mullis

[[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] histogram() with Date class?

[Ola Caster]

Dear help list,

Is it possible to draw lattice histograms (i.e. use the histogram() function
and not the hist() function) with objects of class Date?

I've tried solutions like

histogram(~date, data=my.data, breaks=months)

but it doesn't seem to work.

Any suggestions are welcome.

Many thanks
Ola Caster

[[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] Microseconds for a zoo object?

[Dirk Eddelbuettel]

On 8 May 2008 at 14:58, Creighton, Sean  wrote:
| Hello
| 
| I have a string which contains microseconds, can anyone help on
| constructing this in to a time object, with the microseconds, that I can
| take to a ZOO file?

Easy, just read the docs:

i)  you need %0S instead of %S to parse sub-seconds components

ii) you need to tell R to print sub-second components.

[EMAIL PROTECTED]:~$ r -e 'posixt - strptime(2007-12-11 17:09:53.824123, 
%Y-%m-%d %H:%M:%OS); options(digits.secs=6); print(posixt)'
[1] 2007-12-11 17:09:53.824123

This used littler, but it's the same for R.  I set 
 options(digits.secs=6)
in my ~/.Rprofile/

Hth, D.

| 
| Thanks
| Sean
| 
| 
| 
|  UK[1,3] 
| [1] 17:09:53.824
|  UK[1,1]
| [1] 2007-12-11 00:00:00
|  mydates - paste( substr(UK[,1], 1, 10), UK[,3])  
|  mydates[1]
| [1] 2007-12-11 17:09:53.824
|  is(mydates)
| [1] character vector   
|  dt1 - as.POSIXct(strptime(as.character(mydates),%Y-%m-%d
| %H:%M:%S,tz=GMT))
|  dt1[1]
| [1] 2007-12-11 17:09:53 GMT
|  is(dt1)
| [1] POSIXt   oldClass POSIXct  POSIXlt 
|  
|  
| 
| __
| R-help@r-project.org mailing list
| https://stat.ethz.ch/mailman/listinfo/r-help
| PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
| and provide commented, minimal, self-contained, reproducible code.

-- 
Three out of two people have difficulties with fractions.

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

[kate]

I have a data with 236 observations. After plotting the histogram, I found that 
it looks like non-central t distribution. I would like to get MLE for mu and 
df. 

I found an example to find MLE for gamma distribution from fitting 
distributions with R:

library(stats4) ## loading package stats4
ll-function(lambda,alfa) {n-200
x-x.gam
-n*alfa*log(lambda)+n*log(gamma(alfa))-(alfa-
1)*sum(log(x))+lambda*sum(x)} ## -log-likelihood function
est-mle(minuslog=ll, start=list(lambda=2,alfa=1))

Is anyone how how to write down -log-likelihood function for noncentral t 
distribution?

Thanks a lot!!

Kate 
[[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] hashmap in R

[Earl F. Glynn]

Qiang Li (Jonathan) [EMAIL PROTECTED] wrote in message 
news:[EMAIL PROTECTED]...
 Hi friends on R list,

 Have people tried to implement a hashmap in R? What is the generic way to 
 implement a lookup table in R?

Does this help?

 x - rnorm(4)
 names(x) - c(a, b, c, d)
 x
 a  b  c  d
-1.4122868  1.3588267 -0.5499391 -0.3581889
 x[d]
 d
-0.3581889


efg
Earl F Glynn
Bioinformatics
Stowers Institute for Medical Research

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

Re: [R] help with updating to R2.7

[ravi]

Hi,
Ouch! That really hurt. But I get the point.
Here's what I did now. I copied all the package folders from R2.6, except for 
the R.css file, and copied them into the R2.7 folder. 
In the process, I overwrote the common files that came with the installation of 
R2.7.
Here's the output that I obtained in R2.7 :
 library()
Error in unique.default(ans) : 
  2 arguments passed to .Internal(unique) which requires 3
 update.packages(checkBuilt=TRUE, ask=FALSE)
Error: could not find function update.packages

My intention with my previous mail was to check and obtain confirmation on 
details like what I should do with the R.css file.
But I realise that I should have then restricted myself to that single 
question. 
That way, I would not be wasting the valuable time of the R list experts, from 
whom I have received a lot of help previously.
After reading the FAQ once again, I am guessing that the problem could be a 
case of the internet downloading functions failing.
As suggested in the FAQ, I tried the following commands :
 path_to_R\bin\Rgui.exe http_proxy=http://gannet/ http_proxy_user=ask
Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy
 path_to_R\bin\Rgui.exe http_proxy=http://user:[EMAIL PROTECTED]:80/
Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy
Should I run the above commands in R, or in the windows command window?

Thanking You,
Ravi

- Original Message 
From: Duncan Murdoch [EMAIL PROTECTED]
To: ravi [EMAIL PROTECTED]
Cc: r-help@r-project.org
Sent: Thursday, 8 May, 2008 2:57:07 PM
Subject: Re: [R] help with updating to R2.7

On 5/8/2008 8:30 AM, ravi wrote:
 I know that it would be best if I reproduced the exact error messages, but I 
 have tried so many different things now 
 and have lost track of the exact error messages at each stage.

This is not a reasonable request.  Rather than trying those two methods 
one more time, and recording what goes wrong, you expect someone on the 
list to do it for you?  If the advice was in the FAQ, presumably it 
works for others, so the reason it didn't work for you is likely that 
you misinterpreted part of it, or your system has a unique problem, etc. 
  The *only* way someone could diagnose that would be to see the error 
messages.

You need to use some common sense when you're asking for help.

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] help with updating to R2.7

[Duncan Murdoch]

On 5/8/2008 10:40 AM, ravi wrote:

Hi,
Ouch! That really hurt. But I get the point.
Here's what I did now. I copied all the package folders from R2.6, except for the R.css file, and copied them into the R2.7 folder. 


That's your problem.  You've hosed the 2.7 libraries.

You need to reinstall 2.7.  Only copy package folders from 2.6 if they 
don't exist in 2.7.  If they already exist there, they've already been 
updated.


The older packages may not work in 2.7, but once you run the 
update.packages() function, you should get the latest versions.


Duncan Murdoch


In the process, I overwrote the common files that came with the installation of 
R2.7.
Here's the output that I obtained in R2.7 :

library()
Error in unique.default(ans) : 
  2 arguments passed to .Internal(unique) which requires 3

update.packages(checkBuilt=TRUE, ask=FALSE)

Error: could not find function update.packages

My intention with my previous mail was to check and obtain confirmation on 
details like what I should do with the R.css file.
But I realise that I should have then restricted myself to that single question. 
That way, I would not be wasting the valuable time of the R list experts, from whom I have received a lot of help previously.

After reading the FAQ once again, I am guessing that the problem could be a 
case of the internet downloading functions failing.
As suggested in the FAQ, I tried the following commands :

path_to_R\bin\Rgui.exe http_proxy=http://gannet/ http_proxy_user=ask

Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy

path_to_R\bin\Rgui.exe http_proxy=http://user:[EMAIL PROTECTED]:80/

Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy
Should I run the above commands in R, or in the windows command window?

Thanking You,
Ravi

- Original Message 
From: Duncan Murdoch [EMAIL PROTECTED]
To: ravi [EMAIL PROTECTED]
Cc: r-help@r-project.org
Sent: Thursday, 8 May, 2008 2:57:07 PM
Subject: Re: [R] help with updating to R2.7

On 5/8/2008 8:30 AM, ravi wrote:
I know that it would be best if I reproduced the exact error messages, but I have tried so many different things now 
and have lost track of the exact error messages at each stage.


This is not a reasonable request.  Rather than trying those two methods 
one more time, and recording what goes wrong, you expect someone on the 
list to do it for you?  If the advice was in the FAQ, presumably it 
works for others, so the reason it didn't work for you is likely that 
you misinterpreted part of it, or your system has a unique problem, etc. 
  The *only* way someone could diagnose that would be to see the error 
messages.


You need to use some common sense when you're asking for help.

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] MLE for noncentral t distribution

[Duncan Murdoch]

On 5/8/2008 10:34 AM, kate wrote:
I have a data with 236 observations. After plotting the histogram, I found that it looks like non-central t distribution. I would like to get MLE for mu and df. 


I found an example to find MLE for gamma distribution from fitting distributions 
with R:

library(stats4) ## loading package stats4
ll-function(lambda,alfa) {n-200
x-x.gam
-n*alfa*log(lambda)+n*log(gamma(alfa))-(alfa-
1)*sum(log(x))+lambda*sum(x)} ## -log-likelihood function
est-mle(minuslog=ll, start=list(lambda=2,alfa=1))

Is anyone how how to write down -log-likelihood function for noncentral t 
distribution?



dt() has a non-centrality parameter and a log parameter, so it would 
simply be


ll - function(x, ncp, df) sum(dt(x, ncp=ncp, df=df, log=TRUE))

Make sure you convert mu into the ncp properly; the man page says how 
ncp is interpreted.


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] help with updating to R2.7

[Gabor Grothendieck]

That's not what movedir.bat and copydir.bat in batchfiles do.  They
will not overwrite any files.


On Thu, May 8, 2008 at 10:40 AM, ravi [EMAIL PROTECTED] wrote:
 Hi,
  Ouch! That really hurt. But I get the point.
  Here's what I did now. I copied all the package folders from R2.6, except 
 for the R.css file, and copied them into the R2.7 folder.
  In the process, I overwrote the common files that came with the installation 
 of R2.7.
  Here's the output that I obtained in R2.7 :
   library()
  Error in unique.default(ans) :
2 arguments passed to .Internal(unique) which requires 3

  update.packages(checkBuilt=TRUE, ask=FALSE)
  Error: could not find function update.packages

  My intention with my previous mail was to check and obtain confirmation on 
 details like what I should do with the R.css file.
  But I realise that I should have then restricted myself to that single 
 question.
  That way, I would not be wasting the valuable time of the R list experts, 
 from whom I have received a lot of help previously.
  After reading the FAQ once again, I am guessing that the problem could be a 
 case of the internet downloading functions failing.
  As suggested in the FAQ, I tried the following commands :
   path_to_R\bin\Rgui.exe http_proxy=http://gannet/ http_proxy_user=ask
  Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy
   path_to_R\bin\Rgui.exe http_proxy=http://user:[EMAIL PROTECTED]:80/
  Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy
  Should I run the above commands in R, or in the windows command window?

  Thanking You,
  Ravi


  - Original Message 
  From: Duncan Murdoch [EMAIL PROTECTED]
  To: ravi [EMAIL PROTECTED]
  Cc: r-help@r-project.org

 Sent: Thursday, 8 May, 2008 2:57:07 PM
  Subject: Re: [R] help with updating to R2.7



 On 5/8/2008 8:30 AM, ravi wrote:
   I know that it would be best if I reproduced the exact error messages, but 
 I have tried so many different things now
   and have lost track of the exact error messages at each stage.

  This is not a reasonable request.  Rather than trying those two methods
  one more time, and recording what goes wrong, you expect someone on the
  list to do it for you?  If the advice was in the FAQ, presumably it
  works for others, so the reason it didn't work for you is likely that
  you misinterpreted part of it, or your system has a unique problem, etc.
The *only* way someone could diagnose that would be to see the error
  messages.

  You need to use some common sense when you're asking for help.

  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.


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

Re: [R] MLE for noncentral t distribution

[Prof Brian Ripley]

On Thu, 8 May 2008, kate wrote:

I have a data with 236 observations. After plotting the histogram, I 
found that it looks like non-central t distribution. I would like to get 
MLE for mu and df.


So you mean 'non-central'?  See ?dt.


I found an example to find MLE for gamma distribution from fitting distributions 
with R:

library(stats4) ## loading package stats4
ll-function(lambda,alfa) {n-200
x-x.gam
-n*alfa*log(lambda)+n*log(gamma(alfa))-(alfa-
1)*sum(log(x))+lambda*sum(x)} ## -log-likelihood function
est-mle(minuslog=ll, start=list(lambda=2,alfa=1))

Is anyone how how to write down -log-likelihood function for noncentral t 
distribution?


Just use dt. E.g.


library(MASS)
?fitdistr


shows you a worked example for location, scale and df, but note the 
comments.  You could fit a non-central t, but it would be unusual to do 
so.




Thanks a lot!!

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



--
Brian D. Ripley,  [EMAIL PROTECTED]
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] MLE for noncentral t distribution

[kate]

Thanks for your quick reply.

I try the command as follows,

library(stats4) ## loading package stats4
ll - function(change, ncp, df) {-sum(dt(x, ncp=ncp, df=df, 
log=TRUE))}#-log-likelihood function

est-mle(minuslog=ll, start=list(ncp=-0.3,df=2))

But the warnings appears as follows,

invalid class mle object: invalid object for slot fullcoef in class 
mle: got class list, should be or extend class numeric


When I typed warnings(), I get

In dt(x, ncp = ncp, df = df, log = TRUE) :
 full precision was not achieved in 'pnt'

Does anyone know how to solve it?

Thanks,

Kate


- Original Message - 
From: Duncan Murdoch [EMAIL PROTECTED]

To: kate [EMAIL PROTECTED]
Cc: r-help@r-project.org
Sent: Thursday, May 08, 2008 9:46 AM
Subject: Re: [R] MLE for noncentral t distribution



On 5/8/2008 10:34 AM, kate wrote:
I have a data with 236 observations. After plotting the histogram, I 
found that it looks like non-central t distribution. I would like to get 
MLE for mu and df. I found an example to find MLE for gamma distribution 
from fitting distributions with R:


library(stats4) ## loading package stats4
ll-function(lambda,alfa) {n-200
x-x.gam
-n*alfa*log(lambda)+n*log(gamma(alfa))-(alfa-
1)*sum(log(x))+lambda*sum(x)} ## -log-likelihood function
est-mle(minuslog=ll, start=list(lambda=2,alfa=1))

Is anyone how how to write down -log-likelihood function for noncentral t 
distribution?



dt() has a non-centrality parameter and a log parameter, so it would 
simply be


ll - function(x, ncp, df) sum(dt(x, ncp=ncp, df=df, log=TRUE))

Make sure you convert mu into the ncp properly; the man page says how ncp 
is interpreted.


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] LaTeX in system()

[Erik Iverson]

Are you running R from the shell, or R.app?

I don't own or use a Mac, but I've seen something like this happen to 
people running R through ESS on some Emacs on a Mac.


Apologies for lack of precision here in terminology, but it had 
something to do with the PATH getting set through a shell initialization 
file (e.g., .cshrc, .bashrc).  Thus, starting in a shell, the path was 
set properly.  But the path was not set in the 'Mac GUI environment', so 
things started 'from there' (i.e., by clicking something), did not know 
about the path.


If it is a situation like the above, it really has nothing to do with R 
directly.  If that is the problem, there is probably some way to set 
environment variables so that the Mac GUI shell knows about them.


HTH,
Erik Iverson

Christoph Heibl wrote:

Dear list,

I want to run latex from an R script:

system(latex mysource.tex)

or:

texi2dvi(mysource.tex, pdf = TRUE, clean = FALSE, quiet = 
TRUE, texi2dvi = latex)


but latex does not seem to be on the search path:

/bin/sh: line 1: latex: command not found.

Although 'printenv PATH' tells me that the usr/texbin is looked for 
executables:


/Library/Frameworks/Python.framework/Versions/Current/bin:/sw/bin:/sw/sbin:/bin:/sbin:/usr/bin:/usr/sbin:/usr/local/bin:/usr/texbin:/usr/X11R6/bin 



Or am I wrong here?

In any case this is strange because if I call latex from the Terminal 
shell it runs without problems.  On the other hand texi2dvi from the 
tools package also does not work (for the same reason?)


I use R version 2.6.2 on Intel Mac OS X 10.4.11

Why is there a difference between the way the call to latex behaves 
directly in the shell or via the R system () command?


Thanks in advance for any advice!
Christoph












Christoph Heibl

Systematic Botany
Ludwig-Maximilians-Universität München
Menzinger Str. 67
D-80638 München
GERMANY

phone: +49-(0)89-17861-251
e-mail:[EMAIL PROTECTED]

http://www.botanik.biologie.uni-muenchen.de/botsyst/heibl/ch-home.html

SAVE PAPER - THINK BEFORE YOU PRINT

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] poisson regression with robust error variance ('eyestudy

[Paul Johnson]

On Thu, May 8, 2008 at 8:38 AM, Ted Harding
[EMAIL PROTECTED] wrote:
 The below is an old thread:

  On 02-Jun-04 10:52:29, Lutz Ph. Breitling wrote:
   Dear all,
  
   i am trying to redo the 'eyestudy' analysis presented on the site
   http://www.ats.ucla.edu/stat/stata/faq/relative_risk.htm
   with R (1.9.0), with special interest in the section on relative
   risk estimation by poisson regression with robust error variance.
  
   so i guess rlm is the function to use. but what is its equivalent
   to the glm's argument family to indicate 'poisson'? or am i
   somehow totally wrong and this is not applicable here?
There have been several questions about getting robust standard errors
in glm lately.  I went and read that UCLA website on the RR eye study
and the Zou article that uses a glm with robust standard errors.

I don't think rlm is the right way to go because that gives
different parameter estimates.  You need to estimate with glm and then
get standard errors that are adjusted for heteroskedasticity. Well,
you may wish to use rlm for other reasons, but to replicate that
eyestudy project, you need to take the ordinary usual estimates of the
b's and adjust the standard errors.

The Zou article does not give the mathematical formulae used to
estimate those robust errors, it rather gives a code snippit on using
SAS proc GENMOD to get those estimates.  Presumably, if we had access
to the SAS formula, we could easily get the calculations we need with
R. It is a little irksome to me that people think saying use SAS proc
GENMOD is informative.  Rather, we really need to know which TYPE of
robust standard error is being considered, since there are about 10
competing formulations.

I started checking various R packages for calculating sandwich
variance estimates.  There is a package sandwich for that purpose,
and in the car package there is a function hccm that can do so.
The hccm in car refuses to take glm objects, but the function vcovHC
in sandwich will do so. The discussion for sandwich's vcovHC
function refers to linear models, and lm objects are used in the
examples, but the vignette sandwich distributed with the package
states The HAC estimators are already available for generalized
linear models (fitted by glm) and robust regression (fitted by rlm in
package MASS). .

As a result, one can get a var/covar matrix from the routines in the
sandwich package, but I'm not entirely sure they are match the ones
SAS gives.  The vcovHC help page has a nice explanation of the
differences across sandwich estimators.  It says they are all variants
on

(X'X)^{-1} X' Omega X (X'X)^{-1}

With different stipulations about Omega.  If we knew the stipulation
about OMEGA used in the SAS routine, we could try it.

I tried to get the eyestudy data, but the link on the UCLA website is
no longer valid.

So I generated some phony 0-1 data:

y - as.numeric(rnorm(1000)  0)
x - rnorm(1000)
 glm1 - glm(y~x, family=binomial)
 hccm(glm1)
Error in hccm.lm(glm1) : requires unweighted lm
 vcovH
vcovHAC  vcovHC
 glm1 - glm(y~x, family=poisson(link=log))
 vcovHC(glm1)
  (Intercept) x
(Intercept)  1.017588e-03 -3.722456e-05
x   -3.722456e-05  9.492517e-04


The default type of estimation the method called HC3,  which is
recommended by authors of some Monte Carlo research studies.  vcovHC
calculates White's basic correction if you run:


 vcovHC(glm1, type=HC0)
  (Intercept) x
(Intercept)  1.013508e-03 -3.691381e-05
x   -3.691381e-05  9.417839e-04


Compare against the non-robust glm var/covar matrix.

 vcov(glm1)
  (Intercept) x
(Intercept)  0.0020152998 -0.778422
x   -0.778422  0.0018721903


In conclusion, use glm followed by vcovHC and I believe you will find
estimates like the ones provided by SAS or Stata.  But, without access
to their data, I can't be entirely sure.



-- 
Paul E. Johnson
Professor, Political Science
1541 Lilac Lane, Room 504
University of Kansas

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

[elnano]

I've basically solved the problem using the nls.lm function from the
minpack.lm (thanks Katharine) with some modifications for ignoring residuals
above a given percentile. This is to avoid the strong influence of points
which push my modeled vs. measured values away from the 1:1 line.
I based it on the example given for nls.lm. Here it is:

R   # soil respiration data
ST - ST [!is.na(R)] # soil temeprature data. Had to remove na to make
nls.lm work
SM - SM [!is.na(R)]# soil moisture data
R - R [!is.na(R)]
q - 0.95# quantile
   
p - c(a=-0.003, b=0.13, c=0.50, E=400)# model parameters with
estimated values

Rf - function(ST, SM, a, b, c, E)
{
expr -
expression((a*SM^I(2)+b*SM+c)*exp(E*((1/(283.15-227.13))-(1/(ST+273.15-227.13)
eval(expr)
}

optim.f - function(p, ST, SM, R, Rfcall)
{
res - R - do.call(Rfcall, c(list(ST = ST, SM = SM), as.list(p)))   #
the nls.lm example divides this by sqrt(R), I don't know why. I removed
that.
abserr - abs(res)
qnum - quantile(abserr, probs=q, na.rm=T)# calculate the q
quantile of the absolute errors
res[abserr  qnum]=0  # convert
residuals above qnum to  0
res
}

Rmodel- nls.lm(par = p, fn = optim.f, Rfcall = Rf, ST = ST, SM = SM, R
= R)

summary(Rmodel)

The only error I still get is when using lower q values. A q of around 0.95
or less (depending on the dataset) gives me completely wrong parameter
estimates resulting in negative predicted values. Maybe someone has a
suggestion here.

Fernando



Katharine Mullen wrote:
 
 The error message means that the gradient (first derivative of residual
 vector with respect to the parameter vector) is not possible to work with;
 calling the function qr on the gradient multiplied by the square root of
 the weight vector .swts (in your case all 1's) fails.
 
 If you want concrete advice it would be helpful to provide the commented,
 minimal, self-contained, reproducible code that the posting guide
 requests.  what are the values of ST04, SM08b, ch2no, and tower?
 
 General comments:  If your goal is to minimize sum( (observed -
 predicted)^2), the function you give nls to minimize (optim.fun in your
 case) should return the vector (observed - predicted), not the scalar sum(
 (observed - predicted)^2).  You may want to see the nls.lm function in
 package minpack.lm for another way to deal with the problem.
 
 On Wed, 7 May 2008, Fernando Moyano wrote:
 
 Greetings R users, maybe there is someone who can help
 me with this problem:

 I define a function optim.fun and want as output the
 sum of squared errors between predicted and measured
 values, as follows:

 optim.fun - function (ST04, SM08b, ch2no, a, b, d, E)
 {
 predR -
 (a*SM08b^I(2)+b*SM08b+d)*exp(E*((1/(283.15-227.13))-(1/(ST04+273.15-227.13
 abserr - abs(ch2no-predR)
 qnum - quantile(abserr, probs=0.95, na.rm=T)

 is.na(abserr) - (abserr  qnum)
 errsq - sum(abserr^2, na.rm=T)
 errsq
 }

 Then I want to optimize parameters a,b,d and E as to
 minimize the function output with the following:

 optim.model-nls(nulldat ~ optim.fun(ST04, SM08b,
 ch2no, a, b, d, E), data=tower,
 start=list(a=-0.003,b=0.13,d=0.50, E=400), na.action =
 na.exclude )

 were nulldat=0
 At this point I get the following error message:

 Error in qr.default(.swts * attr(rhs, gradient)) :
   NA/NaN/Inf in foreign function call (arg 1)
 Warning messages:
 1: In if (na.rm) x - x[!is.na(x)] else if
 (any(is.na(x))) stop(missing values and NaN's not
 allowed if 'na.rm' is FALSE) ... :
   the condition has length  1 and only the first
 element will be used
 (this warning is repeated 12 times)

 Question: what does the error mean? What am I doing
 wrong?
 Thanks a bunch.
 Nano
 Jen, Germany
 Max Planck for Biogeochemistry




  
 

 [[elided Yahoo spam]]

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

-- 
View this message in context: 
http://www.nabble.com/function-in-nls-argument-tp17108100p17127514.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 

Re: [R] help with updating to R2.7

[ravi]

Hi,
Thanks for the help. I have now solved the problem of installing the old 
packages in R2.7. A preliminary check showed that they seemed to work.
But I had the following problem with updating :
 update.packages(checkBuilt=TRUE, ask=FALSE)
--- Please select a CRAN mirror for use in this session ---
Warning: unable to access index for repository 
http://cran.uk.r-project.org/bin/windows/contrib/2.7
Warning: unable to access index for repository 
http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.7
Warning message:
In open.connection(con, r) :
  unable to connect to 'cran.r-project.org' on port 80.
This must be due to a firewall problem. Is there a way of solving this? 
Maybe, I should just repeat this process in another pc and internet connection 
which does not have this firewall problem. 
I can then copy over the resulting files into my present pc.
I have now downloaded the batchfiles folder with the movedir.bat and 
copydir.bat commands. 
This appears to offer a much better method of updating libraries.
But I am not very sure how I can get started with it.
It does not appear to be a regular package to be installed in R.
It seems that I need to run it from the windows command window. 
Can you please show me a small example of how I can get started with it.
Thanking You,
Ravi


- Original Message 
From: Gabor Grothendieck [EMAIL PROTECTED]
To: ravi [EMAIL PROTECTED]
Cc: Duncan Murdoch [EMAIL PROTECTED]; r-help@r-project.org
Sent: Thursday, 8 May, 2008 4:51:47 PM
Subject: Re: [R] help with updating to R2.7

That's not what movedir.bat and copydir.bat in batchfiles do.  They
will not overwrite any files.


On Thu, May 8, 2008 at 10:40 AM, ravi [EMAIL PROTECTED] wrote:
 Hi,
  Ouch! That really hurt. But I get the point.
  Here's what I did now. I copied all the package folders from R2.6, except 
for the R.css file, and copied them into the R2.7 folder.
  In the process, I overwrote the common files that came with the installation 
of R2.7.
  Here's the output that I obtained in R2.7 :
   library()
  Error in unique.default(ans) :
    2 arguments passed to .Internal(unique) which requires 3

  update.packages(checkBuilt=TRUE, ask=FALSE)
  Error: could not find function update.packages

  My intention with my previous mail was to check and obtain confirmation on 
details like what I should do with the R.css file.
  But I realise that I should have then restricted myself to that single 
question.
  That way, I would not be wasting the valuable time of the R list experts, 
from whom I have received a lot of help previously.
  After reading the FAQ once again, I am guessing that the problem could be a 
case of the internet downloading functions failing.
  As suggested in the FAQ, I tried the following commands :
   path_to_R\bin\Rgui.exe http_proxy=http://gannet/ http_proxy_user=ask
  Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy
   path_to_R\bin\Rgui.exe http_proxy=http://user:[EMAIL PROTECTED]:80/
  Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy
  Should I run the above commands in R, or in the windows command window?

  Thanking You,
  Ravi


  - Original Message 
  From: Duncan Murdoch [EMAIL PROTECTED]
  To: ravi [EMAIL PROTECTED]
  Cc: r-help@r-project.org

 Sent: Thursday, 8 May, 2008 2:57:07 PM
  Subject: Re: [R] help with updating to R2.7



 On 5/8/2008 8:30 AM, ravi wrote:
   I know that it would be best if I reproduced the exact error messages, but 
I have tried so many different things now
   and have lost track of the exact error messages at each stage.

  This is not a reasonable request.  Rather than trying those two methods
  one more time, and recording what goes wrong, you expect someone on the
  list to do it for you?  If the advice was in the FAQ, presumably it
  works for others, so the reason it didn't work for you is likely that
  you misinterpreted part of it, or your system has a unique problem, etc.
    The *only* way someone could diagnose that would be to see the error
  messages.

  You need to use some common sense when you're asking for help.

  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.



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

[Federico Calboli]

Hi everyone,

I am confused on how to specify some nesting and interaction terma with lme().

I have a dataset where some flies where selected for accessory gland size, made 
to mate in presence/absence of another male and the level of some protein 
measured. Now the complex stuff.


The selection has been replicated twice, so that the selection term has got two 
levels (large and small) with replicates large1/large2 and small1/small2.


A second complication comes from the fact the experiment has been replicated 
three times at three different months, in two blocks for each months.


In allthen I have

selection (fixed) with 2 levels
line%in%selection (random) with 4 levels (2 large, 2 small)
number of males (fixed) and continuous
replica (random) with 3 levels
block%in%replica (random) with 6 levels (2 for each month).

The easiest model ignores the nested random effects and uses just selection, 
males and replica and the relative interactions. The model


lme(y ~ selection * males, random = ~1|replica/selection/males, mydata)

gives and anova table

numDF denDF  F-value p-value
(Intercept) 1   228 870.5669  .0001
selection   1 2   0.2393  0.6731
males   1 4   0.0228  0.8874
selection:males 1 4   0.0941  0.7744

where the denDF for males and selection:males are wrong (it should be 2 in both 
cases). So my model is wrongly specified.


To sum it up, how would I model the 3 possible interactions between replica, 
selection and males?


and if I had the masochistic desire to add to the model line%in%selection and 
block%in%replica, how could I model that?


Chers,

Federico





--
Federico C. F. Calboli
Department of Epidemiology and Public Health
Imperial College, St Mary's Campus
Norfolk Place, London W2 1PG

Tel  +44 (0)20 7594 1602 Fax (+44) 020 7594 3193

f.calboli [.a.t] imperial.ac.uk
f.calboli [.a.t] gmail.com

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

Re: [R] LaTeX in system()

[Charilaos Skiadas]

This is exactly the problem: apps launched through the Finder do not  
go through the usual shell initialization process, where the PATH is  
typically set up. The two solutions would be to either use the full  
path to the command, or else start R.app from the Terminal, via the  
command:


open -a R

If you start R.app this way, it will get the same path as your shell.

Haris Skiadas
Department of Mathematics and Computer Science
Hanover College

On May 8, 2008, at 11:25 AM, Erik Iverson wrote:


Are you running R from the shell, or R.app?

I don't own or use a Mac, but I've seen something like this happen  
to people running R through ESS on some Emacs on a Mac.


Apologies for lack of precision here in terminology, but it had  
something to do with the PATH getting set through a shell  
initialization file (e.g., .cshrc, .bashrc).  Thus, starting in a  
shell, the path was set properly.  But the path was not set in the  
'Mac GUI environment', so things started 'from there' (i.e., by  
clicking something), did not know about the path.


If it is a situation like the above, it really has nothing to do  
with R directly.  If that is the problem, there is probably some  
way to set environment variables so that the Mac GUI shell knows  
about them.


HTH,
Erik Iverson

Christoph Heibl wrote:

Dear list,
I want to run latex from an R script:
system(latex mysource.tex)
or:
texi2dvi(mysource.tex, pdf = TRUE, clean = FALSE, quiet  
= TRUE, texi2dvi = latex)

but latex does not seem to be on the search path:
/bin/sh: line 1: latex: command not found.
Although 'printenv PATH' tells me that the usr/texbin is looked  
for executables:
/Library/Frameworks/Python.framework/Versions/Current/bin:/sw/bin:/ 
sw/sbin:/bin:/sbin:/usr/bin:/usr/sbin:/usr/local/bin:/usr/texbin:/ 
usr/X11R6/bin Or am I wrong here?
In any case this is strange because if I call latex from the  
Terminal shell it runs without problems.  On the other hand  
texi2dvi from the tools package also does not work (for the same  
reason?)

I use R version 2.6.2 on Intel Mac OS X 10.4.11
Why is there a difference between the way the call to latex  
behaves directly in the shell or via the R system () command?

Thanks in advance for any advice!
Christoph

Christoph Heibl
Systematic Botany
Ludwig-Maximilians-Universität München
Menzinger Str. 67
D-80638 München
GERMANY
phone: +49-(0)89-17861-251
e-mail:[EMAIL PROTECTED]
http://www.botanik.biologie.uni-muenchen.de/botsyst/heibl/ch- 
home.html

SAVE PAPER - THINK BEFORE YOU PRINT


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

Re: [R] acf function

[Christofer]

It should be ok

-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On
Behalf Of Irene Mantzouni
Sent: Thursday, May 08, 2008 5:45 PM
To: [EMAIL PROTECTED]
Subject: [R] acf function

Dear all, 

 

I have an annual time-series of population numbers and I would like to
estimate the auto-correlation. Can I use acf() function and judge
whether auto-correlation is significant by the plots? The acf array, eg:

 

Autocorrelations of series 'x$log.s.r', by lag

 

 0  1  2  3  4  5  6  7  8  9
10 11 12 

 1.000  0.031 -0.171 -0.451  0.021  0.070  0.238 -0.079 -0.046  0.006
0.188 -0.134 -0.016 

13 14 15 

-0.153  0.146  0.096

 

 gives the auto-correlation at lags 1, 2 Is that right? 

 

Thank you!


[[alternative HTML version deleted]]

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

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

[Tasneem Zaihra]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] stuck wid bootstrap (parametric)

[Tasneem Zaihra]

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

[Fernando Moyano]

I solved the problem arising from using certain quantile values simply by 
printing the iterations with the argument nprint. This gave me correct 
estimates. I have no idea why.

- Original Message 
From: elnano [EMAIL PROTECTED]
To: r-help@r-project.org
Sent: Thursday, May 8, 2008 5:43:31 PM
Subject: Re: [R] function in nls argument


I've basically solved the problem using the nls.lm function from the
minpack.lm (thanks Katharine) with some modifications for ignoring residuals
above a given percentile. This is to avoid the strong influence of points
which push my modeled vs. measured values away from the 1:1 line.
I based it on the example given for nls.lm. Here it is:

R   # soil respiration data
ST - ST [!is.na(R)] # soil temeprature data. Had to remove na to make
nls.lm work
SM - SM [!is.na(R)]# soil moisture data
R - R [!is.na(R)]
q - 0.95# quantile
  
p - c(a=-0.003, b=0.13, c=0.50, E=400)# model parameters with
estimated values

Rf - function(ST, SM, a, b, c, E)
{
expr -
expression((a*SM^I(2)+b*SM+c)*exp(E*((1/(283.15-227.13))-(1/(ST+273.15-227.13)
eval(expr)
}

optim.f - function(p, ST, SM, R, Rfcall)
{
res - R - do.call(Rfcall, c(list(ST = ST, SM = SM), as.list(p)))   #
the nls.lm example divides this by sqrt(R), I don't know why. I removed
that.
abserr - abs(res)
qnum - quantile(abserr, probs=q, na.rm=T)# calculate the q
quantile of the absolute errors
res[abserr  qnum]=0  # convert
residuals above qnum to  0
res
}

Rmodel- nls.lm(par = p, fn = optim.f, Rfcall = Rf, ST = ST, SM = SM, R
= R)

summary(Rmodel)

The only error I still get is when using lower q values. A q of around 0.95
or less (depending on the dataset) gives me completely wrong parameter
estimates resulting in negative predicted values. Maybe someone has a
suggestion here.

Fernando



Katharine Mullen wrote:
 
 The error message means that the gradient (first derivative of residual
 vector with respect to the parameter vector) is not possible to work with;
 calling the function qr on the gradient multiplied by the square root of
 the weight vector .swts (in your case all 1's) fails.
 
 If you want concrete advice it would be helpful to provide the commented,
 minimal, self-contained, reproducible code that the posting guide
 requests.  what are the values of ST04, SM08b, ch2no, and tower?
 
 General comments:  If your goal is to minimize sum( (observed -
 predicted)^2), the function you give nls to minimize (optim.fun in your
 case) should return the vector (observed - predicted), not the scalar sum(
 (observed - predicted)^2).  You may want to see the nls.lm function in
 package minpack.lm for another way to deal with the problem.
 
 On Wed, 7 May 2008, Fernando Moyano wrote:
 
 Greetings R users, maybe there is someone who can help
 me with this problem:

 I define a function optim.fun and want as output the
 sum of squared errors between predicted and measured
 values, as follows:

 optim.fun - function (ST04, SM08b, ch2no, a, b, d, E)
 {
 predR -
 (a*SM08b^I(2)+b*SM08b+d)*exp(E*((1/(283.15-227.13))-(1/(ST04+273.15-227.13
 abserr - abs(ch2no-predR)
 qnum - quantile(abserr, probs=0.95, na.rm=T)

is.na(abserr) - (abserr  qnum)
 errsq - sum(abserr^2, na.rm=T)
 errsq
 }

 Then I want to optimize parameters a,b,d and E as to
 minimize the function output with the following:

 optim.model-nls(nulldat ~ optim.fun(ST04, SM08b,
 ch2no, a, b, d, E), data=tower,
 start=list(a=-0.003,b=0.13,d=0.50, E=400), na.action =
 na.exclude )

 were nulldat=0
 At this point I get the following error message:

 Error in qr.default(.swts * attr(rhs, gradient)) :
   NA/NaN/Inf in foreign function call (arg 1)
 Warning messages:
 1: In if (na.rm) x - x[!is.na(x)] else if
 (any(is.na(x))) stop(missing values and NaN's not
 allowed if 'na.rm' is FALSE) ... :
   the condition has length  1 and only the first
 element will be used
 (this warning is repeated 12 times)

 Question: what does the error mean? What am I doing
 wrong?
 Thanks a bunch.
 Nano
 Jen, Germany
 Max Planck for Biogeochemistry




  
 

 [[elided Yahoo spam]]

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

-- 
View this message in context: 

[R] R strucchange question -- robust regression

[Kittler, Richard]

Is it possible to use some form of robust regression with the
breakpoints routine so that it is less sensitive to outliers? 

--Rich

Richard Kittler
Advanced Micro Devices, Inc.  
Sunnyvale, CA 

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

[statmobile]

Hey All,

I was wondering if I could solicit a little input on what I'm trying to 
do here.  I have a list of matrices, and I want to set their dimnames, 
but all I can come up with is this:


x - matrix(1:4,2)
y - matrix(5:8,2)

z - list(x,y)
nm - c(a,b)
nms - list(nm,nm)

z - lapply(z,function(x)dimnames(x)-nms)

As you can see, this doesn't quite work.  Could anybody help me, and 
bonus points if you can help me do this efficiently.  Please cc my email 
address, because I only get summaries of the list.


Thanks in advance!

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

[R] miniscule font size on R console output

[juanita choo]

Hi,

I am having a situation where I cannot change the output size of the R
console. I have played around with the font format menu but  the changes are
only reflected to the script that I type in but not to the output. Everytime
I run a script, I have to go back to font format to increase the output
script, which is currently showing up as small as the dust on my computer
screen. I have mac by the way.
Thanks for your suggestions.

juanita

[[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] Replicating Rows

[Tony Plate]

Another way is using straightforward indexing:

 x - cbind(trips=c(1,3,2), y=1:3, z=4:6)
 x
 trips y z
[1,] 1 1 4
[2,] 3 2 5
[3,] 2 3 6
 # generate row indices with the appropriate
 # number of repeats
 ii - rep(seq(len=nrow(x)), x[,1])
[1] 1 2 2 2 3 3
 # use these indices to select data rows
 x[ii, -1]
 y z
[1,] 1 4
[2,] 2 5
[3,] 2 5
[4,] 2 5
[5,] 3 6
[6,] 3 6


Jorge Ivan Velez wrote:

Hi Marion,

Try this:

set.seed(123)
mydf=data.frame(trips=rpois(10,5), matrix(rnorm(10*5),ncol=5))
mydf
sapply(mydf[,-1],rep,mydf[,1])


HTH,

Jorge


On Wed, May 7, 2008 at 11:41 PM, [EMAIL PROTECTED] wrote:



Hi,

 I have a data matrix in which there are 1000 rows by 30 columns. The
first
column of each row is a numeric indicating the number of trips taken to a
particular location with location attributes in the following column
entries for
that row.

I want to repeat each row based on the number of trips taken (as indicated
by
the number in the first column)...i.e., if 1,1 indicates 4 trips, I want
to
replicate row 1 four times, and do this for each entry of column 1.

I have played with rep command with little luck. Can anyone help? This is
probably very simple.

Thank you,

mw


Marion Wittmann, Ph.D. candidate
Environmental Science and Management
University of California
Santa Barbara, CA 93106-5131

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

[Achim Zeileis]

On Thu, 8 May 2008, Kittler, Richard wrote:

 Is it possible to use some form of robust regression with the
 breakpoints routine so that it is less sensitive to outliers?

Conceptually, it is possible to use the underlying dynamic programming
algorithm for other objective functions than the residual sum of squares.
However, the implementation of breakpoints() exploits the special
structure of OLS to speed up computations.

For package fxregime, I've written a more general (and hence even
slower) object-oriented implementation. Because it is still a bit instable
it's currently hidden in the namespace and essentially undocumented.

It could be used for what you want to do if:
  - Your robust regression is available in R in a function, fit() say,
with a formula interface
  fm - fit(formula, data)
  - Your robust regression has an objective function which is additive
in the observations and accessible in an extractor, objfun() say:
  objfun(fm)

For OLS these would be lm() and deviance():

  ## OLS-optimized interface
  bp1 - breakpoints(Nile ~ 1)

  ## object-oriented interface
  nile - data.frame(Nile = Nile)
  bp2 - fxregime:::gbreakpoints(Nile ~ 1, data = nile,
fit = lm, objfun = deviance)

  ## compare results
  summary(bp1)$breakpoints
  summary(bp2)$breakpoints

So if you've got something sensible as an alternative to lm() and
deviance(), you could plug that in. The downsides are: (1) This can be
terribly slow. (2) Information criteria are not automatically available
for selecting the number of breakpoints.

I hope that helps,
Z

 --Rich

 Richard Kittler
 Advanced Micro Devices, Inc.
 Sunnyvale, CA

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

[jim holtman]

I think what you want is this -- you have to return 'x' from the lapply:

x - matrix(1:4,2)
y - matrix(5:8,2)

z - list(x,y)
nm - c(a,b)
nms - list(nm,nm)

z - lapply(z,function(x){
dimnames(x)-nms
x
})

On Thu, May 8, 2008 at 1:33 PM, statmobile [EMAIL PROTECTED] wrote:
 Hey All,

 I was wondering if I could solicit a little input on what I'm trying to do
 here.  I have a list of matrices, and I want to set their dimnames, but all
 I can come up with is this:

 x - matrix(1:4,2)
 y - matrix(5:8,2)

 z - list(x,y)
 nm - c(a,b)
 nms - list(nm,nm)

 z - lapply(z,function(x)dimnames(x)-nms)

 As you can see, this doesn't quite work.  Could anybody help me, and bonus
 points if you can help me do this efficiently.  Please cc my email address,
 because I only get summaries of the list.

 Thanks in advance!

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




-- 
Jim Holtman
Cincinnati, OH
+1 513 646 9390

What is the problem you are trying to solve?

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

Re: [R] Setting matrix dimnames in a list

[Marc Schwartz]

on 05/08/2008 12:33 PM statmobile wrote:

Hey All,

I was wondering if I could solicit a little input on what I'm trying
to do here.  I have a list of matrices, and I want to set their
dimnames, but all I can come up with is this:

x - matrix(1:4,2) y - matrix(5:8,2)

z - list(x,y) nm - c(a,b) nms - list(nm,nm)

z - lapply(z,function(x)dimnames(x)-nms)

As you can see, this doesn't quite work.  Could anybody help me, and
 bonus points if you can help me do this efficiently.  Please cc my
email address, because I only get summaries of the list.

Thanks in advance!


The problem is that dimnames() does not return the matrix, so you need 
to explicitly return 'x':


 lapply(z, function(x) {dimnames(x) - nms; x})
[[1]]
  a b
a 1 3
b 2 4

[[2]]
  a b
a 5 7
b 6 8


HTH,

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] significance threshold in CCF

[E C]

Hi everyone,

When the CCF between two series of observations is plotted in R, a line 
indicating (presumably) the significance threshold appears across the plot. 
Does anyone know how this threshold is determined (it is different for each set 
of series) and how its value can be extracted from R? I've tried saving the CCF 
into an object and unclassing the object, but there's nothing there to indicate 
this.
Some sample code to show what I mean:
x - c(1,2,3,4,5,6,7,8,9)
y - c(3,4,5,6,7,8,9,10,11)
ccf(x, y, plot=T)
Thanks in advance!

_
If you like crossword puzzles, then you'll love Flexicon, a game which 
comb[[elided Hotmail spam]]

[[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] Reading multiple tables from file

[jim holtman]

Will this do it for you:

 x - readLines(textConnection(1
+ Pietje
+ I1 I2 Value
+ 1  1  0.11
+ 1  2  0.12
+ 2  1  0.21
+
+ 2
+ Jantje
+ I1 I2 I3 Value
+ 1  1  1  0.111
+ 3  3  3  0.333))
 closeAllConnections()
 start - grep(^[[:digit:]]+$, x)
 mark - vector('integer', length(x))
 mark[start] - 1
 # determine limits of each table
 mark - cumsum(mark)
 # split the data for reading
 df - lapply(split(x, mark), function(.data){
+ .input - read.table(textConnection(.data), skip=2, header=TRUE)
+ attr(.input, 'name') - .data[2]  # save the name
+ .input
+ })
 # rename the list
 names(df) - sapply(df, attr, 'name')
 df
$Pietje
  I1 I2 Value
1  1  1  0.11
2  1  2  0.12
3  2  1  0.21

$Jantje
  I1 I2 I3 Value
1  1  1  1 0.111
2  3  3  3 0.333



On Thu, May 8, 2008 at 7:40 AM, G. Draisma [EMAIL PROTECTED] wrote:
 Dear R-users,
 I have output files having a variable number of tables
 in the following format:
 -
 1
 Pietje
 I1 I2 Value
 1  1  0.11
 1  2  0.12
 2  1  0.21

 2
 Jantje
 I1 I2 I3 Value
 1  1  1  0.111
 3  3  3  0.333
 ...
 -

 Would there be an easy way
 of turning this into (a list of) data.frames
 with names Pietje, Jantje
 and variables I1,I2,...Value?
 (I1,I2 are string or categorical,
 Value is double).

 I used an sed script to extract the tables
 from the file into separte file,
 but would rather be able to process
 the output files directly.

 Thanks,
 Gerrit.

 --
 Gerrit Draisma
 Department of Public Health
 Erasmus MC, University Medical Center Rotterdam
 Room AE-103
 P.O. Box 2040 3000 CA  Rotterdam The Netherlands
 Phone: +31 10 7043124 Fax: +31 10 010-7038474
 http://mgzlx4.erasmusmc.nl/pwp/?gdraisma

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




-- 
Jim Holtman
Cincinnati, OH
+1 513 646 9390

What is the problem you are trying to solve?

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

Re: [R] help with the unique function

[Claudia Beleites]

Ravi,

if you have a large data.frame you might want to have a look at the count.rows 
function I collected from older threads and put into the wiki  
(http://wiki.r-project.org/rwiki/doku.php?id=tips:data-frames:count_and_extract_unique_rows)

WIth table I run into memory trouble - just as with aggregate.

Claudia

-- 
Claudia Beleites
Dipartimento dei Materiali e delle Risorse Naturali
Università degli Studi di Trieste
Via Alfonso Valerio 2
I-34127 Trieste

phone: +39 (0 40) 5 88-34 47
email: [EMAIL PROTECTED]

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

Re: [R] significance threshold in CCF

[Ted Harding]

On 08-May-08 18:23:27, E C wrote:
 Hi everyone,
 
 When the CCF between two series of observations is plotted in R, a line
 indicating (presumably) the significance threshold appears across the
 plot. Does anyone know how this threshold is determined (it is
 different for each set of series) and how its value can be extracted
 from R? I've tried saving the CCF into an object and unclassing the
 object, but there's nothing there to indicate this.
 Some sample code to show what I mean:
 x - c(1,2,3,4,5,6,7,8,9)
 y - c(3,4,5,6,7,8,9,10,11)
 ccf(x, y, plot=T)
 Thanks in advance!

Have a look at ?plot.acf (also used for plotting ccf objects)
Ted.


E-Mail: (Ted Harding) [EMAIL PROTECTED]
Fax-to-email: +44 (0)870 094 0861
Date: 08-May-08   Time: 19:54:20
-- XFMail --

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

[R] problem with caretNWS on linux

[Peter Tait]

Hi,
I am using caretNWS on a RHEL x86_64 system and I am getting an error 
message that is nearly identical to the one occuring in

http://www.r-project.org/nosvn/R.check/r-release-macosx-ix86/caretNWS-00check.txt

Error in socketConnection(serverHost, port = port, open = a+b, blocking = 
TRUE) :

 unable to open connection
Calls: system.time ... .local - tryCatch - tryCatchList - 
socketConnection

In addition: Warning message:
In socketConnection(serverHost, port = port, open = a+b, blocking = TRUE) 
:

 penguin:8765 cannot be opened

The software versions I am using are the following:
Python 2.3.4
twistd (the Twisted daemon) 1.3.0rc1
nwsserver-1.5.2
R version 2.6.2 (2008-02-08)
nws 1.6.3

Thank you for the help.
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.

Re: [R] help with updating to R2.7

[John C Frain]

When I upgrade in Windows from say 2.6.2 tp 2.7.0 I do the following

1. Install 2.7.0 in a new directory

2 Rename the library subdirectory in the new version from library to library2

3 Copy the library subdirectory in 2.6.2 to 2.7.0

4 Copy the contents of library2 to the transferred library.  This
re-installs the updated base packages.

5 Within R run update.packages(checkBuilt = TRUE) and you will be
presented with a list of packages and be given the option of updating
each package.

This procedure is easy to implement.  I have used it for several
updates and it has always worked well

Best Regards

John.

2008/5/8 Duncan Murdoch [EMAIL PROTECTED]:
 On 5/8/2008 10:40 AM, ravi wrote:

 Hi,
 Ouch! That really hurt. But I get the point.
 Here's what I did now. I copied all the package folders from R2.6, except
 for the R.css file, and copied them into the R2.7 folder.

 That's your problem.  You've hosed the 2.7 libraries.

 You need to reinstall 2.7.  Only copy package folders from 2.6 if they don't
 exist in 2.7.  If they already exist there, they've already been updated.

 The older packages may not work in 2.7, but once you run the
 update.packages() function, you should get the latest versions.

 Duncan Murdoch

 In the process, I overwrote the common files that came with the
 installation of R2.7.
 Here's the output that I obtained in R2.7 :

 library()

 Error in unique.default(ans) :  2 arguments passed to .Internal(unique)
 which requires 3

 update.packages(checkBuilt=TRUE, ask=FALSE)

 Error: could not find function update.packages

 My intention with my previous mail was to check and obtain confirmation on
 details like what I should do with the R.css file.
 But I realise that I should have then restricted myself to that single
 question. That way, I would not be wasting the valuable time of the R list
 experts, from whom I have received a lot of help previously.
 After reading the FAQ once again, I am guessing that the problem could be
 a case of the internet downloading functions failing.
 As suggested in the FAQ, I tried the following commands :

 path_to_R\bin\Rgui.exe http_proxy=http://gannet/ http_proxy_user=ask

 Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy

 path_to_R\bin\Rgui.exe http_proxy=http://user:[EMAIL PROTECTED]:80/

 Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy
 Should I run the above commands in R, or in the windows command window?

 Thanking You,
 Ravi

 - Original Message 
 From: Duncan Murdoch [EMAIL PROTECTED]
 To: ravi [EMAIL PROTECTED]
 Cc: r-help@r-project.org
 Sent: Thursday, 8 May, 2008 2:57:07 PM
 Subject: Re: [R] help with updating to R2.7

 On 5/8/2008 8:30 AM, ravi wrote:

 I know that it would be best if I reproduced the exact error messages,
 but I have tried so many different things now and have lost track of the
 exact error messages at each stage.

 This is not a reasonable request.  Rather than trying those two methods
 one more time, and recording what goes wrong, you expect someone on the list
 to do it for you?  If the advice was in the FAQ, presumably it works for
 others, so the reason it didn't work for you is likely that you
 misinterpreted part of it, or your system has a unique problem, etc.  The
 *only* way someone could diagnose that would be to see the error messages.

 You need to use some common sense when you're asking for help.

 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.




-- 
John C Frain
Trinity College Dublin
Dublin 2
Ireland
www.tcd.ie/Economics/staff/frainj/home.html
mailto:[EMAIL PROTECTED]
mailto:[EMAIL PROTECTED]

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] applying cor.test to a (m, n) matrix

[Monica Pisica]

Hi everybody,

I would like to apply cor.test to a matrix with m rows and n columns and get 
the results in a list of matrices , one matrix for p.val, one for the 
statistic, one for the correlation and 2 for upper and lower confidence 
intervals, something analog with cor() applied to a matrix.

I have done my own function to get a matrix of p.values and i suppose i can 
build similar functions for all the others. But i have used for loops and i 
wonder if there is any way to actually use one of the functions from the 
apply family to do this in a quicker way.

Here is my little function:

cor.pval - function(x, method = c(pearson, kendal, spearman), digit=8) {
n - dim(x)[2]
pval - matrix(paste(rep(c, n*n), seq(1,n*n), sep = ), n, n, byrow = T)
for (i in 1:n) {
for (j in 1:n){
pval[i, j] - cor.test(x[,i], x[,j], method = method)$p.value
  }
   }
pval - matrix(round(as.numeric(pval),digit), n, n, byrow = T)
rownames(pval) - colnames(x)
colnames(pval) - colnames(x)
return(pval)
}

Thanks for any input,

Monica







_


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

[Katharine Mullen]

To make your example reproducible you have to provide the data somehow;
I am pretty sure nprint doesn't effect the solution, but if it does this
would be a bug and I would appreciate a reproducible report.

The example in nls.lm is a little complicated in order to show how to use
an analytical expression for the gradient (possible to provide in the
argument jac); since you don't need this, note that your residual function
can be simplified into something along the lines of

p - list(a=-0.003, b=0.13, c=0.50, E=400)

optim.f - function(p, ST, SM, R, q) {
 res - R 
-(p$a*SM^I(2)+p$b*SM+p$c)*exp(p$E*((1/(283.15-227.13))-(1/(ST+273.15-227.13))
 abserr - abs(res)
 qnum - quantile(abserr, probs=q, na.rm=T)
 res[abserr  qnum] - 0
 res
}

Rmodel - nls.lm(par = p, fn = optim.f, ST = ST, SM = SM, R = R, q = q)

On Thu, 8 May 2008, Fernando Moyano wrote:

 I solved the problem arising from using certain quantile values simply
 by printing the iterations with the argument nprint. This gave me
 correct estimates. I have no idea why.

 - Original Message 
 From: elnano [EMAIL PROTECTED]
 To: r-help@r-project.org
 Sent: Thursday, May 8, 2008 5:43:31 PM
 Subject: Re: [R] function in nls argument


 I've basically solved the problem using the nls.lm function from the
 minpack.lm (thanks Katharine) with some modifications for ignoring residuals
 above a given percentile. This is to avoid the strong influence of points
 which push my modeled vs. measured values away from the 1:1 line.
 I based it on the example given for nls.lm. Here it is:

 R   # soil respiration data
 ST - ST [!is.na(R)] # soil temeprature data. Had to remove na to make
 nls.lm work
 SM - SM [!is.na(R)]# soil moisture data
 R - R [!is.na(R)]
 q - 0.95# quantile

 p - c(a=-0.003, b=0.13, c=0.50, E=400)# model parameters with
 estimated values

 Rf - function(ST, SM, a, b, c, E)
 {
 expr -
 expression((a*SM^I(2)+b*SM+c)*exp(E*((1/(283.15-227.13))-(1/(ST+273.15-227.13)
 eval(expr)
 }

 optim.f - function(p, ST, SM, R, Rfcall)
 {
 res - R - do.call(Rfcall, c(list(ST = ST, SM = SM), as.list(p)))   #
 the nls.lm example divides this by sqrt(R), I don't know why. I removed
 that.
 abserr - abs(res)
 qnum - quantile(abserr, probs=q, na.rm=T)# calculate the q
 quantile of the absolute errors
 res[abserr  qnum]=0  # convert
 residuals above qnum to  0
 res
 }

 Rmodel- nls.lm(par = p, fn = optim.f, Rfcall = Rf, ST = ST, SM = SM, R
 = R)

 summary(Rmodel)

 The only error I still get is when using lower q values. A q of around 0.95
 or less (depending on the dataset) gives me completely wrong parameter
 estimates resulting in negative predicted values. Maybe someone has a
 suggestion here.

 Fernando



 Katharine Mullen wrote:
 
  The error message means that the gradient (first derivative of residual
  vector with respect to the parameter vector) is not possible to work with;
  calling the function qr on the gradient multiplied by the square root of
  the weight vector .swts (in your case all 1's) fails.
 
  If you want concrete advice it would be helpful to provide the commented,
  minimal, self-contained, reproducible code that the posting guide
  requests.  what are the values of ST04, SM08b, ch2no, and tower?
 
  General comments:  If your goal is to minimize sum( (observed -
  predicted)^2), the function you give nls to minimize (optim.fun in your
  case) should return the vector (observed - predicted), not the scalar sum(
  (observed - predicted)^2).  You may want to see the nls.lm function in
  package minpack.lm for another way to deal with the problem.
 
  On Wed, 7 May 2008, Fernando Moyano wrote:
 
  Greetings R users, maybe there is someone who can help
  me with this problem:
 
  I define a function optim.fun and want as output the
  sum of squared errors between predicted and measured
  values, as follows:
 
  optim.fun - function (ST04, SM08b, ch2no, a, b, d, E)
  {
  predR -
  (a*SM08b^I(2)+b*SM08b+d)*exp(E*((1/(283.15-227.13))-(1/(ST04+273.15-227.13
  abserr - abs(ch2no-predR)
  qnum - quantile(abserr, probs=0.95, na.rm=T)
 
 is.na(abserr) - (abserr  qnum)
  errsq - sum(abserr^2, na.rm=T)
  errsq
  }
 
  Then I want to optimize parameters a,b,d and E as to
  minimize the function output with the following:
 
  optim.model-nls(nulldat ~ optim.fun(ST04, SM08b,
  ch2no, a, b, d, E), data=tower,
  start=list(a=-0.003,b=0.13,d=0.50, E=400), na.action =
  na.exclude )
 
  were nulldat=0
  At this point I get the following error message:
 
  Error in qr.default(.swts * attr(rhs, gradient)) :
NA/NaN/Inf in foreign function call (arg 1)
  Warning messages:
  1: In if (na.rm) x - x[!is.na(x)] else if
  (any(is.na(x))) stop(missing values and NaN's not
  allowed if 'na.rm' is FALSE) ... :
the condition has length  1 and 

Re: [R] Intervention/Impact analysis in time series

[Jill E. List]

Hi R users: 
Does any one know about a R library to deal with 
intervention/impact analysis in time series (eg. Box-Tiao et. al. theory?). 


Thank you for your help.

 
Jill List

(408)892-5742

[[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] MANCOVA in R

[Bruno Estigarribia]

Hello,

I have subjects in 4 groups: X1, X2, X3, X4. There are 33 subjects in 
group X1, 35 in X2, 31 in X3, and 46 in group X4. I have 7 continuous 
response variables (actually integers, approximately normal) measured 
for each subject: Y1 to Y7, and two continuous covariates C1, C2 (they 
are both integers).
I want to perform all pairwise comparisons for each response variable 
between groups. I have searched for a way to do a MANCOVA in R to no 
avail. I am familiar with summary.manova, and with Venables  Ripley 
Modern Applied Statistics With S and Everitt's An R and S-Plus 
Companion to Multivariate Analysis. However, I am neither a 
statistician nor a programmer so I am finding it hard to figure this 
out. Can summary.manova be adapted to use covariates? What is the impact 
of the unbalanced design? Can I adjust for multiple comparisons?

Thank you

Bruno Estigarribia
UNC Chapel Hill

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Heat map on large sample size

[Anh Tran]

Hey,
I'm trying to generate a heat map of 30,000 fragments from probably 5-10
samples. Windows complains about memory shortage. Should I resort to Unix
system?

Also, if I only plot 1000 fragments out, they can finish it rather fast.
5000 would take more than 10 minutes. I don't know what to expect for
30,000...

And on the side note, it seems like R only uses up to 50% of CPU while doing
number crunching. Is there anyway to utilize 100% of CPU for R? I'm using R
2.6.2 on Windows XP SP2, average config.

Thanks.

UCLA Neurology Research Lab

-- 
Regards,
Anh Tran

[[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] poisson regression with robust error variance ('eyestudy

[Paul Johnson]

Ted Harding said:
 I can get the estimated RRs from

 RRs - exp(summary(GLM)$coef[,1])

 but do not see how to implement confidence intervals based
 on robust error variances using the output in GLM.


Thanks for the link to the data.  Here's my best guess.   If you use
the following approach, with the HC0 type of robust standard errors in
the sandwich package (thanks to Achim Zeileis), you get almost the
same numbers as that Stata output gives. The estimated b's from the
glm match exactly, but the robust standard errors are a bit off.



### Paul Johnson 2008-05-08
### sandwichGLM.R
system(wget http://www.ats.ucla.edu/stat/stata/faq/eyestudy.dta;)

library(foreign)

dat - read.dta(eyestudy.dta)

### Ach, stata codes factor contrasts backwards

dat$carrot0 - ifelse(dat$carrot==0,1,0)
dat$gender1 - ifelse(dat$gender==1,1,0)

glm1 - glm(lenses~carrot0,  data=dat, family=poisson(link=log))

summary(glm1)

library(sandwich)

vcovHC(glm1)


sqrt(diag(vcovHC(glm1)))


sqrt(diag(vcovHC(glm1, type=HC0)))

### Result:
#  summary(glm1)
# Call:
# glm(formula = lenses ~ carrot0, family = poisson(link = log),
#data = dat)

# Deviance Residuals:
# Min   1Q   Median   3Q  Max
# -1.1429  -0.9075   0.3979   0.3979   0.7734

# Coefficients:
#Estimate Std. Error z value Pr(|z|)
# (Intercept)  -0.8873 0.2182  -4.066 4.78e-05 ***
# carrot0   0.4612 0.2808   1.6420.101
# ---
# Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

# (Dispersion parameter for poisson family taken to be 1)

# Null deviance: 67.297  on 99  degrees of freedom
# Residual deviance: 64.536  on 98  degrees of freedom
# AIC: 174.54

# Number of Fisher Scoring iterations: 5

#  sqrt(diag(vcovHC(glm1, type=HC0)))
# (Intercept) carrot0
#  0.1673655   0.1971117
# 

### Compare against
## http://www.ats.ucla.edu/stat/stata/faq/relative_risk.htm
## robust standard errors are:
## .1682086  .1981048


glm2 - glm(lenses~carrot0 +gender1 +latitude, data=dat,
family=poisson(link=log))


sqrt(diag(vcovHC(glm2, type=HC0)))

### Result:


#  summary(glm2)

#Call:
# glm(formula = lenses ~ carrot0 + gender1 + latitude, family =
poisson(link = log),
# data = dat)

# Deviance Residuals:
# Min   1Q   Median   3Q  Max
# -1.2332  -0.9316   0.2437   0.5470   0.9466

# Coefficients:
#Estimate Std. Error z value Pr(|z|)
# (Intercept) -0.652120.69820  -0.934   0.3503
# carrot0  0.483220.28310   1.707   0.0878 .
# gender1  0.205200.27807   0.738   0.4605
# latitude-0.010010.01898  -0.527   0.5980
# ---
# Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

# (Dispersion parameter for poisson family taken to be 1)

#Null deviance: 67.297  on 99  degrees of freedom
# Residual deviance: 63.762  on 96  degrees of freedom
# AIC: 177.76

# Number of Fisher Scoring iterations: 5

#  sqrt(diag(vcovHC(glm2, type=HC0)))
# (Intercept) carrot0 gender1latitude
# 0.49044963  0.19537743  0.18481067  0.01275001

### UCLA site has:

## .4929193   .1963616  .1857414  .0128142


So, either there is some rounding error or Stata is not using
exactly the same algorithm as vcovHC in sandwich.

-- 
Paul E. Johnson
Professor, Political Science
1541 Lilac Lane, Room 504
University of Kansas

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

Re: [R] Setting matrix dimnames in a list

[statmobile]

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

Re: [R] miniscule font size on R console output

[Yasir Kaheil]

from the upper menu, go to edit/gui preferences 
then you can increase the font size. i used windows so im not sure if this
works on mac.
thanks


juanita choo wrote:
 
 Hi,
 
 I am having a situation where I cannot change the output size of the R
 console. I have played around with the font format menu but  the changes
 are
 only reflected to the script that I type in but not to the output.
 Everytime
 I run a script, I have to go back to font format to increase the output
 script, which is currently showing up as small as the dust on my computer
 screen. I have mac by the way.
 Thanks for your suggestions.
 
 juanita
 
   [[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.
 
 


-
Yasir H. Kaheil, Ph.D.
Catchment Research Facility
The University of Western Ontario 

-- 
View this message in context: 
http://www.nabble.com/miniscule-font-size-on-R-console-output-tp17132821p17136574.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] speeding up a special product of three arrays

[Giuseppe Paleologo]

I am struggling with R code optimization, a recurrent topic on this list.

I have three arrays, say A, B and C, all having the same number of columns.
I need to compute an array D whose generic element is

D[i, j, k] - sum_n A[i, n]*B[j, n]*C[k, n]

Cycling over the three indices and subsetting the columns won't do. Is there
any way to implement this efficiently in R or should I resign to do this in
C?

Thanks,

Giuseppe

[[alternative HTML version deleted]]

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

Re: [R] help with regsubsets

[Yasir Kaheil]

is it possible that regsubsets divides the subset dataset into training and
testing set.. so when it calculates R^2 it's not the same as what you'd get
with lm. You probably need to look into the source code to know if this is
true or not



Edwards, David J wrote:
 
 Hi,
  
 I'm new to R and this mailing list, so I will attempt to state my question
 as appropriately as possible. 
  
 I am running R version 2.7 with Windows XP and have recently been
 exploring the use of the function regsubsets in the leaps package in order
 to perform all-subsets regression. 
  
 So, I'm calling the function as:
  
 sub=regsubsets(X,Y,nbest=5,nvmax=maxnv,method=exhaustive,really.big=F)
  
 I'm dealing with 14 observations and 23 variables and I'm aware of the
 linear dependencies that arise in these situations.  I also know that
 regsubsets performs a branch-and-bound algorithm and that at times it is
 necessary to reorder the variables. 
  
 However, anytime I get the message: Reordering variables and trying
 again the R^2 for a model chosen by regsubsets does not match the R^2 I
 get from fitting the same model using the usual linear model function
 (lm).  If I don't get the Reordering variables and trying again message,
 I don't have any problems with R^2 (i.e. regubsets and lm match). 
  
 Just for example, suppose the best 1-variable model chosen by regsubsets
 gives me an R^2 of 0.665. Fitting that same model using lm gives me an R^2
 of 0.6317.  
  
 Does anyone have any insight on why this may be happening? 
  
 
   [[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.
 
 


-
Yasir H. Kaheil, Ph.D.
Catchment Research Facility
The University of Western Ontario 

-- 
View this message in context: 
http://www.nabble.com/help-with-regsubsets-tp17115128p17136650.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] problem with caretNWS on linux

[Patrick Shields]

Peter,
Are you running the NWS server on the same machine as the R session (ie 
the machine running 'twistd -y /etc/nws.tac')?


Pat

Peter Tait wrote:

Hi,
I am using caretNWS on a RHEL x86_64 system and I am getting an error 
message that is nearly identical to the one occuring in
http://www.r-project.org/nosvn/R.check/r-release-macosx-ix86/caretNWS-00check.txt 



Error in socketConnection(serverHost, port = port, open = a+b, 
blocking = TRUE) :

 unable to open connection
Calls: system.time ... .local - tryCatch - tryCatchList - 
socketConnection

In addition: Warning message:
In socketConnection(serverHost, port = port, open = a+b, blocking = 
TRUE) :

 penguin:8765 cannot be opened

The software versions I am using are the following:
Python 2.3.4
twistd (the Twisted daemon) 1.3.0rc1
nwsserver-1.5.2
R version 2.6.2 (2008-02-08)
nws 1.6.3

Thank you for the help.
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.


--
Pat Shields
REvolution Computing
(203) 777-7442 x250
[EMAIL PROTECTED]

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

Re: [R] help with updating to R2.7

[Vincent Goulet]

Disclaimer: I don't use R on Windows much myself. But: I can't help  
noticing such threads occurring once in a while on r-help.


Why not simply define some folder to host one's personal library  
*outside* of the main R library? As explained in ?Startup, one can put  
something like


R_LIBS=path_to_personal_library:${R_LIBS}

in a .Renviron file located in his/her HOME directory (see http://cran.r-project.org/bin/windows/base/rw-FAQ.html#What-are-HOME-and-working-directories_003f 
 for information on this). Then, install.packages() will, well,  
install packages in this personal folder and R updates do not require  
shuffling  around of the library. At most an update.packages() once in  
a while.


I have been using this scheme for a long time with much success (and  
peace of mind) on Linux and OS X.


Hope this helps.

Le jeu. 8 mai à 15:02, John C Frain a écrit :


When I upgrade in Windows from say 2.6.2 tp 2.7.0 I do the following

1. Install 2.7.0 in a new directory

2 Rename the library subdirectory in the new version from library to  
library2


3 Copy the library subdirectory in 2.6.2 to 2.7.0

4 Copy the contents of library2 to the transferred library.  This
re-installs the updated base packages.

5 Within R run update.packages(checkBuilt = TRUE) and you will be
presented with a list of packages and be given the option of updating
each package.

This procedure is easy to implement.  I have used it for several
updates and it has always worked well

Best Regards

John.

2008/5/8 Duncan Murdoch [EMAIL PROTECTED]:

On 5/8/2008 10:40 AM, ravi wrote:


Hi,
Ouch! That really hurt. But I get the point.
Here's what I did now. I copied all the package folders from R2.6,  
except

for the R.css file, and copied them into the R2.7 folder.


That's your problem.  You've hosed the 2.7 libraries.

You need to reinstall 2.7.  Only copy package folders from 2.6 if  
they don't
exist in 2.7.  If they already exist there, they've already been  
updated.


The older packages may not work in 2.7, but once you run the
update.packages() function, you should get the latest versions.

Duncan Murdoch


In the process, I overwrote the common files that came with the
installation of R2.7.
Here's the output that I obtained in R2.7 :


library()


Error in unique.default(ans) :  2 arguments passed  
to .Internal(unique)

which requires 3


update.packages(checkBuilt=TRUE, ask=FALSE)


Error: could not find function update.packages

My intention with my previous mail was to check and obtain  
confirmation on

details like what I should do with the R.css file.
But I realise that I should have then restricted myself to that  
single
question. That way, I would not be wasting the valuable time of  
the R list

experts, from whom I have received a lot of help previously.
After reading the FAQ once again, I am guessing that the problem  
could be

a case of the internet downloading functions failing.
As suggested in the FAQ, I tried the following commands :


path_to_R\bin\Rgui.exe http_proxy=http://gannet/  
http_proxy_user=ask


Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy


path_to_R\bin\Rgui.exe http_proxy=http://user:[EMAIL PROTECTED]:80/


Error: unexpected symbol in path_to_R\bin\Rgui.exe http_proxy
Should I run the above commands in R, or in the windows command  
window?


Thanking You,
Ravi

- Original Message 
From: Duncan Murdoch [EMAIL PROTECTED]
To: ravi [EMAIL PROTECTED]
Cc: r-help@r-project.org
Sent: Thursday, 8 May, 2008 2:57:07 PM
Subject: Re: [R] help with updating to R2.7

On 5/8/2008 8:30 AM, ravi wrote:


I know that it would be best if I reproduced the exact error  
messages,
but I have tried so many different things now and have lost track  
of the

exact error messages at each stage.


This is not a reasonable request.  Rather than trying those two  
methods
one more time, and recording what goes wrong, you expect someone  
on the list
to do it for you?  If the advice was in the FAQ, presumably it  
works for

others, so the reason it didn't work for you is likely that you
misinterpreted part of it, or your system has a unique problem,  
etc.  The
*only* way someone could diagnose that would be to see the error  
messages.


You need to use some common sense when you're asking for help.

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.





--
John C Frain
Trinity College Dublin
Dublin 2
Ireland
www.tcd.ie/Economics/staff/frainj/home.html
mailto:[EMAIL PROTECTED]
mailto:[EMAIL PROTECTED]

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

Re: [R] poisson regression with robust error variance ('eyestudy

[Achim . Zeileis]

Paul  Ted:

  I can get the estimated RRs from

  RRs - exp(summary(GLM)$coef[,1])

  but do not see how to implement confidence intervals based
  on robust error variances using the output in GLM.


 Thanks for the link to the data.  Here's my best guess.   If you use
 the following approach, with the HC0 type of robust standard errors in
 the sandwich package (thanks to Achim Zeileis), you get almost the
 same numbers as that Stata output gives. The estimated b's from the
 glm match exactly, but the robust standard errors are a bit off.

Thanks for posting the code reproducing this example!

 ### Paul Johnson 2008-05-08
 ### sandwichGLM.R
 system(wget http://www.ats.ucla.edu/stat/stata/faq/eyestudy.dta;)

 library(foreign)

 dat - read.dta(eyestudy.dta)

 ### Ach, stata codes factor contrasts backwards
 dat$carrot0 - ifelse(dat$carrot==0,1,0)
 dat$gender1 - ifelse(dat$gender==1,1,0)

 glm1 - glm(lenses~carrot0,  data=dat, family=poisson(link=log))

 library(sandwich)
 sqrt(diag(vcovHC(glm1, type=HC0)))
 # (Intercept) carrot0
 #  0.1673655   0.1971117

 ### Compare against
 ## http://www.ats.ucla.edu/stat/stata/faq/relative_risk.htm
 ## robust standard errors are:
 ## .1682086  .1981048

I have the solution. Note that the ratio of both standard errors to those
from sandwich is almost constant which suggests a scaling difference. In
sandwich I have implemented two scaling strategies: divide by n
(number of observations) or by n-k (residual degrees of freedom). This
leads to

R sqrt(diag(sandwich(glm1)))
(Intercept) carrot0
  0.1673655   0.1971117
R sqrt(diag(sandwich(glm1, adjust = TRUE)))
(Intercept) carrot0
  0.1690647   0.1991129

(Equivalently, you could youse vcovHC() with types HC0 and HC1,
respectively.)

Their standard errors are between those, so I thought that they used a
different scaling. Here, n = 100 and n-k = 98 which only left

R sqrt(diag(sandwich(glm1) * 100/99))
(Intercept) carrot0
  0.1682086   0.1981048

Bingo! So they scaled with n-1 rather than n or n-k. This is, of course,
also consistent (if the other estimators are), but I wouldn't know of a
good theoretical underpinning for it.

 glm2 - glm(lenses~carrot0 +gender1 +latitude, data=dat,
 family=poisson(link=log))

 ### UCLA site has:
 ## .4929193   .1963616  .1857414  .0128142

R round(sqrt(diag(sandwich(glm2) * 100/99)), digits = 7)
(Intercept) carrot0 gender1latitude
  0.4929204   0.1963617   0.1857417   0.0128142

Only the constant is a bit off, but that is not unusual in replication
exercises.

Some further advertising: If you want to get the full table of
coefficients, you can use coeftest() from lmtest

R library(lmtest)
R coeftest(glm2, sandwich(glm2) * 100/99)

z test of coefficients:

 Estimate Std. Error z value Pr(|z|)
(Intercept) -0.652122   0.492920 -1.3230  0.18584
carrot0  0.483220   0.196362  2.4609  0.01386 *
gender1  0.205201   0.185742  1.1048  0.26926
latitude-0.010009   0.012814 -0.7811  0.43474
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

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

Re: [R] poisson regression with robust error variance ('eyestudy

[Ted Harding]

Once again, Paul, many thanks for your thorough examination
of this question! And for spelling out your approach!!!

It certainly looks as though you're very close to target
(or even spot-on).

I've only one comment -- see at end.

On 08-May-08 20:35:38, Paul Johnson wrote:
 Ted Harding said:
 I can get the estimated RRs from
 
 RRs - exp(summary(GLM)$coef[,1])
 
 but do not see how to implement confidence intervals based
 on robust error variances using the output in GLM.
 
 Thanks for the link to the data. Here's my best guess. If you
 use the following approach, with the HC0 type of robust standard
 errors in the sandwich package (thanks to Achim Zeileis),
 you get almost the same numbers as that Stata output gives.
 The estimated b's from the glm match exactly, but the robust
 standard errors are a bit off.
 
### Paul Johnson 2008-05-08
### sandwichGLM.R
 system(wget http://www.ats.ucla.edu/stat/stata/faq/eyestudy.dta;)
 library(foreign)
 
 dat - read.dta(eyestudy.dta)
### Ach, stata codes factor contrasts backwards
 
 dat$carrot0 - ifelse(dat$carrot==0,1,0)
 dat$gender1 - ifelse(dat$gender==1,1,0)
 
 glm1 - glm(lenses~carrot0,  data=dat, family=poisson(link=log))
 summary(glm1)
 
 library(sandwich)
 vcovHC(glm1)
 sqrt(diag(vcovHC(glm1)))
 sqrt(diag(vcovHC(glm1, type=HC0)))
 
### Result:
#  summary(glm1)
# Call:
# glm(formula = lenses ~ carrot0, family = poisson(link = log),
#data = dat)
 
# Deviance Residuals:
# Min   1Q   Median   3Q  Max
# -1.1429  -0.9075   0.3979   0.3979   0.7734
 
# Coefficients:
#Estimate Std. Error z value Pr(|z|)
# (Intercept)  -0.8873 0.2182  -4.066 4.78e-05 ***
# carrot0   0.4612 0.2808   1.6420.101
# ---
# Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
 
# (Dispersion parameter for poisson family taken to be 1)
 
# Null deviance: 67.297  on 99  degrees of freedom
# Residual deviance: 64.536  on 98  degrees of freedom
# AIC: 174.54
 
# Number of Fisher Scoring iterations: 5
 
#  sqrt(diag(vcovHC(glm1, type=HC0)))
# (Intercept) carrot0
#  0.1673655   0.1971117
 
### Compare against
## http://www.ats.ucla.edu/stat/stata/faq/relative_risk.htm
## robust standard errors are:
## .1682086  .1981048
 
 glm2 - glm(lenses~carrot0 +gender1 +latitude, data=dat,
 family=poisson(link=log))
 
 sqrt(diag(vcovHC(glm2, type=HC0)))
### Result:
#  summary(glm2)
 
#Call:
# glm(formula = lenses ~ carrot0 + gender1 + latitude, family =
 poisson(link = log),
# data = dat)
 
# Deviance Residuals:
# Min   1Q   Median   3Q  Max
# -1.2332  -0.9316   0.2437   0.5470   0.9466
 
# Coefficients:
#Estimate Std. Error z value Pr(|z|)
# (Intercept) -0.652120.69820  -0.934   0.3503
# carrot0  0.483220.28310   1.707   0.0878 .
# gender1  0.205200.27807   0.738   0.4605
# latitude-0.010010.01898  -0.527   0.5980
# ---
# Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
 
# (Dispersion parameter for poisson family taken to be 1)
 
#Null deviance: 67.297  on 99  degrees of freedom
# Residual deviance: 63.762  on 96  degrees of freedom
# AIC: 177.76
 
# Number of Fisher Scoring iterations: 5
 
#  sqrt(diag(vcovHC(glm2, type=HC0)))
# (Intercept) carrot0 gender1latitude
# 0.49044963  0.19537743  0.18481067  0.01275001
 
### UCLA site has:
 
## .4929193   .1963616  .1857414  .0128142
 
 
 So, either there is some rounding error or Stata is not using
 exactly the same algorithm as vcovHC in sandwich.

The above differences look somewhat systematic (though very
small). The percentage differences (vcovHC relative to STATA)
for the two cases you analyse above are

vcovHC HC0:  0.1673655   0.1971117
STATA   :  0.1682086   0.1981048
-
% Difference: -0.5012229  -0.5013003

vcovHC HC0:  0.49044963  0.19537743  0.18481067  0.01275001
STATA: 0.4929193   0.1963616   0.1857414.0128142
-
% Difference: -0.5010293  -0.5012029  -0.5010891  -0.5009287

so, in all cases, very close to -0.501%, despite varying absolute
values. This suggests a very subtle difference in algorithm,
rather than rounding error.
Though there is conceivably a convergence issue in the background.

Does ANYONE here know exactly how STATA does it?

Best wishes,
Ted.


E-Mail: (Ted Harding) [EMAIL PROTECTED]
Fax-to-email: +44 (0)870 094 0861
Date: 08-May-08   Time: 22:28:36
-- XFMail --

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

Re: [R] MLE for noncentral t distribution

[David Scott]

QRMlib has routines for fitting t distributions. Have a look at that 
package. Also sn has routines for skew-t distributions


David Scott


On Thu, 8 May 2008, kate wrote:


I have a data with 236 observations. After plotting the histogram, I found that 
it looks like non-central t distribution. I would like to get MLE for mu and df.

I found an example to find MLE for gamma distribution from fitting distributions 
with R:

library(stats4) ## loading package stats4
ll-function(lambda,alfa) {n-200
x-x.gam
-n*alfa*log(lambda)+n*log(gamma(alfa))-(alfa-
1)*sum(log(x))+lambda*sum(x)} ## -log-likelihood function
est-mle(minuslog=ll, start=list(lambda=2,alfa=1))

Is anyone how how to write down -log-likelihood function for noncentral t 
distribution?

Thanks a lot!!

Kate
[[alternative HTML version deleted]]

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



_
David Scott Department of Statistics, Tamaki Campus
The University of Auckland, PB 92019
Auckland 1142,NEW ZEALAND
Phone: +64 9 373 7599 ext 86830 Fax: +64 9 373 7000
Email:  [EMAIL PROTECTED]

Graduate Officer, Department of Statistics
Director of Consulting, Department of Statistics

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

[Jack Tuomikoski]

Hello,

I have a question regarding the boot function.  I am non-parametrically
bootstrapping a function that I've written to estimate survival for
animals following some closed form estimators (i.e. no optimization
needed).  I'm using the boot function for this and it performs well when
using inputs with say less than 30,000 animals (each record = one animal
with 4 columns of info).  With larger datasets (~50,000) animals, and
1000 bootstraps, the program errors out based on a lack of memory.  The
problem seems to be that the boot function creates the entire 1000
bootstrap databases at once; is the an option for the boot function
where one can go through one iteration at a time (and save the output
from that iteration which contains MUCH less information than the
input)?  In these sorts of studies, samples sizes of 50,000 animals are
often required.

Thanks for any help,

Jack 


[[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] variable of mathematics string not interpreted

[delpacho]

Hi everybody,

my goal is to display symbols on the x-axis of a barplot.
I read some mathematics strings in a file and convert them to an expression
as follows:

tt- scan(file = fstr ,'what' ='character', sep = );

for (iaa in 1:length(tt)) { 
tt[iaa]-do.call(expression, lapply(tt[iaa], as.name));
}

I obtained the following result:

tt =

expression(`alpha[tr]`, beta, widthD, `S[ts]^-5([165,275]Hz)`, 
`S[ts]^2([165,275]Hz)`, `S[ts]^-5([275,460]Hz)`,
`S[r]^2([165,275]Hz)`, 
`S[r]^-5([275,460]Hz)`, `S[r]^3([100,165]Hz)`,
`T[ts]^4([165,275]Hz)`, 
`T[ts]^7([21,30]Hz)`, `T[ts]^8([13,21]Hz)`, `T[r]^7([21,30]Hz)`, 
`T[r]^8([13,21]Hz)`, `T[r]^9([8,13]Hz)`, `E([100,165]Hz)`, 
`E([165,275]Hz)`, `E([60,100]Hz)`, `E([8,13]Hz)`,
`E([30,60]Hz)`, 
`E([21,30]Hz)`)

only the strings without ` are interpreted properly when I put the
property names.arg =tt in barplot (ex:beta displays the greek letter,
whereas `alpha[tr]` displays alpha[tr]), could you explain me why these
symbol ` intervenes when I convert the strings to an expression and/or how 
can I have the symbols associated with the strings displayed properly?

Thanks you very much.
-- 
View this message in context: 
http://www.nabble.com/variable-of-mathematics-string-not-interpreted-tp17138472p17138472.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] variable of mathematics string not interpreted

[Duncan Murdoch]

delpacho wrote:

Hi everybody,

my goal is to display symbols on the x-axis of a barplot.
I read some mathematics strings in a file and convert them to an expression
as follows:

tt- scan(file = fstr ,'what' ='character', sep = );

for (iaa in 1:length(tt)) { 
tt[iaa]-do.call(expression, lapply(tt[iaa], as.name));

}

I obtained the following result:

tt =

expression(`alpha[tr]`, beta, widthD, `S[ts]^-5([165,275]Hz)`, 
`S[ts]^2([165,275]Hz)`, `S[ts]^-5([275,460]Hz)`,
`S[r]^2([165,275]Hz)`, 
`S[r]^-5([275,460]Hz)`, `S[r]^3([100,165]Hz)`,
`T[ts]^4([165,275]Hz)`, 
`T[ts]^7([21,30]Hz)`, `T[ts]^8([13,21]Hz)`, `T[r]^7([21,30]Hz)`, 
`T[r]^8([13,21]Hz)`, `T[r]^9([8,13]Hz)`, `E([100,165]Hz)`, 
`E([165,275]Hz)`, `E([60,100]Hz)`, `E([8,13]Hz)`,
`E([30,60]Hz)`, 
`E([21,30]Hz)`)


only the strings without ` are interpreted properly when I put the
property names.arg =tt in barplot (ex:beta displays the greek letter,
whereas `alpha[tr]` displays alpha[tr]), could you explain me why these
symbol ` intervenes when I convert the strings to an expression and/or how 
can I have the symbols associated with the strings displayed properly?
  


You used as.name() to convert the strings to expressions, so it 
converted them to names.  You really want to parse them, not just 
convert them to names.  For example,


 tt - c(alpha[tr], beta, S[ts]^-5(\[165,275]Hz\))
 barplot(1:3, names.arg=parse(text=tt))

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.

[R] scrime Package simulatedSNP function

[Claire_6700]

Hello,

I need some help with the simulatedSNPs function from scrime package.

I am trying to simulate some genotype of a case/control disease locus. The
allele frequence are cases/controls

Sample cases  controls
2000 .5.10
1500 .6.40

In each of the row, i need to simulate 100 snp and calculate the pvalue

##Download Scrime Package###
library(scrime)
n.obs-1000
n.snp-100
vec.ia-1
simulateSNPs(n.obs, n.snp, vec.ia, prop.explain = 1,
list.ia.val = NULL, vec.ia.num = NULL, maf = c(0.1, 0.12),
prob.val = rep(1/3, 3), list.equal = NULL, prob.equal = 0.8,
rm.redundancy = TRUE, shuffle = FALSE, shuffle.obs = FALSE, rand = NA)

What is the right parameter. I am pretty new with R.

wrong result.
  Interaction Cases Controls
1   SNP1 == 2   5000

any help will be appreciated.

thanks
Claire
-- 
View this message in context: 
http://www.nabble.com/scrime-Package-simulatedSNP-function-tp17138820p17138820.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] poisson regression with robust error variance ('eyestudy

[Ted Harding]

I'd like to thank Paul Johnson and Achim Zeileis heartily
for their thorough and accurate responses to my query.

I think that the details og how to use the procedure, and
of its variants, which they have sent to the list should
be definitive -- and very helpfully usable -- for folks
like myself who may in future grope in the archives
concerning this question.

And, just to confirm, it all worked perfectly for me
in the end.

Best wishes,
Ted.


E-Mail: (Ted Harding) [EMAIL PROTECTED]
Fax-to-email: +44 (0)870 094 0861
Date: 09-May-08   Time: 00:56:14
-- XFMail --

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

[R] histogram

[Roslina Zakaria]

Dear R-expert,
For histogram function, can we get the table of bin and frequency like in 
excel, together with the histogram?
Therefore, we can check the number of data included.
Thank you so much for your attention and help.


  


[[elided Yahoo spam]]

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

Re: [R] histogram

[Rolf Turner]

On 9/05/2008, at 12:04 PM, Roslina Zakaria wrote:


Dear R-expert,
For histogram function, can we get the table of bin and frequency  
like in excel, together with the histogram?

Therefore, we can check the number of data included.
Thank you so much for your attention and help.


?hist

Look at the ``Value'' material.

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] histogram

[David Scott]

On Thu, 8 May 2008, Roslina Zakaria wrote:


Dear R-expert,
For histogram function, can we get the table of bin and frequency like in 
excel, together with the histogram?
Therefore, we can check the number of data included.
Thank you so much for your attention and help.



Easy one: just use the argument plot=FALSE and hist will produce the 
counts etc for you


David Scott

_
David Scott Department of Statistics, Tamaki Campus
The University of Auckland, PB 92019
Auckland 1142,NEW ZEALAND
Phone: +64 9 373 7599 ext 86830 Fax: +64 9 373 7000
Email:  [EMAIL PROTECTED]

Graduate Officer, Department of Statistics
Director of Consulting, Department of Statistics

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

Re: [R] histogram

[Ted Harding]

On 09-May-08 00:12:46, David Scott wrote:
 On Thu, 8 May 2008, Roslina Zakaria wrote:
 
 Dear R-expert,
 For histogram function, can we get the table of bin and
 frequency like in excel, together with the histogram?
 Therefore, we can check the number of data included.
 Thank you so much for your attention and help.
 
 Easy one: just use the argument plot=FALSE and hist will produce the 
 counts etc for you
 
 David Scott

Well, there's a slightly nicer way than that! (Don't forget that
we're competing with Excel here ... ).

Since the value of hist() includes breakpoints and counts,
you can get the plot and the table as follows:

  set.seed(54321)
  X - 5*rnorm(500)
  H-hist(X)### Here you get the histogram drawn
  B-H$breaks ; C-H$counts ; k-length(B)
  cbind(From=B[1:(k-1)],To=B[2:k],Count=C) ### Here bins and frequencies
# From  To Count
#[1,]  -15 -1016
#[2,]  -10  -569
#[3,]   -5   0   181
#[4,]0   5   149
#[5,]5  1071
#[6,]   10  1514

So we can indeed get the table of bin and frequency like
in excel, together with the histogram.

Best wishes to all,
Ted.


E-Mail: (Ted Harding) [EMAIL PROTECTED]
Fax-to-email: +44 (0)870 094 0861
Date: 09-May-08   Time: 01:26:35
-- XFMail --

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

Re: [R] histogram

[Jorge Ivan Velez]

Hi all,

Why not only?

set.seed(54321)
 X - 5*rnorm(500)
hist(X,label=TRUE,ylim=c(0,200))

Thanks,

Jorge



On Thu, May 8, 2008 at 8:26 PM, Ted Harding [EMAIL PROTECTED]
wrote:

 On 09-May-08 00:12:46, David Scott wrote:
  On Thu, 8 May 2008, Roslina Zakaria wrote:
 
  Dear R-expert,
  For histogram function, can we get the table of bin and
  frequency like in excel, together with the histogram?
  Therefore, we can check the number of data included.
  Thank you so much for your attention and help.
 
  Easy one: just use the argument plot=FALSE and hist will produce the
  counts etc for you
 
  David Scott

 Well, there's a slightly nicer way than that! (Don't forget that
 we're competing with Excel here ... ).

 Since the value of hist() includes breakpoints and counts,
 you can get the plot and the table as follows:

  set.seed(54321)
  X - 5*rnorm(500)
  H-hist(X)### Here you get the histogram drawn
  B-H$breaks ; C-H$counts ; k-length(B)
  cbind(From=B[1:(k-1)],To=B[2:k],Count=C) ### Here bins and frequencies
 # From  To Count
 #[1,]  -15 -1016
 #[2,]  -10  -569
 #[3,]   -5   0   181
 #[4,]0   5   149
 #[5,]5  1071
 #[6,]   10  1514

 So we can indeed get the table of bin and frequency like
 in excel, together with the histogram.

 Best wishes to all,
 Ted.

 
 E-Mail: (Ted Harding) [EMAIL PROTECTED]
 Fax-to-email: +44 (0)870 094 0861
 Date: 09-May-08   Time: 01:26:35
 -- XFMail --

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


[[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] Is it possible to do fancy area plots?

[Sean Carmody]

Does anyone have any ideas about how you could use R to produce a fancy area
plot like this one in the NY Times? http://tinyurl.com/6rr22g
Regards,
Sean,

[[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] Is it possible to do fancy area plots?

[Charilaos Skiadas]

On May 8, 2008, at 9:11 PM, Sean Carmody wrote:

Does anyone have any ideas about how you could use R to produce a  
fancy area

plot like this one in the NY Times? http://tinyurl.com/6rr22g


I certainly hope not, I wouldn't want my favorite statistics program  
to produce an area graph where the sizes of the areas are not  
proportional to the percentages represented, and where the shapes of  
the areas are so irregular as to make effective area comparisons near  
impossible...


Unless my eyes are seriously deceiving me, which is quite possible.

It would be interesting however, to find out how such a graph is  
constructed. I wonder if hyperbolic geometry enters the picture.



Regards,
Sean,


Haris Skiadas
Department of Mathematics and Computer Science
Hanover College

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

Re: [R] Is it possible to do fancy area plots?

[Sean Carmody]

After a bit more searching, I've discovered that this chart is a variant of
the treemap or map of the market. I'll play around with the sample code
posted here https://stat.ethz.ch/pipermail/r-sig-finance/2006q2/000880.html,
but if anyone's taken that further, I'd be keen to know. I'm happy to use
the more conventional box approach than circles and irregular areas.
Sean.

On Fri, May 9, 2008 at 11:35 AM, Charilaos Skiadas [EMAIL PROTECTED]
wrote:

 On May 8, 2008, at 9:11 PM, Sean Carmody wrote:

  Does anyone have any ideas about how you could use R to produce a fancy
 area
 plot like this one in the NY Times? http://tinyurl.com/6rr22g


 I certainly hope not, I wouldn't want my favorite statistics program to
 produce an area graph where the sizes of the areas are not proportional to
 the percentages represented, and where the shapes of the areas are so
 irregular as to make effective area comparisons near impossible...

 Unless my eyes are seriously deceiving me, which is quite possible.

 It would be interesting however, to find out how such a graph is
 constructed. I wonder if hyperbolic geometry enters the picture.

  Regards,
 Sean,


 Haris Skiadas
 Department of Mathematics and Computer Science
 Hanover College






[[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] Defining factors for graphical presentations and incorporation of a smoother overlay line

[Sparkes, Earl]

Hello,
I have some photosynthesis data that is grouped by SPECIES and REPLICATION
I have tried to develop a grouped data set with the following:
photo.gd - groupedData(Photo~PARi|Species, 
 data = photo.frm,  
 labels = list(x = PPFD (expression(paste)µ mol photons m^-2 s^-1), 
   y = Net photosynthetic rate (expression(paste)µ mol CO[2]m^-2 
s^-1)) ,
 units = list(y = (std temp))
 )
plot(photo.gd,aspect=1)

How do I add REPLICATION as an ordered factor to the grouped data setup?
How can I incorporate a mean (distance weight least squares) line as part of a 
graph for each replicate set.


Thank you, 
Earl Sparkes 

DISCLAIMER**...{{dropped:15}}

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

Re: [R] significance threshold in CCF

[E C]

Ahh, thanks, that helped! Is the standard error being calculated by 1/sqrt(N-3) 
though? I ask because visually inspecting the plots of the 6 CCFs I've done, 
only 4 have a C.I. line that look about right according to my own 
calculations using this formula. The other 2 are a little below my calculations 
(although probably close enough). However, maybe I caught these two only 
because in these cases, my calculation happened to be a value above a tick mark 
while the line shown was below that tick mark. Maybe I'm not using the exactly 
right formula? Is there a way of extracting the value used? Thanks again!



 Date: Thu, 8 May 2008 19:54:25 +0100
 From: [EMAIL PROTECTED]
 To: [EMAIL PROTECTED]
 Subject: RE: [R] significance threshold in CCF
 CC: r-help@r-project.org
 
 On 08-May-08 18:23:27, E C wrote:
  Hi everyone,
  
  When the CCF between two series of observations is plotted in R, a line
  indicating (presumably) the significance threshold appears across the
  plot. Does anyone know how this threshold is determined (it is
  different for each set of series) and how its value can be extracted
  from R? I've tried saving the CCF into an object and unclassing the
  object, but there's nothing there to indicate this.
  Some sample code to show what I mean:
  x - c(1,2,3,4,5,6,7,8,9)
  y - c(3,4,5,6,7,8,9,10,11)
  ccf(x, y, plot=T)
  Thanks in advance!
 
 Have a look at ?plot.acf (also used for plotting ccf objects)
 Ted.
 
 
 E-Mail: (Ted Harding) [EMAIL PROTECTED]
 Fax-to-email: +44 (0)870 094 0861
 Date: 08-May-08   Time: 19:54:20
 -- XFMail --

_
If you like crossword puzzles, then you'll love Flexicon, a game which 
comb[[elided Hotmail spam]]

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