Re: [R] How can I fitting this function from values.

[Frede Aakmann Tøgersen]

R is case sensitive ;-) so in your list of start values you must have B=3 and 
not b=3. Perhaps you have an object named b of a different length than your 
data somewhere giving rise to your error message.

Also brute force wont work as Joerg point out. We have 

A*1/(1+t/B)*1/sqrt(1+t/C) = D/(B+t)/sqrt(C+t),

where D = A*B*sqrt(C).

We notice that we need B > -t and C > -t.


Now one can use the "port" algorithme of nls() with constraints on parameters:


foo <- function(t,ccc){
  D <- ccc[1]
  B <- ccc[2]
  C <- ccc[3]
  D/(B+t)/sqrt(C+t)
}

ttt <- seq(1,10,length=100)
yyy <- foo(ttt,c(5*3*sqrt(2),3,2))
set.seed(42)

dat.test <- data.frame(t=ttt,y=yyy+rnorm(100,0,.1))


fit.test <- nls(y~D/(B+t)/sqrt(C+t),data=dat.test,
start=c(D=5*3*sqrt(2),B=3,C=2),
lower=c(-Inf,-1,-1),alg="port")
# -Inf could be 0, -1 could be -min(dat.test$t)

with(dat.test,plot(t,y))
with(dat.test,lines(t,foo(t,c(21.21,3,2)),col="green")) # Looks cool.
lines(ttt,fitted(fit.test),col="red")


Best regards

Frede Aakmann Tøgersen
Scientist


UNIVERSITY OF AARHUS
Faculty of Agricultural Sciences
Dept. of Genetics and Biotechnology
Blichers Allé 20, P.O. BOX 50
DK-8830 Tjele

Phone:   +45 8999 1900
Direct:  +45 8999 1878

E-mail:  [EMAIL PROTECTED]
Web:   http://www.agrsci.org

This email may contain information that is confidential.
Any use or publication of this email without written permission from Faculty of 
Agricultural Sciences is not allowed.
If you are not the intended recipient, please notify Faculty of Agricultural 
Sciences immediately and delete this email.


 

> -Oprindelig meddelelse-
> Fra: [EMAIL PROTECTED] 
> [mailto:[EMAIL PROTECTED] På vegne af Rolf Turner
> Sendt: 12. november 2007 21:16
> Til: Horacio Castellini
> Cc: r-help@r-project.org
> Emne: Re: [R] How can I fitting this function from values.
> 
> 
> On 13/11/2007, at 8:27 AM, Horacio Castellini wrote:
> 
> > Hi all, I'd like fit this function:
> >
> > G(t)=A*1/(1+t/B)*1/sqrt(1+t/C)
> >
> > where A, B and C are fitting constants that one like search. I have 
> > got a fcs<-(t,G(t)) value table in a data frame.
> >
> > Any one can help me? Tahnks Horacio.
> 
> 
> I ***thought*** that nls would solve your problem in a trice.
> But I tried a toy example before replying to you, and it messed up
> mightily:
> 
> foo <- function(t,ccc){
> ccc[1]*1/(1+t/ccc[2])*1/sqrt(1+1/ccc[3])
> }
> ttt <- seq(1,10,length=100)
> yyy <- foo(ttt,c(5,3,2))
> set.seed(42)
> dat.test <- data.frame(t=ttt,y=yyy+rnorm(100,0,0.1))
> with(dat.test,plot(t,y))
> with(dat.test,lines(t,foo(t,c(5,3,2 # Looks cool.
> fit.test <- nls(y~A*1/(1+t/B)*1/sqrt(1+t/C),start=list
> (A=5,b=3,C=2),data=dat.test)
> # Gives error; no idea what the implications are.
> 
> The error message:
> 
> Error in `[[<-.data.frame`(`*tmp*`, var, value = 
> c(0.70934153524875, 0.761288405463172,  :
>replacement has 65 rows, data has 100 In addition: Warning message:
> In t/B : longer object length is not a multiple of shorter 
> object length
> 
> I don't think I have *ever* had a good experience with nls(); 
> it is always pretty flaky.
> (Well, maybe it worked for me once ..., can't remember for sure. :-) )
> 
>   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.
> 

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

[Dörte Salecker]

dear r-helpers

i've got a table that in extracts looks like this:

   V1   V2  V3V4  V5
1  01/01/1975 00:00:00 125.837 3.691 296.618
2  01/01/1975 01:00:00 124.799 3.679 281.307
3  01/01/1975 02:00:00 111.607 3.536 281.307
4  02/24/1976 11:00:00  21.602 2.555  93.893
5  02/24/1976 12:00:00  27.804 2.623  93.893
6  02/24/1976 13:00:00  26.105 2.604 114.716
7  10/18/1977 00:00:00  18.850 2.525 330.203
8  10/18/1977 01:00:00  21.561 2.555 320.465
9  10/18/1977 02:00:00  21.941 2.559 314.505
10 11/05/1978 01:00:00  24.981 2.592 310.464
11 11/05/1978 02:00:00  19.066 2.528 310.464
12 11/05/1978 03:00:00  26.198 2.605 305.749
13 06/29/1979 19:00:00  78.829 3.179 306.199
14 06/29/1979 20:00:00  92.221 3.325 306.199
15 06/29/1979 21:00:00  72.364 3.108 303.049

i need to find the annual maximun for each year in V3, V4 and V5 So that the 
result would be a table like this:

   V1  V4
1  01/01/1975   3.691 
2  02/24/1976   2.623
3  10/18/1977   2.559
.
.
.and so on.

i hope you can help me with this,

thank you in advance,
doerte

--

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

Re: [R] Query an Access database based on a date attribute

[Prof Brian Ripley]

This is a question about SQL, or more precisely, Microsoft's peculiar 
dialect of SQL.  You haven't even mentioned (let alone credited) package 
RODBC which you appear to be using.

In SQL queries you need to quote numeric values if you want them to be 
treated as character.  Why did you quote 'alt-bos' and not '2007-11-20' ? 
As I recall, some DBMSs need date values quoted and some do not, so try 
quoting.

On Mon, 12 Nov 2007, Tudor Bodea wrote:

> Dear useRs,
>
> I would like to query an Access database through R based on a date 
> attribute but, unfortunately, I fail to do so. For example, the table 
> test_table of the test.mdb looks like:
>
>   ID cd  market competitor dd price
> 1   1 2007-11-20 atl-bos  delta 2007-11-20   210
> 2   2 2007-11-20 atl-bos  delta 2007-11-21   190
> 3   3 2007-11-20 atl-bos  delta 2007-11-22   180
> 4   4 2007-11-20 atl-bos united 2007-11-20   205
> 5   5 2007-11-20 atl-bos united 2007-11-21   195
> 6   6 2007-11-20 atl-bos united 2007-11-22   175
> 7   7 2007-11-20 atl-sfa  delta 2007-11-20   350
> 8   8 2007-11-20 atl-sfa  delta 2007-11-21   320
> 9   9 2007-11-20 atl-sfa  delta 2007-11-22   300
> 10 10 2007-11-20 atl-sfa aa 2007-11-20   340
> 11 11 2007-11-20 atl-sfa aa 2007-11-21   320
> 12 12 2007-11-20 atl-sfa aa 2007-11-22   300
>
> In this context, I try to get all the records for which market is 
> atl-bos, competitor is delta and dd is 2007-11-20 (first record above). 
> To do this I used
>
>> # channel <- odbcConnectAccess("test.mdb")
>> res <- sqlQuery(channel, "select * from test_table where market = 'atl-bos'
> and competitor = 'delta' and dd = 2007-11-20")
>
> but the result seems to be incorrect.
>
>> res
> [1] ID cd market competitor dd price
> <0 rows> (or 0-length row.names)
>
> Does any of you know what it is going on? I suspect that the last 
> expression in the select statement is not executed/interpreted correctly 
> and, as a consequence, the query leads to the above empty set.
>
> Thank you.
>
> Tudor
>
> --
> Tudor Dan Bodea
> Georgia Institute of Technology
> School of Civil and Environmental Engineering
> Web: http://www.prism.gatech.edu/~gtg757i

-- 
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] Using lme (nlme) to find the conditional variance of therandom effects

[Rick Bilonick]

On Tue, 2007-11-13 at 01:03 -0500, Rick Bilonick wrote:

> 
> Is there some way to get ranef with postVar=TRUE to show what the
> variances are, or what the lower and upper bounds are? qqmath makes nice
> plots but I need to obtain the numerical values.
> 
> Rick B.
> 
I found a way:

attr(ranef(lmer.13,postVar=TRUE)[[2]],"postVar")

But I still don't understand why it's not OK to access the bVar slot
directly. 

Also, the code I originally showed and the results from ranef are very
similar with a correlation of 0.9983 (it varies very slightly from
subject to subject):

>
round(data.frame(a=as.numeric([EMAIL 
PROTECTED](attr(VarCorr(lmer.13),"sc")^2)[1]),
  b=as.numeric(attr(ranef(lmer.13,postVar=TRUE)[[2]],"postVar"))),10)
  ab
1  5.41e-08 5.44e-08
2  4.77e-08 4.83e-08
3  6.24e-08 6.25e-08
4  4.44e-08 4.52e-08
5  6.50e-08 6.50e-08
6  2.67e-08 2.92e-08
7  5.07e-08 5.12e-08
8  6.43e-08 6.43e-08
9  3.64e-08 3.79e-08
10 4.86e-08 4.92e-08
11 6.33e-08 6.33e-08
12 3.44e-08 3.60e-08
13 4.16e-08 4.26e-08
14 3.69e-08 3.83e-08
15 5.96e-08 5.97e-08
16 6.46e-08 6.46e-08
17 3.28e-08 3.46e-08
18 4.71e-08 4.77e-08
19 5.18e-08 5.22e-08
20 2.81e-08 3.04e-08
21 3.97e-08 4.09e-08
22 5.70e-08 5.72e-08
23 6.06e-08 6.07e-08
24 3.23e-08 3.42e-08
25 4.94e-08 4.99e-08
26 5.35e-08 5.38e-08
27 3.86e-08 3.98e-08
28 6.73e-08 6.73e-08
29 4.68e-08 4.74e-08
30 6.15e-08 6.16e-08
31 4.67e-08 4.74e-08
32 2.04e-08 2.37e-08
33 3.45e-08 3.61e-08
34 6.28e-08 6.29e-08
35 5.53e-08 5.55e-08

Not sure why they are not exactly the same.

Rick B.

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

Re: [R] a repetition of simulation

[sigalit mangut-leiba]

Thank you,
Sigalit.


On 11/13/07, Moshe Olshansky <[EMAIL PROTECTED]> wrote:
>
> Hi,
>
> Look at names(log_v) in your notation to see what you
> really need.
> Then you could do something like:
>
> results <- list(600)
> for (ij in 1:600) {
> do what you did
> results[[ij]] <- log_v
> }
>
> So now you have a list of 600 results and can do
> anything you need (look at lapply/sapply).
>
> Regards,
>
> Moshe.
>
> --- sigalit mangut-leiba <[EMAIL PROTECTED]> wrote:
>
> > Hello,
> > I have a simple (?) simulation problem.
> > I'm doing a simulation with logistic model and I
> > want to reapet it 600
> > times.
> > The simulation looks like this:
> >
> > z <- 0
> > x <- 0
> > y <- 0
> > aps <- 0
> > tiss <- 0
> > for (i in 1:500){
> > z[i] <- rbinom(1, 1, .6)
> > x[i] <- rbinom(1, 1, .95)
> > y[i] <- z[i]*x[i]
> > if (y[i]==1) aps[i] <- rnorm(1,mean=13.4, sd=7.09)
> > else aps[i] <-
> > rnorm(1,mean=12.67, sd=6.82)
> > if (y[i]==1) tiss[i] <-
> > rnorm(1,mean=20.731,sd=9.751) else  tiss[i] <-
> > rnorm(1,mean=18.531,sd=9.499)
> > }
> > v <- data.frame(y, aps, tiss)
> > log_v <- glm(y~., family=binomial, data=v)
> > summary(log_v)
> >
> > I want to do a repetition of this 600 times (I want
> > to have 600 logistic
> > models), and see all the coefficients of the
> > covariates aps & tiss.
> > Thanks in advance,
> > Sigalit.
> >
> >   [[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] cor.test() and binary sequences

[andrew collier]

peter, thanks for your help with my questions regarding cor.test(). i have 
another question though: does this function make any assumptions about the 
underlying distribution of the two sequences? does it assume that they have a 
gaussian distribution?

i ask because the data that i am working with is two binary sequences. just 
series of 0 and 1. will the confidence intervals and p-values generated by 
cor.test() still apply?

cheers,
andrew.

--
Get a free email account with anti spam protection.
http://www.bluebottle.com/tag/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] percent inclusion in cph

[sushi4u]

Dear R-help,

I have to calculate the percent inclusion of each variable in a bootstrap 
validation of a cox proportional hazards model(described in Sauerbrei and 
Schumacher, Stat Med 11:1093, 1992).

First I need to get a bootstrap sample from my dataset, which I did with the 
sample function. Then I tried to run a cph model and looked which covariates 
are significant. This I would repeat 200 times and at the end calculate how 
many percent the covariates were included. 

This is what I entered:
boot1 <- sample(Dataset, 300, replace=T)
cph1 <- cph(Surv(months,status) ~ cov1 + cov2 + cov3 + cov4, data=boot1)

Unfortunately, I get exactly the same results (coefficient, SE, p-value) then 
when I would fit a Cox model without drawing a bootstrap sample before.

How do I do it right? Or is there another way to calculate the percentage?

Sorry for my bad English. 

Thanks.

Dott. Mario Rossi
University of Foggia, Italy

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

[sushi4u]

Dear R-help,

I have to calculate the percent inclusion of each variable in a bootstrap 
validation of a cox proportional hazards model(described in Sauerbrei and 
Schumacher, Stat Med 11:1093, 1992).

First I need to get a bootstrap sample from my dataset, which I did with the 
sample function. Then I tried to run a cph model and looked which covariates 
are significant. This I would repeat 200 times and at the end calculate how 
many percent the covariates were included. 

This is what I entered:
boot1 <- sample(Dataset, 300, replace=T)
cph1 <- cph(Surv(months,status) ~ cov1 + cov2 + cov3 + cov4, data=boot1)

Unfortunately, I get exactly the same results (coefficient, SE, p-value) then 
when I would fit a Cox model without drawing a bootstrap sample before.

How do I do it right? Or is there another way to calculate the percentage?

Sorry for my bad English. 

Thanks.

Dott. Mario Rossi
University of Foggia, Italy

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

Re: [R] Using lme (nlme) to find the conditional variance of therandom effects

[Rick Bilonick]

On Mon, 2007-11-12 at 16:45 -0500, Doran, Harold wrote:
> No, don't reach into the bVar slot. Use the proper extractor function
> ranef() with postVar=T. There is no similar function for lme() 
> 
> > -Original Message-
> > From: [EMAIL PROTECTED] 
> > [mailto:[EMAIL PROTECTED] On Behalf Of Rick Bilonick
> > Sent: Monday, November 12, 2007 4:40 PM
> > To: R Help
> > Subject: [R] Using lme (nlme) to find the conditional 
> > variance of therandom effects
> > 
> > Using lmer in the lme4 package, you can compute the 
> > conditional variance-covariance matrix of the random effects 
> > using the bVar slot:
> > 
> > bVar: A list of the diagonal inner blocks (upper triangles 
> > only) of the positive-definite matrices on the diagonal of 
> > the inverse of ZtZ+Omega.
> > With the appropriate scale factor (and conversion to a 
> > symmetric matrix) these are the conditional 
> > variance-covariance matrices of the random effects.
> > 
> > Is there anything similar in the nlme package using the lme function?
> > 
> > Rick B.
> > 
> > __


Is there some way to get ranef with postVar=TRUE to show what the
variances are, or what the lower and upper bounds are? qqmath makes nice
plots but I need to obtain the numerical values.

Rick B.

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

Re: [R] a repetition of simulation

[Stephen Weigand]

Dear Sigalit,

On Nov 12, 2007 2:18 PM, sigalit mangut-leiba <[EMAIL PROTECTED]> wrote:
> Hello,
> I have a simple (?) simulation problem.
> I'm doing a simulation with logistic model and I want to reapet it 600
> times.
> The simulation looks like this:
>
> z <- 0
> x <- 0
> y <- 0
> aps <- 0
> tiss <- 0
> for (i in 1:500){
> z[i] <- rbinom(1, 1, .6)
> x[i] <- rbinom(1, 1, .95)
> y[i] <- z[i]*x[i]
> if (y[i]==1) aps[i] <- rnorm(1,mean=13.4, sd=7.09) else aps[i] <-
> rnorm(1,mean=12.67, sd=6.82)
> if (y[i]==1) tiss[i] <- rnorm(1,mean=20.731,sd=9.751) else  tiss[i] <-
> rnorm(1,mean=18.531,sd=9.499)
> }
> v <- data.frame(y, aps, tiss)
> log_v <- glm(y~., family=binomial, data=v)
> summary(log_v)
>
> I want to do a repetition of this 600 times (I want to have 600 logistic
> models), and see all the coefficients of the covariates aps & tiss.
> Thanks in advance,
> Sigalit.
>

I won't comment on your general approach to simulation but does this
answer your question?

### Create an empty matrix to hold the results as they
### are created
Results <- matrix(NA, ncol  = 3, nrow = 600)

for(i in 1:600) {
  generate your data
  Results[i, ] <- coef(glm(y~., family=binomial, data=v))
  }

summary(Results)

Stephen

-- 
Rochester, Minn. USA

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

[Rolf Turner]

I am bewildered by the behaviour of tickmarks as demonstrated by the
following code.  (What I'm trying to do is draw a single tick mark  
extending
from the axis all the way down to the tick label, which is two lines  
from
the axis to make sure it doesn't overlap with the ``ordinary'' tick  
labels.)

# Try 1:
# Gap between tickmark and label.
plot(1:10)
axis(side = 1, at = 2.35, labels = "2.35", tick = TRUE,
  tcl=1,line=1,col="red", col.axis = "red")
readline('Go? ')

# Try 2:
# Gap between tickmark and axis.
plot(1:10)
axis(side = 1, at = 2.35, labels = "2.35", tick = TRUE,
  tcl=-1,line=1,col="red", col.axis = "red")
readline('Go? ')

# Try 3:
# Gap between tickmark and label not filled.
plot(1:10)
axis(side = 1, at = 2.35, labels = "2.35", tick = TRUE,
tcl=1,line=1,col="red", col.axis = "red")
axis(side = 1, at = 2.35, labels ="", tick = TRUE, tcl=-1,col="red")
readline('Go? ')

# Try 4:
# Gap between tickmark and label ***is*** filled.
# Gives what I want, but ***why***?
plot(1:10)
axis(side = 1, at = 2.35, labels = "2.35", tick = TRUE,
  tcl=-1,line=1,col="red", col.axis = "red")
axis(side = 1, at = 2.35, labels ="", tick = TRUE, tcl=-1,col="red")

What I can't understand is why the tcl value in the second call to  
axis() seems
to have the opposite effect of its value in the first call to axis().

Can anyone enlighten me?

BTW setting tcl=2 and tcl=-2 have impacts which also bewilder me and  
do not
have the desired effect.

cheers,

Rolf Turner

P. S. Version info:

 > version
_
platform   i386-apple-darwin8.10.1
arch   i386
os darwin8.10.1
system i386, darwin8.10.1
status
major  2
minor  6.0
year   2007
month  10
day03
svn rev43063
language   R
version.string R version 2.6.0 (2007-10-03)


##
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] help on drawing a tree with "ape"?

[Hua Li]

Thank you very much for the help, Ben!

As a follow up, is there a way to specify the labels,
through the way the text is written, rather than
reading the edge positions from the graph? For
example, 

mytree =
"((A:51.78,(C:24.6,D:24.6):27.18):40.06,B:91.84):0.0;"
plot(read.tree(text = mytree))
edgelabels(c("a-b","b","b-d","d","d","a"))

would give the graph I want. You may notice that if
set a = 91.84; b= 51.78; d= 24.6, the labels I put on
each edge is the corresponding lengths in variable
form.  The way I do it right now need me to first draw
the tree and then manually set the labels, and I
suppose there should be an easier way to get them by
reading the way the tree is written
("((A:51.78,(C:24.6,D:24.6):27.18):40.06,B:91.84):0.0;")?


Sorry for my low level questions. Just started using
"ape" and research on phylogenetic trees.

Best wishes,

Hua





--- Ben Bolker <[EMAIL PROTECTED]> wrote:

> 
> 
> 
> Hua Li wrote:
> > 
> > Dear Simon and everyone,
> > 
> > Thanks for response. 
> > 
> > Specifically, I have a tree with species "A, B, C"
> and
> > can be written as
> > "((A:45.15,C:45.15):46.19,B:91.34):0.0;" . If I
> use
> > the R command
> >> mytreeABC
> > 
> > [1] "((A:45.15,C:45.15):46.19,B:91.34):0.0;"
> > 
> >> plot(read.tree(text = mytreeABC))
> > 
> > --- Simon Blomberg <[EMAIL PROTECTED]> wrote:
> > 
> >> ?edgelabels
> >> 
> > 
> 
>Did you look at the examples in the file
> suggested by Simon?
> Unfortunately, the PS file but not the JPG came
> through in the
> version I read.
> 
>   It seems to me you may want _node_ labels rather
> than _edge_ labels?See below ...
> 
> mytreeABC = "((A:45.15,C:45.15):46.19,B:91.34):0.0;"
> t1 = read.tree(text=mytreeABC)
> library(ape)
> plot(t1)
> edgelabels(1:4)
> nodelabels(c("a","b"))
> 
> -- 
> View this message in context:
>
http://www.nabble.com/help-on-drawing-a-tree-with-%22ape%22--tf4789546.html#a13711995
> 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.
> 



  

Never miss a thing.  Make Yahoo your home page.

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

Re: [R] Multivariate time series

[Paul Gilbert]

Package dse does do linear, multivariate time series with multivariate 
exogenous inputs, using VAR, multivariate ARMA, and state-space models, 
just like you are describing. You can specify the model or have it 
automatically determined. Perhaps you could give a few more details 
about what you are looking for. (I assume you have looked at the dse 
Users Guide distributed with the bundle, or some of the help, for 
example ?estVARXls, ?bft, and ?forecast. Sections 6 and 7 of the guide 
are about the prediction problem, very much like you describe.)

Paul

Giusy wrote:
> Hello to everyone!
> I have a question for you..I need to predict multivariate time series, for
> example sales of 2 products related one to the other, having the 2 prices
> like inputs..
> Is there in R a function to do it? I saw dse package but I didn't find what
> a I'm looking for..
> Could anyone help me?
> Thank you very much
> Giusy


La version française suit le texte anglais.



This email may contain privileged and/or confidential in...{{dropped:26}}

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

Re: [R] connection diagram

[elw]

hi,

You should probably be looking at the functions in the following packages:

sna
network(s)
graph
dynamicgraph
mathgraph
igraph
Matrix

and a few others ;)

what you're describing sounds like, to my ear, a restricted social network 
diagram;  many of the problems you describe are typical of such problems, 
and are accounted for in the packages described above.

The most difficult part is likely to be the plots;  handling an annotated, 
weighted, complex network is fairly straightforward in terms of data 
handling and analytic tools (e.g. regressions on networks are common...).

--elijah




On Tue, 13 Nov 2007, Manal Helal wrote:

> Date: Tue, 13 Nov 2007 12:44:16 +1100
> From: Manal Helal <[EMAIL PROTECTED]>
> To: r-help@r-project.org
> Subject: [R] connection diagram
> 
> Hi
>
> I am practically new to R, and need to construct connection diagrams,
> I have a table of data, of nodes in vertical rows, and horizontally
> the number of outgoing connections to other nodes, and the indices of
> these nodes, each in a column, so some columns are used, and some are
> not, based on how many connections I have
>
> the node is identified by these variables (dimension, wave number,
> partition number, index)
>
> the number of incoming and outgoing connections to each nodes varies,
> but obviously there is a relationship
>
> First: I need to draw diagrams of these connections
> Second: I need to apply regression analysis on this data, to predict a
> closed formula of how the 4 variables above decide how many incoming
> connections, and outgoing connections, and from/to which node(s)
>
> Am I making sense? if so, is this doable in R? or do I need to use
> other tools? If R can do it, I really need to find a tutorial or a
> starting link that I can follow to learn more how I can do these,
>
> sorry for being that ignorant about R, but I think I will need it a
> lot if it does what I need,
>
> thank you very much for your help 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] connection diagram

[Manal Helal]

Hi

I am practically new to R, and need to construct connection diagrams,
I have a table of data, of nodes in vertical rows, and horizontally
the number of outgoing connections to other nodes, and the indices of
these nodes, each in a column, so some columns are used, and some are
not, based on how many connections I have

the node is identified by these variables (dimension, wave number,
partition number, index)

the number of incoming and outgoing connections to each nodes varies,
but obviously there is a relationship

First: I need to draw diagrams of these connections
Second: I need to apply regression analysis on this data, to predict a
closed formula of how the 4 variables above decide how many incoming
connections, and outgoing connections, and from/to which node(s)

Am I making sense? if so, is this doable in R? or do I need to use
other tools? If R can do it, I really need to find a tutorial or a
starting link that I can follow to learn more how I can do these,

sorry for being that ignorant about R, but I think I will need it a
lot if it does what I need,

thank you very much for your help in advance,


-- 
Kind Regards,

Manal Helal

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

Re: [R] The "reverse" of do.call("rbind",) on a list

[Sebastian P. Luque]

On Tue, 13 Nov 2007 02:07:39 +0100,
Søren Højsgaard <[EMAIL PROTECTED]> wrote:

> Dear List, I want to turn the matrix
>> xm
>  [,1] [,2] [1,] "a" "b" [2,] "d" "e"

> into a list "by rows" as: [[1]] [1] "a" "b" [[2]] [1] "d" "e"

> A (bad?) way of doing this is as
>> unlist(apply(xm,1, list), recursive=F)
> [[1]] [1] "a" "b" [[2]] [1] "d" "e"

> - but there must be a more elegant way. Can anyone help on this?

How about:

lapply(seq(nrow(xm)), function(k) xm[k, ])

?

-- 
Seb

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

Re: [R] The "reverse" of do.call("rbind",) on a list

[Christos Hatzis]

Try split:

> x <- matrix(letters[1:4], 2)
> x
 [,1] [,2]
[1,] "a"  "c" 
[2,] "b"  "d" 

> split(x, row(x))
$`1`
[1] "a" "c"

$`2`
[1] "b" "d"
 
-Christos

> -Original Message-
> From: [EMAIL PROTECTED] 
> [mailto:[EMAIL PROTECTED] On Behalf Of Søren Højsgaard
> Sent: Monday, November 12, 2007 8:08 PM
> To: [EMAIL PROTECTED]
> Subject: [R] The "reverse" of do.call("rbind",) on a list
> 
> Dear List,
> I want to turn the matrix
> > xm
>  [,1] [,2]
> [1,] "a"  "b" 
> [2,] "d"  "e" 
> 
> into a list "by rows" as:
> [[1]]
> [1] "a" "b"
> [[2]]
> [1] "d" "e"
> 
> A (bad?) way of doing this is as
> > unlist(apply(xm,1, list), recursive=F)
> [[1]]
> [1] "a" "b"
> [[2]]
> [1] "d" "e"
>  
> - but there must be a more elegant way. Can anyone help on 
> this? Thanks!
> Søren
>  
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide 
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 

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

Re: [R] The "reverse" of do.call("rbind",) on a list

[Gabor Grothendieck]

Try:

split(xm, row(xm))

On Nov 12, 2007 8:07 PM, Søren Højsgaard <[EMAIL PROTECTED]> wrote:
> Dear List,
> I want to turn the matrix
> > xm
> [,1] [,2]
> [1,] "a"  "b"
> [2,] "d"  "e"
>
> into a list "by rows" as:
> [[1]]
> [1] "a" "b"
> [[2]]
> [1] "d" "e"
>
> A (bad?) way of doing this is as
> > unlist(apply(xm,1, list), recursive=F)
> [[1]]
> [1] "a" "b"
> [[2]]
> [1] "d" "e"
>
> - but there must be a more elegant way. Can anyone help on this? Thanks!
> Søren
>
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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

Re: [R] Resample a matrix (table)

[Moshe Olshansky]

Hi Luca,

One (not very elegant) way of doing this is to use
double loop as below:
Let A be your matrix (100x50), W your weights matrix
(10x10) and you want to produce matrix B (10X5)
You can do:

B <- matrix(0,nrow=10,ncol=5)
for (i in 1:10) for (j in 1:5) B[i,j] <-
sum(W*A[(10*(i-1)+1):10*i,(10*(j-1)+1):10*j])

Regards,

Moshe.


--- Luca Penasa <[EMAIL PROTECTED]> wrote:

> I have a matrix of values... for exemple the matrix
> could be 100 rows
> and 50 columns...
> I would like to resample the matrix to obtain a new
> matrix -> a 10x5 matrix
> So i need to calculate a new value for every cell, i
> would like to use a
> weightened mean to do this. I explain this thing
> better:
> In the exemple the original matrix would be divided
> into 50
> partial-matrix, for every partial-matrix should be
> extracted a value
> obtained from a weightened mean: the cell at center
> of the
> partial-matrix should have the max importance into
> the mean, the fhurter
> cell should have the lower importance (something
> like a gaussian shape,
> where i can decide the sigma used) In this way
> i'll obtain a
> resampled matrix with the values calculated using a
> weightened mean of
> the single partial-matrix...
> 
> Is there anybody who can suggest me a way to do
> this??? is there a
> function that can do this?
> 
> PS: sorry for my english... and... im a newbie,
> please try to explain in
> a simple way.. thank you
> 
> Luca Penasa,
> geology student,
> University of Padua, Italy
> 
> 
>  
>  
>  --
>  Email.it, the professional e-mail, gratis per te:
> http://www.email.it/f
>  
>  Sponsor:
>  In REGALO 'All the Good Thing' di NELLY FURTADO
>  Clicca qui:
>
http://adv.email.it/cgi-bin/foclick.cgi?mid=6617&d=12-11
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] The "reverse" of do.call("rbind",) on a list

[Søren Højsgaard]

Dear List,
I want to turn the matrix
> xm
 [,1] [,2]
[1,] "a"  "b" 
[2,] "d"  "e" 

into a list "by rows" as:
[[1]]
[1] "a" "b"
[[2]]
[1] "d" "e"

A (bad?) way of doing this is as
> unlist(apply(xm,1, list), recursive=F)
[[1]]
[1] "a" "b"
[[2]]
[1] "d" "e"
 
- but there must be a more elegant way. Can anyone help on this? Thanks!
Søren
 

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

Re: [R] update matrix with subset of it where only row names match

[jim holtman]

Here is one way of doing it that uses the row and column names:

> # create test data
> mat1 <- matrix(0, nrow=10, ncol=3)
> dimnames(mat1) <- list(paste('row', 1:10, sep=''), LETTERS[1:3])
> mat2 <- matrix(1:3, ncol=1, dimnames=list(c('row3', 'row7', 'row5'), "B"))
> mat2
 B
row3 1
row7 2
row5 3
> # create indexing matrix
> indx <- cbind(match(rownames(mat2), rownames(mat1)), match(colnames(mat2), 
> colnames(mat1)))
> indx
 [,1] [,2]
[1,]32
[2,]72
[3,]52
> mat1[indx] <- mat2
> mat1
  A B C
row1  0 0 0
row2  0 0 0
row3  0 1 0
row4  0 0 0
row5  0 3 0
row6  0 0 0
row7  0 2 0
row8  0 0 0
row9  0 0 0
row10 0 0 0


On Nov 12, 2007 4:54 PM, Martin Waller <[EMAIL PROTECTED]> wrote:
> I guess this has a simple solution:
>
> I have matrix 'mat1' which has row and column names, e.g.:
>
>A   B   C
> row10   0   0
> row20   0   0
> 
> rown0   0   0
>
> I have a another matrix 'mat2', essentially a subset of 'mat1' where the
> rownames are all in 'mat1' e.g.:
>
>B
> row35
> row86
> row54   7
>
>
> I want to insert the values of matrix mat2 for column B (in reality it
> could be some or all of column names A, B or C, etc.) (same name in both
> matrices if that matters - rownames of mat2 guaranteed to be in mat1)
> into matrix mat1 where the rownames match, so final desired result is:
>
> matrix mat1:
>A   B   C
> row10   0   0
> row20   0   0
> row30   5   0
> ...
> row80   6   0
> ...
> row54   0   7   0
> ..
> rown0   0   0
>
> My solution was (along the lines of):
>
> mat1[rownames(mat2)%in%rownames(mat1),"B"]=mat2[,"B"]
>
> Is there a better way? It doesn't 'feel' right?
>
> Thanks - hope I explained it right (its late and I had a little drink
> about an hour ago,etc).
>
>
> Martin
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
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] Tick mark puzzle.

[Bert Gunter]

Well, let me have a crack at it ... (inline below) 


Bert Gunter
Genentech Nonclinical Statistics


-Original Message-
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] On
Behalf Of Rolf Turner
Sent: Monday, November 12, 2007 2:16 PM
To: r-help@r-project.org
Subject: [R] Tick mark puzzle.


I am bewildered by the behaviour of tickmarks as demonstrated by the
following code.  (What I'm trying to do is draw a single tick mark  
extending
from the axis all the way down to the tick label, which is two lines  
from
the axis to make sure it doesn't overlap with the ``ordinary'' tick  
labels.)

# Try 1:
# Gap between tickmark and label.
plot(1:10)
axis(side = 1, at = 2.35, labels = "2.35", tick = TRUE,
  tcl=1,line=1,col="red", col.axis = "red")
readline('Go? ')

## From the first line position, below the plotting region (line =1) the
tick mark is
## drawn up 1 line height of text to the 0 position at the edge of the
plotting region
## The label is drawn below the line =1 position, leaving the gap

# Try 2:
# Gap between tickmark and axis.
plot(1:10)
axis(side = 1, at = 2.35, labels = "2.35", tick = TRUE,
  tcl=-1,line=1,col="red", col.axis = "red")
readline('Go? ')

## From the first line position one line height below the graph the tick
mark is drawn down to where
## the label is located, leaving the gap between tick mark and graph
boundary.

# Try 3:
# Gap between tickmark and label not filled.
plot(1:10)
axis(side = 1, at = 2.35, labels = "2.35", tick = TRUE,
tcl=1,line=1,col="red", col.axis = "red")
axis(side = 1, at = 2.35, labels ="", tick = TRUE, tcl=-1,col="red")

readline('Go? ')

## First axis command as in Try 1. Second command uses default line = 0 =
edge of plot region
## and so redraws tick the "other way" (down instead of up; try the second
axis it in a different ## color and you'll see it overplot the tick in the
different color)

# Try 4:
# Gap between tickmark and label ***is*** filled.
# Gives what I want, but ***why***?
plot(1:10)
axis(side = 1, at = 2.35, labels = "2.35", tick = TRUE,
  tcl=-1,line=1,col="red", col.axis = "red")
axis(side = 1, at = 2.35, labels ="", tick = TRUE, tcl=-1,col="red")

## Well you should be able to figure this one out now, right ... ;-)

Cheers,
Bert

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

Re: [R] how to assign a group mean to individual cases?

[Moshe Olshansky]

You can do the following:

x <- 1:10
g <- rep(3:5,len=10)
df <- data.frame(g=g,x=x)
y <- aggregate(df$x,list(df$g))
z <- sapply(df$g,function(x) which(y[,1]==x))
df1 <- data.frame(df,group.mean=y[z,2])

--- Casey Klofstad <[EMAIL PROTECTED]> wrote:

> I need advice on how to create a variable that is
> the group mean of
> another variable.
> 
> For example, I have a variable called x for which
> each row in the data
> set has a value. I also have a nominal variable
> called g that
> indicates which of 100 different groups each row
> belongs to.
> 
> So, I want to create a new variable called w, which
> is the group mean
> of x for which ever group the row belongs to (as
> indicated by variable
> g). Ideally, I'd also like to take out each row's
> value of x before
> calculating the group mean assigned to that row.
> 
> I've already tried the aggregate command. That gives
> me the group
> means, but does not allow me to assign them to each
> row in the data
> set.
> 
> THANKS!
> -- 
> Casey A. Klofstad
> University of Miami
> Department of Political Science
> Coral Gables, FL
> 
> [EMAIL PROTECTED]
> http://moya.bus.miami.edu/~cklofstad
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] time plotting problem

[John Kane]

Thanks to Gabor and Jim. I am not sure if the first
entry year = 2009 is all the problem I'm getting but
it is certainly seems like the worst of it.  

My stupidity:  Someone sent me the data set in Excel 
and I didn't do the basic data checks on. I _KNEW_ the
data went from 2006 to 2007. 

--- Gabor Grothendieck <[EMAIL PROTECTED]>
wrote:

> In your examples the first line of your data refers
> to the
> year 2009 and Oct 1st is repeated.  Is that really
> what
> you meant?
> 
> I can't tell what your problem is from your
> description
> other than the data problems cited but there are
> lots of
> examples of plotting with zoo in the following which
> may
> help you:
> vignette("zoo")
> vignette("zoo-quickref")
> ?plot.zoo
> ?xyplot.zoo
> 
> Note that zoo series must be time series, i.e. they
> must
> have unique times.
> 
> On Nov 12, 2007 1:47 PM, John Kane
> <[EMAIL PROTECTED]> wrote:
> > I am completely misunderstanding how to handle
> dates.
> > I want to plot a couple of data series against
> some
> > dates.  Simple example 1 below works fine.
> > Unfortunately I have multiple observations per day
> (no
> > time breakdowns) and observations across years.
> > (example 2 very simplistic version )
> >
> > Can anyone suggest a quick fix or point me to
> > something to read?  I thought that zoo might do it
> but
> > I seem to be missing something there too.
> >
> > Any suggestions gratefully recieved.
> >
> >
> > Example 1 consecutive dates same year.
> > =
> > x <- "days
> > 9/26/09
> > 9/27/06
> > 9/28/06
> > 9/29/06
> > 9/29/06
> > 9/29/06
> > 10/1/06
> > 10/1/06
> > 10/2/06
> > 10/3/06"
> >
> > mydata <- read.table(textConnection(x),
> header=TRUE,
> > as.is=TRUE); mydata
> >
> > mydates <- as.Date(mydata[,1], "%m/%d/%y");
> mydates
> > mynums <- rnorm(10)
> > plot(mydates, mynums)
> > 
> > Example 2 (things go blooy!)
> > non-consecutive dates different years.
> >
> > =
> > x <- "days
> > 9/26/09
> > 9/27/06
> > 9/28/06
> > 9/29/06
> > 9/29/06
> > 9/29/06
> > 10/1/07  # <- year changes
> > 10/1/07
> > 10/2/07
> > 10/3/07"
> >
> > mydata <- read.table(textConnection(x),
> header=TRUE,
> > as.is=TRUE); mydata
> >
> > mydates <- as.Date(mydata[,1], "%m/%d/%y");
> mydates
> > mynums <- rnorm(10)
> > plot(mydates, mynums)
> >
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/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] time plotting problem

[Gabor Grothendieck]

In your examples the first line of your data refers to the
year 2009 and Oct 1st is repeated.  Is that really what
you meant?

I can't tell what your problem is from your description
other than the data problems cited but there are lots of
examples of plotting with zoo in the following which may
help you:
vignette("zoo")
vignette("zoo-quickref")
?plot.zoo
?xyplot.zoo

Note that zoo series must be time series, i.e. they must
have unique times.

On Nov 12, 2007 1:47 PM, John Kane <[EMAIL PROTECTED]> wrote:
> I am completely misunderstanding how to handle dates.
> I want to plot a couple of data series against some
> dates.  Simple example 1 below works fine.
> Unfortunately I have multiple observations per day (no
> time breakdowns) and observations across years.
> (example 2 very simplistic version )
>
> Can anyone suggest a quick fix or point me to
> something to read?  I thought that zoo might do it but
> I seem to be missing something there too.
>
> Any suggestions gratefully recieved.
>
>
> Example 1 consecutive dates same year.
> =
> x <- "days
> 9/26/09
> 9/27/06
> 9/28/06
> 9/29/06
> 9/29/06
> 9/29/06
> 10/1/06
> 10/1/06
> 10/2/06
> 10/3/06"
>
> mydata <- read.table(textConnection(x), header=TRUE,
> as.is=TRUE); mydata
>
> mydates <- as.Date(mydata[,1], "%m/%d/%y"); mydates
> mynums <- rnorm(10)
> plot(mydates, mynums)
> 
> Example 2 (things go blooy!)
> non-consecutive dates different years.
>
> =
> x <- "days
> 9/26/09
> 9/27/06
> 9/28/06
> 9/29/06
> 9/29/06
> 9/29/06
> 10/1/07  # <- year changes
> 10/1/07
> 10/2/07
> 10/3/07"
>
> mydata <- read.table(textConnection(x), header=TRUE,
> as.is=TRUE); mydata
>
> mydates <- as.Date(mydata[,1], "%m/%d/%y"); mydates
> mynums <- rnorm(10)
> plot(mydates, mynums)
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] how to assign a group mean to individual cases?

[Gabor Grothendieck]

Look at:

?ave

On Nov 12, 2007 5:01 PM, Casey Klofstad <[EMAIL PROTECTED]> wrote:
> I need advice on how to create a variable that is the group mean of
> another variable.
>
> For example, I have a variable called x for which each row in the data
> set has a value. I also have a nominal variable called g that
> indicates which of 100 different groups each row belongs to.
>
> So, I want to create a new variable called w, which is the group mean
> of x for which ever group the row belongs to (as indicated by variable
> g). Ideally, I'd also like to take out each row's value of x before
> calculating the group mean assigned to that row.
>
> I've already tried the aggregate command. That gives me the group
> means, but does not allow me to assign them to each row in the data
> set.
>
> THANKS!
> --
> Casey A. Klofstad
> University of Miami
> Department of Political Science
> Coral Gables, FL
>
> [EMAIL PROTECTED]
> http://moya.bus.miami.edu/~cklofstad
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] strange `nls' behaviour

[Duncan Murdoch]

On 12/11/2007 2:56 PM, Joerg van den Hoff wrote:
> On Mon, Nov 12, 2007 at 11:09:21AM -0500, Duncan Murdoch wrote:
>> On 11/12/2007 9:14 AM, Joerg van den Hoff wrote:
>>> On Mon, Nov 12, 2007 at 07:36:34AM -0500, Duncan Murdoch wrote:
 On 11/12/2007 6:51 AM, Joerg van den Hoff wrote:
> I initially thought, this should better be posted to r-devel
> but alas! no response. 
 I think the reason there was no response is that your example is too 
 complicated.  You're doing a lot of strange things (fitfunc as a result 
 of deriv, using as.name, as.call, as.formula, etc.)  You should simplify 
>>> thanks for the feedback.
>>>
>>> concerning  "lot  of  strange  things":  OK.  I  thought the
>>> context might be important ("why, for heaven's sake  do  you
>>> want  to  do  this!?"), but, then, maybe not. so the easiest
>>> way to trigger a similar (not the  identical)  behaviour  is
>>> something like
>>>
>>> f <- function (n) {
>>>   #-
>>>   #define n data points for a (hardcoded) model:
>>>   #---
>>>   x <- seq(0, 5, length  = n)
>>>   y <- 2 * exp(-1*x) + 2; 
>>>   y <- rnorm(y,y, 0.01*y)
>>>
>>>   #the model (same as the above hardcoded one):
>>>   model <- y ~ a * exp (-b*x) + c
>>>
>>>   #call nls as usual:
>>>   res1 <- try(nls(model, start = c(a=2, b=1, c=2)))
>>>
>>>   #call it a bit differently:
>>>   res2 <- nls(y ~ eval(model[[3]]), start = c(a=2, b=1, c=2))
>>>
>>>   list(res1 = res1, res2 = res2)
>>>   #-
>>> }
>> I'd say the problem is relying on the default for the envir parameter to 
>> eval.  It's reasonable to expect nls to set things up so that terms in 
>> the model formula are taken from the right place, but when your model 
>> formula is playing with evaluation, you should be very careful to make 
>> sure the evaluation takes place in the right context.
> 
> agreed.
> 
>> The default for envir is "parent.frame()", and that can't be right: that 
>> will see local variables in whatever function called it, so if one of 
>> them has a variable called "model", you'll subscript it.
> 
> for one, in my actual application, I _do_ specify envir (and
> remember the above was really  not  the  original  situation
> where the subsetting of `model' is actually done only in the
> `deriv' call).   second, if I'm not missing  something,  doing
> the  evaluation  in  parent.frame()  is  fine  in  the above
> example and does not explain  the  error,  whose  triggering
> solely depends on the number of data points used.

You are missing something.  parent.frame() will be evaluated by eval(), 
so it refers to whatever function calls eval().  You don't know what 
function that will be, because it's buried deep within nls().

Perhaps you think that the formula

y ~ eval(model[[3]])

is the same as the original one?  It's not.  Try printing it:

 > y ~ eval(model[[3]])
y ~ eval(model[[3]])

The eval doesn't get called at this point, it gets called for every step 
of the least squares minimization.  Who knows what parent.frame() means 
then?

>> If I were trying to do what you're doing, I would construct the formula 
>> before the call to nls, and pass that.  I.e. something like the 
>> following silly code:
>>
>> model2 <- model
>> model2[[3]] <- model[[3]] # The eval() is implicit
>> res2 <- nls(model2, start = c(a=2, b=1, c=2))
> 
> 
> I  was  trying that (more or less) in my original example, 
> I think, but I will reexamine my application and see, whether I can  bypass
> the problem somehow along these lines.
> 
>> If you really want to put eval() in a formula, then I can't see how you 
>> could avoid an explicit specification of the envir parameter.  So I'd 
>> call this a bug in your code.
> 
> 
> as I said, this seems right in principle (still, if the call
> does happen at some well  defined  place  such  as  a  small
> function local to a user visible one, the eval without envir
> might be quite OK) but not w.r.t. explaining the nls  crash.

No, it's not okay.  Your formula is looking in places it shouldn't.  It 
makes absolutely no sense for your formula to depend on what function 
evaluates it.  It should only depend on the data that is available in 
the data argument to nls() or visible in the current environment.

> katherine  mullen  seems  to  have located the exact spot in
> `nls' where something goes wrong. so, for now I think martin
> might be right ("bug in my thinking" by assuming I'm allowed
> to use eval in this place), but otherwise it is  a  property
> (or bug) of `nls'.

She may have discovered where the error call was triggered, but things 
went wrong in this example when you mis-used eval().  If you think your 
original code used it properly, then simplify that down to a minimal 
example that triggers the "bug".  Don't put things in the example unless 
they are necessary to trigger the error.  For example:

  - drop the wrapper function f.
  - f

Re: [R] referring to a data-frame column using a variable

[Don MacQueen]

Or

   dat[,c1]

-Don

At 8:49 AM +1300 11/13/07, Rolf Turner wrote:
>On 13/11/2007, at 8:42 AM, Bernd Jagla wrote:
>
>>  Hi,
>>
>>
>>
>>  I would like to refer to a column in a data frame using a variable. 
>>  How
>>  would I do this?
>>
>>
>>
>>  Something like:
>>
>>
>>
>>  c1<-"column1"
>>
>>  dat${c1}
>>
>>
>>
>>  where names(dat) includes 'column1'
>
>   Use dat[[c1]].
>
>   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.


-- 
--
Don MacQueen
Environmental Protection Department
Lawrence Livermore National Laboratory
Livermore, CA, USA
925-423-1062

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] Query an Access database based on a date attribute

[Tudor Bodea]

Dear useRs,

I would like to query an Access database through R based on a date attribute
but, unfortunately, I fail to do so. For example, the table test_table of the
test.mdb looks like:

   ID cd  market competitor dd price
1   1 2007-11-20 atl-bos  delta 2007-11-20   210
2   2 2007-11-20 atl-bos  delta 2007-11-21   190
3   3 2007-11-20 atl-bos  delta 2007-11-22   180
4   4 2007-11-20 atl-bos united 2007-11-20   205
5   5 2007-11-20 atl-bos united 2007-11-21   195
6   6 2007-11-20 atl-bos united 2007-11-22   175
7   7 2007-11-20 atl-sfa  delta 2007-11-20   350
8   8 2007-11-20 atl-sfa  delta 2007-11-21   320
9   9 2007-11-20 atl-sfa  delta 2007-11-22   300
10 10 2007-11-20 atl-sfa aa 2007-11-20   340
11 11 2007-11-20 atl-sfa aa 2007-11-21   320
12 12 2007-11-20 atl-sfa aa 2007-11-22   300

In this context, I try to get all the records for which market is atl-bos,
competitor is delta and dd is 2007-11-20 (first record above). To do this I used

># channel <- odbcConnectAccess("test.mdb")
>res <- sqlQuery(channel, "select * from test_table where market = 'atl-bos'
and competitor = 'delta' and dd = 2007-11-20")

but the result seems to be incorrect.

> res
[1] ID cd market competitor dd price
<0 rows> (or 0-length row.names)

Does any of you know what it is going on? I suspect that the last expression in
the select statement is not executed/interpreted correctly and, as a
consequence, the query leads to the above empty set.

Thank you.

Tudor

--
Tudor Dan Bodea
Georgia Institute of Technology
School of Civil and Environmental Engineering
Web: http://www.prism.gatech.edu/~gtg757i

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

[Jenifer Larson-Hall]

I figured this out so I wanted to post a response to my own question. Thanks to 
Brian Ripley for his hint, which got me looking in the right direction. The 
plot I had seen before that I wanted was automatic output to the cor.plot 
command in the mvoutlier library. A similar tolerance ellipse can be obtained 
from the covPlot command in the robustbase library, using 
which="tolEllipsePlot" and classic=T.

Dr. Jenifer Larson-Hall
Assistant Professor of Linguistics
University of North Texas
(940)369-8950

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

Re: [R] State-space model estimation with EM

[Giovanni Petris]

Package "dlm" has a function for maximum likelihood estimation of
parameters in general (linear Normal) state space model. The function,
dlmMLE, computes the likelihood based on singular value decomposition
and appears to be fairly robust. 

No EM algorithm, though. 


Giovanni

> Date: Sun, 11 Nov 2007 18:32:32 + (GMT)
> From: [EMAIL PROTECTED]
> Sender: [EMAIL PROTECTED]
> Priority: normal
> Precedence: list
> 
> Hi - I follow some references and now implement my own state-space model 
> estimation. I have a question. In case, my equations are like this:
> 
> y(t) = Ax(t)+Bu(t)+eps(t) # observation eq
> x(t) = Cx(t-1)+Du(t)+eta(t) # state eq
> 
> Using EM, after backward recursion, you will use the smoothed state 
> estimation to update the A, B, C, and D which is chosen so as to maximize the 
> expectation equation. But for example, if my C is always of matrix zero 
> (model specification), during EM I still get the value of estimated C which 
> turns out to be non-zero. How can I resolve this conflict? Or I just ignore 
> my estimation result and keep it as zero? THank you.
> 
> - adschai
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 

-- 

Giovanni Petris  <[EMAIL PROTECTED]>
Associate Professor
Department of Mathematical Sciences
University of Arkansas - Fayetteville, AR 72701
Ph: (479) 575-6324, 575-8630 (fax)
http://definetti.uark.edu/~gpetris/

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

Re: [R] Multivariate time series

[Giovanni Petris]

Giusy,

There is also a package "dlm" that may be useful, but you need to
specify the model you want to use. 

Giovanni

> Date: Sun, 11 Nov 2007 08:40:42 -0800 (PST)
> From: Giusy <[EMAIL PROTECTED]>
> Sender: [EMAIL PROTECTED]
> Precedence: list
> 
> 
> Hello to everyone!
> I have a question for you..I need to predict multivariate time series, for
> example sales of 2 products related one to the other, having the 2 prices
> like inputs..
> Is there in R a function to do it? I saw dse package but I didn't find what
> a I'm looking for..
> Could anyone help me?
> Thank you very much
> Giusy
> -- 
> View this message in context: 
> http://www.nabble.com/Multivariate-time-series-tf4786472.html#a13693139
> 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.
> 
> 

-- 

Giovanni Petris  <[EMAIL PROTECTED]>
Associate Professor
Department of Mathematical Sciences
University of Arkansas - Fayetteville, AR 72701
Ph: (479) 575-6324, 575-8630 (fax)
http://definetti.uark.edu/~gpetris/

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

Re: [R] Microsoft SOAP - Help!!

[Duncan Temple Lang]

-BEGIN PGP SIGNED MESSAGE-
Hash: SHA1

Hi Konstantinos

 I am not able to reproduce the error. I can use different
search queries of different length, etc.

  The error is coming from libcurl and, if I have gotten things the
right way around, it basically means that the write handler callback
didn't handle the entire collection of text it was given.
Can you

curlPerform(url="http://soap.search.msn.com/webservices.asmx?wsdl";,
 httpheader=c(Accept="text/xml", Accept="multipart/*",
  'Content-Type' = "text/xml; charset=utf-8"),
 postfields=body,
 writefunction = function(txt) { print(txt) ; nchar(txt)},
 verbose = FALSE
)


and verify that the we do get the point where we are processing the
return from the server (and not failing in writing our request to the
server).

If the error is in handling the response from the server, then it may
be that we have to look at the string encodings, etc.


 D.

Metaxab wrote:
> Hello,
> 
> I am trying to access Microsoft Live Search Using SOAP through R.
> In R I am using the RCurl packages to make the calls.
> 
> I have the following situation that looks crazy and cannot figure out how to
> solve it:
> 
> #SOAP Request
> library(RCurl)
> 
> h = basicTextGatherer()
> 
> body='
>  SOAP-ENV:encodingStyle=" http://schemas.xmlsoap.org/soap/encoding/";
> xmlns:SOAP-ENV="http://schemas.xmlsoap.org/soap/envelope/";
> xmlns:xsd=" http://www.w3.org/2001/XMLSchema";
> xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance";
> xmlns:SOAP-ENC=" http://schemas.xmlsoap.org/soap/encoding/";>
> http://tempuri.org";>
> MY ID
> "Cancer"
> en-US
> Off
> 
> 
> Web
> 0
> 10
> Title
> 
> '
> 
> curlPerform(url="http://soap.search.msn.com/webservices.asmx?wsdl";,
> httpheader=c(Accept="text/xml", Accept="multipart/*",
>  'Content-Type' = "text/xml; charset=utf-8"),
> postfields=body,
> writefunction = h$update,
> verbose = FALSE
>   )
> 
> ##Put the responce in body.
> body = h$value()
> h$reset()
> 
> An now the crazy stuff. This code works well for some keywords but not for
> all.
> 
> So for the Query:"Cancer" the code works fine, but for the Query: "Tiger" I
> get the following error message: 
>   Error in curlPerform(url =
> "http://soap.search.msn.com/webservices.asmx?wsdl";,  : 
>   Failed writing body
> 
> I cannot figure that out. Please help
> 
> Thanks 
> Konstantinos Vougas
> 
-BEGIN PGP SIGNATURE-
Version: GnuPG v1.4.7 (Darwin)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iD8DBQFHOM7/9p/Jzwa2QP4RAvSGAJ91fdow7OPuihgE9yX+X75BJGiOzgCfapaw
0dni/l4DJjcAkqsEPcRk84E=
=10Tt
-END PGP SIGNATURE-

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

Re: [R] Sending E-Mail from R

[Brahm, David]

Josh Kalish <[EMAIL PROTECTED]> wrote:

> Has anyone used R to send out an email via an SMTP server?

If your system is Unix or Linux and you use a system call to "mail",
you might be interested to know that "uuencode" is a good way to include
attachments.  Create a temp file containing your body text, then append
the output of "uuencode" to the end of it.  Here's a piece of code:

temp <- tempfile("Smail.")
cat(text, file=temp, sep="\n")
for (i in seq(along=attach))
  system(paste("uuencode", attach[i], name[i], ">>", temp))
recip <- paste(recip, collapse=" ")
system(paste("mail -s '", subject, "' ", recip, "<", temp, sep=""),
TRUE)
unlink(temp)

-- David Brahm ([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] how to assign a group mean to individual cases?

[Casey Klofstad]

I need advice on how to create a variable that is the group mean of
another variable.

For example, I have a variable called x for which each row in the data
set has a value. I also have a nominal variable called g that
indicates which of 100 different groups each row belongs to.

So, I want to create a new variable called w, which is the group mean
of x for which ever group the row belongs to (as indicated by variable
g). Ideally, I'd also like to take out each row's value of x before
calculating the group mean assigned to that row.

I've already tried the aggregate command. That gives me the group
means, but does not allow me to assign them to each row in the data
set.

THANKS!
-- 
Casey A. Klofstad
University of Miami
Department of Political Science
Coral Gables, FL

[EMAIL PROTECTED]
http://moya.bus.miami.edu/~cklofstad

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] update matrix with subset of it where only row names match

[Martin Waller]

I guess this has a simple solution:

I have matrix 'mat1' which has row and column names, e.g.:

A   B   C
row10   0   0
row20   0   0

rown0   0   0

I have a another matrix 'mat2', essentially a subset of 'mat1' where the 
rownames are all in 'mat1' e.g.:

B
row35
row86
row54   7


I want to insert the values of matrix mat2 for column B (in reality it 
could be some or all of column names A, B or C, etc.) (same name in both 
matrices if that matters - rownames of mat2 guaranteed to be in mat1) 
into matrix mat1 where the rownames match, so final desired result is:

matrix mat1:
A   B   C
row10   0   0
row20   0   0
row30   5   0
...
row80   6   0
...
row54   0   7   0
..
rown0   0   0

My solution was (along the lines of):

mat1[rownames(mat2)%in%rownames(mat1),"B"]=mat2[,"B"]

Is there a better way? It doesn't 'feel' right?

Thanks - hope I explained it right (its late and I had a little drink 
about an hour ago,etc).


Martin

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

Re: [R] Using lme (nlme) to find the conditional variance of therandom effects

[Doran, Harold]

No, don't reach into the bVar slot. Use the proper extractor function
ranef() with postVar=T. There is no similar function for lme() 

> -Original Message-
> From: [EMAIL PROTECTED] 
> [mailto:[EMAIL PROTECTED] On Behalf Of Rick Bilonick
> Sent: Monday, November 12, 2007 4:40 PM
> To: R Help
> Subject: [R] Using lme (nlme) to find the conditional 
> variance of therandom effects
> 
> Using lmer in the lme4 package, you can compute the 
> conditional variance-covariance matrix of the random effects 
> using the bVar slot:
> 
> bVar: A list of the diagonal inner blocks (upper triangles 
> only) of the positive-definite matrices on the diagonal of 
> the inverse of ZtZ+Omega.
> With the appropriate scale factor (and conversion to a 
> symmetric matrix) these are the conditional 
> variance-covariance matrices of the random effects.
> 
> Is there anything similar in the nlme package using the lme function?
> 
> Rick B.
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] Using lme (nlme) to find the conditional variance of the random effects

[Rick Bilonick]

Using lmer in the lme4 package, you can compute the conditional
variance-covariance matrix of the random effects using the bVar slot:

bVar: A list of the diagonal inner blocks (upper triangles only) of the
positive-definite matrices on the diagonal of the inverse of ZtZ+Omega.
With the appropriate scale factor (and conversion to a symmetric matrix)
these are the conditional variance-covariance matrices of the random
effects.

Is there anything similar in the nlme package using the lme function?

Rick B.

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

Re: [R] How to set up your library location for multiple users

[James W. MacDonald]

Hi Dave,

Not sure what .lib.loc() is supposed to do, but you might want to look 
at section 6.2 of the R Installation and Administration manual, which 
covers library path issues.

Best,

Jim



Dave Jacoby wrote:
> I'm running R 2.5.0 (with the upgrade to 2.6.0 scheduled) on a and
> trying to write a few applications that will be run as myself, as my
> boss and as the web server, on one of many Solaris x86 servers.
> 
> I'd like to have my library sharable, and I thought I had found the way
> by putting this line early in my R scripts:
> 
> .lib.loc="/home/jacoby/My_R_Lib/"
> 
> But I tried that with my permissions and got this error:
> 
> Computing affinities[1] "Checking to see if your internet connection
> works..."
> Error in getCDF(cdfpackagename) : You do not have write access to
> /usr/local/R-2.5.0/lib/R/library
> Please check your permissions or provide a different 'lib' parameter.
> Execution halted
> 
> Any pointers?
> 

-- 
James W. MacDonald, M.S.
Biostatistician
Affymetrix and cDNA Microarray Core
University of Michigan Cancer Center
1500 E. Medical Center Drive
7410 CCGC
Ann Arbor MI 48109
734-647-5623

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

Re: [R] mtry in ctree_control()

[Achim Zeileis]

On Mon, 12 Nov 2007, Lam, Kelvin wrote:

> What is the actual usage of "mtry" in ctree(), or specifically,
> ctree_control() since it's a single tree?  Thanks in advance.

As the associated man page explains: This is useful if you want to build
an ensemble (such as a random forest) based on ctree()s. If you only want
to fit a single tree, you typically want to keep the default value "mtry =
0", i.e., no random selection.
Z

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

[R] Kronecker product in C file

[Gabriel Rodrigues Alves Margarido]

Hello, everyone.
I am writing a C code that uses quite a lot of matrix operations. I have
written my own functions for matrix-matrix multiplication, hadamard
(direct) and kronecker products.
However, I am sure using BLAS will greatly enhance speed. I did quite a
large search and found that I can use 'dgenn' for matrix-matrix
multiplication and 'dh' for hadamard product (will always use double
precision). On the other hand, I could not find a function to do the
kronecker product.

Any hints?

I tried to send this message ro R-devel, but always got a message saying
it was rejected by filter match. Don't know what's wrong...

Thank you in advance,
Gabriel Margarido.

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

[Lam, Kelvin]

Dear Group,

 

What is the actual usage of "mtry" in ctree(), or specifically,
ctree_control() since it's a single tree?  Thanks in advance.

 

Regards,

 

Kelvin Lam, MSc.

Analyst, Programming & Biostatistics

Institute for Clinical Evaluative Sciences (ICES)

2075 Bayview Avenue, G179 Toronto, ON M4N 3M5

(416) 480-4055 Ext. 3057 Fax: (416) 480-6048

email: [EMAIL PROTECTED]   

 


___
This email may contain confidential and/or privileged in...{{dropped:14}}

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

Re: [R] How to pool a group of samples and take the ave.

[affy snp]

Thanks Jim and Ben.

Allen

On Nov 12, 2007 3:43 PM, Benilton Carvalho <[EMAIL PROTECTED]> wrote:
> it's defined in Biobase.
>
>
> On Nov 12, 2007, at 3:30 PM, affy snp wrote:
>
> > Thanks Jim. It is way simple. Great.
> > Is there any function like rowMedians() which
> > could take the median value across samples?
> >
> > Allen
> >
> > On Nov 12, 2007 3:10 PM, jim holtman <[EMAIL PROTECTED]> wrote:
> >> Try something like this:
> >>
> >> myAvg <- rowMeans(A[,48:243])
> >> B <- A[1:47,] / myAvg
> >>
> >>
> >> On Nov 12, 2007 1:37 PM, affy snp <[EMAIL PROTECTED]> wrote:
> >>> Dear list,
> >>>
> >>> Hi! I have a table A, 238304 rows and 243 columns (representing
> >>> samples). First of all, I would like to pool a group of samples
> >>> from 48th column to 243rd column and take the average across
> >>> them and make a single column,saying as the reference column.
> >>>
> >>> Second, I want to use each column of first 47 columns in table
> >>> A divided by the reference column and end up with a new table
> >>> B with 238304 rows and 47 columns.
> >>>
> >>> Is there any simple code which especially could do sth like
> >>>   reference_column<-(A[,48]+A[,49]+...A[,243])/196
> >>> and B<-A[,1:47]/reference_column?
> >>>
> >>> Thank you very much for your help!
> >>>
> >>> Best,
> >>>  Allen
> >>>
> >>> __
> >>> R-help@r-project.org mailing list
> >>> https://stat.ethz.ch/mailman/listinfo/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.
>

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

Re: [R] How to pool a group of samples and take the ave.

[Benilton Carvalho]

it's defined in Biobase.

On Nov 12, 2007, at 3:30 PM, affy snp wrote:

> Thanks Jim. It is way simple. Great.
> Is there any function like rowMedians() which
> could take the median value across samples?
>
> Allen
>
> On Nov 12, 2007 3:10 PM, jim holtman <[EMAIL PROTECTED]> wrote:
>> Try something like this:
>>
>> myAvg <- rowMeans(A[,48:243])
>> B <- A[1:47,] / myAvg
>>
>>
>> On Nov 12, 2007 1:37 PM, affy snp <[EMAIL PROTECTED]> wrote:
>>> Dear list,
>>>
>>> Hi! I have a table A, 238304 rows and 243 columns (representing
>>> samples). First of all, I would like to pool a group of samples
>>> from 48th column to 243rd column and take the average across
>>> them and make a single column,saying as the reference column.
>>>
>>> Second, I want to use each column of first 47 columns in table
>>> A divided by the reference column and end up with a new table
>>> B with 238304 rows and 47 columns.
>>>
>>> Is there any simple code which especially could do sth like
>>>   reference_column<-(A[,48]+A[,49]+...A[,243])/196
>>> and B<-A[,1:47]/reference_column?
>>>
>>> Thank you very much for your help!
>>>
>>> Best,
>>>  Allen
>>>
>>> __
>>> R-help@r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/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.

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

Re: [R] How to pool a group of samples and take the ave.

[jim holtman]

you can use 'apply':

apply(A[, 48:243], 1, median)

On Nov 12, 2007 3:30 PM, affy snp <[EMAIL PROTECTED]> wrote:
> Thanks Jim. It is way simple. Great.
> Is there any function like rowMedians() which
> could take the median value across samples?
>
> Allen
>
> On Nov 12, 2007 3:10 PM, jim holtman <[EMAIL PROTECTED]> wrote:
>
> > Try something like this:
> >
> > myAvg <- rowMeans(A[,48:243])
> > B <- A[1:47,] / myAvg
> >
> >
> > On Nov 12, 2007 1:37 PM, affy snp <[EMAIL PROTECTED]> wrote:
> > > Dear list,
> > >
> > > Hi! I have a table A, 238304 rows and 243 columns (representing
> > > samples). First of all, I would like to pool a group of samples
> > > from 48th column to 243rd column and take the average across
> > > them and make a single column,saying as the reference column.
> > >
> > > Second, I want to use each column of first 47 columns in table
> > > A divided by the reference column and end up with a new table
> > > B with 238304 rows and 47 columns.
> > >
> > > Is there any simple code which especially could do sth like
> > >   reference_column<-(A[,48]+A[,49]+...A[,243])/196
> > > and B<-A[,1:47]/reference_column?
> > >
> > > Thank you very much for your help!
> > >
> > > Best,
> > >  Allen
> > >
> > > __
> > > R-help@r-project.org mailing list
> > > https://stat.ethz.ch/mailman/listinfo/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?
> >
>



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

[R] How to set up your library location for multiple users

[Dave Jacoby]

I'm running R 2.5.0 (with the upgrade to 2.6.0 scheduled) on a and
trying to write a few applications that will be run as myself, as my
boss and as the web server, on one of many Solaris x86 servers.

I'd like to have my library sharable, and I thought I had found the way
by putting this line early in my R scripts:

.lib.loc="/home/jacoby/My_R_Lib/"

But I tried that with my permissions and got this error:

Computing affinities[1] "Checking to see if your internet connection
works..."
Error in getCDF(cdfpackagename) : You do not have write access to
/usr/local/R-2.5.0/lib/R/library
Please check your permissions or provide a different 'lib' parameter.
Execution halted

Any pointers?

-- 
David Jacoby[EMAIL PROTECTED]
  Purdue Genomics Facility
  S049, WSLR building
  Phone: hah!

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

Re: [R] How to pool a group of samples and take the ave.

[affy snp]

Thanks Jim. It is way simple. Great.
Is there any function like rowMedians() which
could take the median value across samples?

Allen

On Nov 12, 2007 3:10 PM, jim holtman <[EMAIL PROTECTED]> wrote:
> Try something like this:
>
> myAvg <- rowMeans(A[,48:243])
> B <- A[1:47,] / myAvg
>
>
> On Nov 12, 2007 1:37 PM, affy snp <[EMAIL PROTECTED]> wrote:
> > Dear list,
> >
> > Hi! I have a table A, 238304 rows and 243 columns (representing
> > samples). First of all, I would like to pool a group of samples
> > from 48th column to 243rd column and take the average across
> > them and make a single column,saying as the reference column.
> >
> > Second, I want to use each column of first 47 columns in table
> > A divided by the reference column and end up with a new table
> > B with 238304 rows and 47 columns.
> >
> > Is there any simple code which especially could do sth like
> >   reference_column<-(A[,48]+A[,49]+...A[,243])/196
> > and B<-A[,1:47]/reference_column?
> >
> > Thank you very much for your help!
> >
> > Best,
> >  Allen
> >
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/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.

[R] a repetition of simulation

[sigalit mangut-leiba]

Hello,
I have a simple (?) simulation problem.
I'm doing a simulation with logistic model and I want to reapet it 600
times.
The simulation looks like this:

z <- 0
x <- 0
y <- 0
aps <- 0
tiss <- 0
for (i in 1:500){
z[i] <- rbinom(1, 1, .6)
x[i] <- rbinom(1, 1, .95)
y[i] <- z[i]*x[i]
if (y[i]==1) aps[i] <- rnorm(1,mean=13.4, sd=7.09) else aps[i] <-
rnorm(1,mean=12.67, sd=6.82)
if (y[i]==1) tiss[i] <- rnorm(1,mean=20.731,sd=9.751) else  tiss[i] <-
rnorm(1,mean=18.531,sd=9.499)
}
v <- data.frame(y, aps, tiss)
log_v <- glm(y~., family=binomial, data=v)
summary(log_v)

I want to do a repetition of this 600 times (I want to have 600 logistic
models), and see all the coefficients of the covariates aps & tiss.
Thanks in advance,
Sigalit.

[[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] time plotting problem

[jim holtman]

Now your first data point is 9/26/09; is it supposed to be 9/26/06?

On Nov 12, 2007 1:47 PM, John Kane <[EMAIL PROTECTED]> wrote:
> I am completely misunderstanding how to handle dates.
> I want to plot a couple of data series against some
> dates.  Simple example 1 below works fine.
> Unfortunately I have multiple observations per day (no
> time breakdowns) and observations across years.
> (example 2 very simplistic version )
>
> Can anyone suggest a quick fix or point me to
> something to read?  I thought that zoo might do it but
> I seem to be missing something there too.
>
> Any suggestions gratefully recieved.
>
>
> Example 1 consecutive dates same year.
> =
> x <- "days
> 9/26/09
> 9/27/06
> 9/28/06
> 9/29/06
> 9/29/06
> 9/29/06
> 10/1/06
> 10/1/06
> 10/2/06
> 10/3/06"
>
> mydata <- read.table(textConnection(x), header=TRUE,
> as.is=TRUE); mydata
>
> mydates <- as.Date(mydata[,1], "%m/%d/%y"); mydates
> mynums <- rnorm(10)
> plot(mydates, mynums)
> 
> Example 2 (things go blooy!)
> non-consecutive dates different years.
>
> =
> x <- "days
> 9/26/09
> 9/27/06
> 9/28/06
> 9/29/06
> 9/29/06
> 9/29/06
> 10/1/07  # <- year changes
> 10/1/07
> 10/2/07
> 10/3/07"
>
> mydata <- read.table(textConnection(x), header=TRUE,
> as.is=TRUE); mydata
>
> mydates <- as.Date(mydata[,1], "%m/%d/%y"); mydates
> mynums <- rnorm(10)
> plot(mydates, mynums)
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] time plotting problem

[jim holtman]

So what is wrong with the plots?  I used your example, and it appears
to plot the data correctly.  What were you expecting?

On Nov 12, 2007 1:47 PM, John Kane <[EMAIL PROTECTED]> wrote:
> I am completely misunderstanding how to handle dates.
> I want to plot a couple of data series against some
> dates.  Simple example 1 below works fine.
> Unfortunately I have multiple observations per day (no
> time breakdowns) and observations across years.
> (example 2 very simplistic version )
>
> Can anyone suggest a quick fix or point me to
> something to read?  I thought that zoo might do it but
> I seem to be missing something there too.
>
> Any suggestions gratefully recieved.
>
>
> Example 1 consecutive dates same year.
> =
> x <- "days
> 9/26/09
> 9/27/06
> 9/28/06
> 9/29/06
> 9/29/06
> 9/29/06
> 10/1/06
> 10/1/06
> 10/2/06
> 10/3/06"
>
> mydata <- read.table(textConnection(x), header=TRUE,
> as.is=TRUE); mydata
>
> mydates <- as.Date(mydata[,1], "%m/%d/%y"); mydates
> mynums <- rnorm(10)
> plot(mydates, mynums)
> 
> Example 2 (things go blooy!)
> non-consecutive dates different years.
>
> =
> x <- "days
> 9/26/09
> 9/27/06
> 9/28/06
> 9/29/06
> 9/29/06
> 9/29/06
> 10/1/07  # <- year changes
> 10/1/07
> 10/2/07
> 10/3/07"
>
> mydata <- read.table(textConnection(x), header=TRUE,
> as.is=TRUE); mydata
>
> mydates <- as.Date(mydata[,1], "%m/%d/%y"); mydates
> mynums <- rnorm(10)
> plot(mydates, mynums)
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] How can I fitting this function from values.

[Rolf Turner]

On 13/11/2007, at 8:27 AM, Horacio Castellini wrote:

> Hi all, I'd like fit this function:
>
> G(t)=A*1/(1+t/B)*1/sqrt(1+t/C)
>
> where A, B and C are fitting constants that one like search. I have
> got a fcs<-(t,G(t)) value table in a data frame.
>
> Any one can help me? Tahnks Horacio.


I ***thought*** that nls would solve your problem in a trice.
But I tried a toy example before replying to you, and it messed up
mightily:

foo <- function(t,ccc){
ccc[1]*1/(1+t/ccc[2])*1/sqrt(1+1/ccc[3])
}
ttt <- seq(1,10,length=100)
yyy <- foo(ttt,c(5,3,2))
set.seed(42)
dat.test <- data.frame(t=ttt,y=yyy+rnorm(100,0,0.1))
with(dat.test,plot(t,y))
with(dat.test,lines(t,foo(t,c(5,3,2 # Looks cool.
fit.test <- nls(y~A*1/(1+t/B)*1/sqrt(1+t/C),start=list 
(A=5,b=3,C=2),data=dat.test)
# Gives error; no idea what the implications are.

The error message:

Error in `[[<-.data.frame`(`*tmp*`, var, value = c(0.70934153524875,  
0.761288405463172,  :
   replacement has 65 rows, data has 100
In addition: Warning message:
In t/B : longer object length is not a multiple of shorter object length

I don't think I have *ever* had a good experience with nls(); it is  
always pretty flaky.
(Well, maybe it worked for me once ..., can't remember for sure. :-) )

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] How to pool a group of samples and take the ave.

[jim holtman]

Was supposed to be:

B <- A[,1:47] / myAvg

On Nov 12, 2007 3:10 PM, jim holtman <[EMAIL PROTECTED]> wrote:
> Try something like this:
>
> myAvg <- rowMeans(A[,48:243])
> B <- A[1:47,] / myAvg
>
>
> On Nov 12, 2007 1:37 PM, affy snp <[EMAIL PROTECTED]> wrote:
> > Dear list,
> >
> > Hi! I have a table A, 238304 rows and 243 columns (representing
> > samples). First of all, I would like to pool a group of samples
> > from 48th column to 243rd column and take the average across
> > them and make a single column,saying as the reference column.
> >
> > Second, I want to use each column of first 47 columns in table
> > A divided by the reference column and end up with a new table
> > B with 238304 rows and 47 columns.
> >
> > Is there any simple code which especially could do sth like
> >   reference_column<-(A[,48]+A[,49]+...A[,243])/196
> > and B<-A[,1:47]/reference_column?
> >
> > Thank you very much for your help!
> >
> > Best,
> >  Allen
> >
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/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?
>



-- 
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] How to pool a group of samples and take the ave.

[jim holtman]

Try something like this:

myAvg <- rowMeans(A[,48:243])
B <- A[1:47,] / myAvg

On Nov 12, 2007 1:37 PM, affy snp <[EMAIL PROTECTED]> wrote:
> Dear list,
>
> Hi! I have a table A, 238304 rows and 243 columns (representing
> samples). First of all, I would like to pool a group of samples
> from 48th column to 243rd column and take the average across
> them and make a single column,saying as the reference column.
>
> Second, I want to use each column of first 47 columns in table
> A divided by the reference column and end up with a new table
> B with 238304 rows and 47 columns.
>
> Is there any simple code which especially could do sth like
>   reference_column<-(A[,48]+A[,49]+...A[,243])/196
> and B<-A[,1:47]/reference_column?
>
> Thank you very much for your help!
>
> Best,
>  Allen
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] How can I fitting this function from values.

[Joerg van den Hoff]

On Mon, Nov 12, 2007 at 04:27:26PM -0300, Horacio Castellini wrote:
> Hi all, I'd like fit this function:
> 
> G(t)=A*1/(1+t/B)*1/sqrt(1+t/C)
> 
> where A, B and C are fitting constants that one like search. I have
> got a fcs<-(t,G(t)) value table in a data frame.
> 
> Any one can help me? Tahnks Horacio.

if your data are in `df' and the columns are labeled `t' and `G'
this seems what you want (cf. `?nls'):

nls(G ~ A*1/(1+t/B)*1/sqrt(1+t/C), start =list(A=1, B=2, C=3), data = df) 

where `start' should be set to reasonable start values.

and I would think, it is better to use p1 = 1/B, p2= 1/C during the fit
(no singularities of G for certain parameter values...).

joerg

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

[jim holtman]

What happens if you have multiple matches in the comparison between
the content_feat and ob_feat?  Why don't you just iterate through the
content_feat and use 'match' to find the corresponding match in
ob_feat?  This should speed it up.  Also why are you using 'as.matrix'
when the values in the 'if' statement are objects of size 1?  Are any
of the objects dataframes?  If so, convert them to matrices for
efficiency.

On Nov 12, 2007 12:09 PM, Mahmudul Haque <[EMAIL PROTECTED]> wrote:
> hi,
>
> please help me out in the following case. seems like it stuck in some 
> where(already 7 hrs passed). what I want is to combine 4 matrix in to one 
> matrix of desired length.
>
> final_matrix<-function(ob_feat,content_feat,link_feat,link_feat_transformed){
> complete_feat<-matrix(rep(-1),nrow=11402,ncol=278)
>
>
> for(i in 1:8944)
>{q<-c(0)
>for(j in 1:11402)
>{
>if(as.matrix(content_feat[i,2])==as.matrix(ob_feat[j,2]))
>{complete_feat[i,1]=as.matrix(ob_feat[j,2])
>complete_feat[i,2:97]=as.matrix(content_feat[i,3:98])
>complete_feat[i,98:99]=as.matrix(ob_feat[j,3:4])
>complete_feat[i,100:140]=as.matrix(link_feat[j,3:43])
>complete_feat[i,141:278]=as.matrix(link_feat_transformed[j,3:140])
>   q<-1}
>if(q==1)
>break;
>}
>
>
>}
>for (i in 8945:11402){
>complete_feat[i,1]=as.matrix(ob_feat[i,2])
>complete_feat[i,98:99]=as.matrix(ob_feat[i,3:4])
>complete_feat[i,100:140]=as.matrix(link_feat[i,3:43])
>complete_feat[i,141:278]=as.matrix(link_feat_transformed[i,3:140])
>}
>
> list(complete_feat=complete_feat)
> }
>
> kind regards,
>
> mahmudul haque
>
>  __
>
>
>
>[[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.
>



-- 
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] strange `nls' behaviour

[Joerg van den Hoff]

On Mon, Nov 12, 2007 at 11:09:21AM -0500, Duncan Murdoch wrote:
> On 11/12/2007 9:14 AM, Joerg van den Hoff wrote:
> >On Mon, Nov 12, 2007 at 07:36:34AM -0500, Duncan Murdoch wrote:
> >>On 11/12/2007 6:51 AM, Joerg van den Hoff wrote:
> >>>I initially thought, this should better be posted to r-devel
> >>>but alas! no response. 
> >>
> >>I think the reason there was no response is that your example is too 
> >>complicated.  You're doing a lot of strange things (fitfunc as a result 
> >>of deriv, using as.name, as.call, as.formula, etc.)  You should simplify 
> >
> >thanks for the feedback.
> >
> >concerning  "lot  of  strange  things":  OK.  I  thought the
> >context might be important ("why, for heaven's sake  do  you
> >want  to  do  this!?"), but, then, maybe not. so the easiest
> >way to trigger a similar (not the  identical)  behaviour  is
> >something like
> >
> >f <- function (n) {
> >   #-
> >   #define n data points for a (hardcoded) model:
> >   #---
> >   x <- seq(0, 5, length  = n)
> >   y <- 2 * exp(-1*x) + 2; 
> >   y <- rnorm(y,y, 0.01*y)
> >
> >   #the model (same as the above hardcoded one):
> >   model <- y ~ a * exp (-b*x) + c
> >
> >   #call nls as usual:
> >   res1 <- try(nls(model, start = c(a=2, b=1, c=2)))
> >
> >   #call it a bit differently:
> >   res2 <- nls(y ~ eval(model[[3]]), start = c(a=2, b=1, c=2))
> >
> >   list(res1 = res1, res2 = res2)
> >   #-
> >}
> 
> I'd say the problem is relying on the default for the envir parameter to 
> eval.  It's reasonable to expect nls to set things up so that terms in 
> the model formula are taken from the right place, but when your model 
> formula is playing with evaluation, you should be very careful to make 
> sure the evaluation takes place in the right context.

agreed.

> 
> The default for envir is "parent.frame()", and that can't be right: that 
> will see local variables in whatever function called it, so if one of 
> them has a variable called "model", you'll subscript it.

for one, in my actual application, I _do_ specify envir (and
remember the above was really  not  the  original  situation
where the subsetting of `model' is actually done only in the
`deriv' call). second, if I'm not missing  something,  doing
the  evaluation  in  parent.frame()  is  fine  in  the above
example and does not explain  the  error,  whose  triggering
solely depends on the number of data points used.

> 
> If I were trying to do what you're doing, I would construct the formula 
> before the call to nls, and pass that.  I.e. something like the 
> following silly code:
> 
> model2 <- model
> model2[[3]] <- model[[3]] # The eval() is implicit
> res2 <- nls(model2, start = c(a=2, b=1, c=2))


I  was  trying that (more or less) in my original example, 
I think, but I will reexamine my application and see, whether I can  bypass
the problem somehow along these lines.

> 
> If you really want to put eval() in a formula, then I can't see how you 
> could avoid an explicit specification of the envir parameter.  So I'd 
> call this a bug in your code.


as I said, this seems right in principle (still, if the call
does happen at some well  defined  place  such  as  a  small
function local to a user visible one, the eval without envir
might be quite OK) but not w.r.t. explaining the nls  crash.
katherine  mullen  seems  to  have located the exact spot in
`nls' where something goes wrong. so, for now I think martin
might be right ("bug in my thinking" by assuming I'm allowed
to use eval in this place), but otherwise it is  a  property
(or bug) of `nls'.

joerg van den hoff

> 
> Duncan Murdoch
> 
> 
> >this is without all the overhead of my first example and now
> >(since not quite the same) the problem  arises  at  n  ==  3
> >where  the  fit  can't  really  procede  (there  are  also 3
> >parameters -- the first example was more  relevant  in  this
> >respect)  but  anyway  the  second nls-call does not procede
> >beyond the parsing phase of `model.frame'.
> >
> >so,  I  can't  see  room for a real bug in my code, but very
> >well a chance that I misuse `nls'  (i.e.  not  understanding
> >what really is tolerable for the `model' argument of `nls').
> >but  if the latter is not the case, I would think there is a
> >bug in `nls'  (or,  actually,  rather  in  `model.frame'  or
> >whatever)  when  parsing  the  nls call.
> 
> 
> >
> >
> >>it down to isolate the bug.  Thats a lot of work, but you're the one in 
> >>the best position to do it.  I'd say there's at least an even chance 
> >>that the bug is in your code rather than in nls().
> >>
> >>And 2.5.0 *is* ancient; please confirm the bug exists in R-patched if it 
> >>turns out to be an R bug.
> >
> >if  need  be,  I'll  do  that  (if  I  get it compiled under
> >macosX). but assuming  that  you  have  R-patched  installed
> >anyway, I would appreciate if you would copy the new example
> >ab

Re: [R] Bug (?) in read.fwf

[david.jessop]

Jim,

That works - thank you; but still doesn't explain why stringsAsFactors
doesn't seem to have an effect. 

Regards,

David
 

-Original Message-
From: jim holtman [mailto:[EMAIL PROTECTED] 
Sent: 08 November 2007 13:41
To: Jessop, David
Cc: r-help@r-project.org
Subject: Re: [R] Bug (?) in read.fwf

have you tried  "as.is=TRUE"

On Nov 8, 2007 6:20 AM,  <[EMAIL PROTECTED]> wrote:
> Hi,
>
> I'm trying to use read.fwf
>
>temp = read.fwf ("Raw data.txt", widths = c (11, 21, 10, rep 
> (16, 6)) ,skip = 2, n = 2, stringsAsFactors = FALSE, strip.white = 
> TRUE)
>
> but no matter what I do the strings are turned into factors.  I 
> believe it's the "n=2" parameter that causes the problem as it seems 
> to work without this.  Am I missing something?
>
> Thanks in advance,
>
> David Jessop
>
>
> Issued by UBS AG or affiliates to professional investors for 
> information only and its accuracy/completeness is not guaranteed.
> All opinions may change without notice and may differ to 
> opinions/recommendations expressed by other business areas of UBS.
> UBS may maintain long/short positions and trade in instruments 
> referred to. Unless stated otherwise, this is not a personal 
> recommendation, offer or solicitation to buy/sell and any 
> prices/quotations are indicative only. UBS may provide investment 
> banking and other services to, and/or its employees may be directors 
> of, companies referred to. To the extent permitted by law, UBS does 
> not accept any liability arising from the use of this communication.
>
>  (c) 2007 UBS.  All rights reserved. Intended for recipient only and 
> not for further distribution without the consent of UBS.
>
> UBS Limited is a company registered in England & Wales under company 
> number 2035362, whose registered office is at 1 Finsbury Avenue, 
> London, EC2M 2PP, United Kingdom.
>
> UBS AG (London Branch) is registered as a branch of a foreign company 
> under number BR004507, whose registered office is at
> 1 Finsbury Avenue, London, EC2M 2PP, United Kingdom.
>
> UBS Clearing and Execution Services Limited is a company registered in

> England & Wales under company number 03123037, whose registered office

> is at 1 Finsbury Avenue, London, EC2M 2PP, United Kingdom.
>
>
>
>
>
>
>
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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?
Issued by UBS AG or affiliates to professional investors for
information only and its accuracy/completeness is not guaranteed.
All opinions may change without notice and may differ to
opinions/recommendations expressed by other business areas of UBS.
UBS may maintain long/short positions and trade in instruments
referred to. Unless stated otherwise, this is not a personal
recommendation, offer or solicitation to buy/sell and any
prices/quotations are indicative only. UBS may provide investment
banking and other services to, and/or its employees may be directors
of, companies referred to. To the extent permitted by law, UBS does
not accept any liability arising from the use of this communication.

 © 2007 UBS.  All rights reserved. Intended for recipient only
and not for further distribution without the consent of UBS.

UBS Limited is a company registered in England & Wales under company
number 2035362, whose registered office is at 1 Finsbury Avenue,
London, EC2M 2PP, United Kingdom.

UBS AG (London Branch) is registered as a branch of a foreign company
under number BR004507, whose registered office is at
1 Finsbury Avenue, London, EC2M 2PP, United Kingdom.

UBS Clearing and Execution Services Limited is a company registered
in England & Wales under company number 03123037, whose registered
office is at 1 Finsbury Avenue, London, EC2M 2PP, United Kingdom.

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

Re: [R] referring to a data-frame column using a variable

[Rolf Turner]

On 13/11/2007, at 8:42 AM, Bernd Jagla wrote:

> Hi,
>
>
>
> I would like to refer to a column in a data frame using a variable.  
> How
> would I do this?
>
>
>
> Something like:
>
>
>
> c1<-"column1"
>
> dat${c1}
>
>
>
> where names(dat) includes 'column1'

Use dat[[c1]].

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.

[R] referring to a data-frame column using a variable

[Bernd Jagla]

Hi, 

 

I would like to refer to a column in a data frame using a variable. How
would I do this?

 

Something like:

 

c1<-"column1"

dat${c1}

 

where names(dat) includes 'column1'

 

Thanks,

 

Bernd


[[alternative HTML version deleted]]

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

[R] How can I fitting this function from values.

[Horacio Castellini]

Hi all, I'd like fit this function:

G(t)=A*1/(1+t/B)*1/sqrt(1+t/C)

where A, B and C are fitting constants that one like search. I have
got a fcs<-(t,G(t)) value table in a data frame.

Any one can help me? Tahnks Horacio.

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
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 on drawing a tree with "ape"?

[Ben Bolker]

Hua Li wrote:
> 
> Dear Simon and everyone,
> 
> Thanks for response. 
> 
> Specifically, I have a tree with species "A, B, C" and
> can be written as
> "((A:45.15,C:45.15):46.19,B:91.34):0.0;" . If I use
> the R command
>> mytreeABC
> 
> [1] "((A:45.15,C:45.15):46.19,B:91.34):0.0;"
> 
>> plot(read.tree(text = mytreeABC))
> 
> --- Simon Blomberg <[EMAIL PROTECTED]> wrote:
> 
>> ?edgelabels
>> 
> 

   Did you look at the examples in the file suggested by Simon?
Unfortunately, the PS file but not the JPG came through in the
version I read.

  It seems to me you may want _node_ labels rather
than _edge_ labels?See below ...

mytreeABC = "((A:45.15,C:45.15):46.19,B:91.34):0.0;"
t1 = read.tree(text=mytreeABC)
library(ape)
plot(t1)
edgelabels(1:4)
nodelabels(c("a","b"))

-- 
View this message in context: 
http://www.nabble.com/help-on-drawing-a-tree-with-%22ape%22--tf4789546.html#a13711995
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] na's and continuous variables

[Tonker]

Hi 
I am performing a lmer on count data with several explanatory variables both
continuous and categorical. All go in fine apart from two of the continuous
variables that have na's in them. They give errors like:

Leading minor of order 5 in downdated X'X is not positive definite

I used the as.factor function to define my categorical factors as
categorical but when i ask is.factor about my two continuous variables that
i didnt use the as.factor for it replies yes.
I have so far come to conclusion that because the data contains na's R is
automatically reading it as categorical not continuous and i was wondering
if there is a way to overcome this?

Thanks
Tonker

ps i do use the na.action=na.omit in my model so if possible i would like to
not have to delete all the records with na in in these two explanatory
factors.
 
-- 
View this message in context: 
http://www.nabble.com/na%27s-and-continuous-variables-tf4792738.html#a13710596
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] time plotting problem

[John Kane]

I am completely misunderstanding how to handle dates. 
I want to plot a couple of data series against some
dates.  Simple example 1 below works fine. 
Unfortunately I have multiple observations per day (no
time breakdowns) and observations across years.
(example 2 very simplistic version )

Can anyone suggest a quick fix or point me to
something to read?  I thought that zoo might do it but
I seem to be missing something there too.  

Any suggestions gratefully recieved.


Example 1 consecutive dates same year.
=
x <- "days
9/26/09
9/27/06
9/28/06
9/29/06
9/29/06
9/29/06
10/1/06
10/1/06
10/2/06
10/3/06"

mydata <- read.table(textConnection(x), header=TRUE,
as.is=TRUE); mydata

mydates <- as.Date(mydata[,1], "%m/%d/%y"); mydates
mynums <- rnorm(10)
plot(mydates, mynums)

Example 2 (things go blooy!)
non-consecutive dates different years.

=
x <- "days
9/26/09
9/27/06
9/28/06
9/29/06
9/29/06
9/29/06
10/1/07  # <- year changes
10/1/07
10/2/07
10/3/07"

mydata <- read.table(textConnection(x), header=TRUE,
as.is=TRUE); mydata

mydates <- as.Date(mydata[,1], "%m/%d/%y"); mydates
mynums <- rnorm(10)
plot(mydates, mynums)

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

[Casey,Richard]

Hello,

We're installing R v.2.6 on HP Proliant 32-bit server and RedHat Enterprise 
Linux v.4 using the binaries and RPM's in CRAN in /bin/linux/redhat/el4/i386.

I could not find instructions on how to use the RPM's.  I assume you issue:

rpm -i R-2.6.0-3.rh4.i386.rpm
rpm -i libRmath-2.6.0-3.rh4.i386.rpm

Is that all there is to it?

I ask the question because we're having trouble installing Bioconductor v.2.1 
in R.  Following the Bioconductor install instructions, the gcc/g++ compilers 
complain about incompatible R libraries, as in:

/usr/bin/ld: skipping incompatible /usr/lib/R/lib/libR.so when searching for -lR

Just wondering if R is installed correctly for Bioconductor.



Richard

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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 pool a group of samples and take the ave.

[affy snp]

Dear list,

Hi! I have a table A, 238304 rows and 243 columns (representing
samples). First of all, I would like to pool a group of samples
from 48th column to 243rd column and take the average across
them and make a single column,saying as the reference column.

Second, I want to use each column of first 47 columns in table
A divided by the reference column and end up with a new table
B with 238304 rows and 47 columns.

Is there any simple code which especially could do sth like
   reference_column<-(A[,48]+A[,49]+...A[,243])/196
and B<-A[,1:47]/reference_column?

Thank you very much for your help!

Best,
  Allen

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

[Metaxab]

Hello,

I am trying to access Microsoft Live Search Using SOAP through R.
In R I am using the RCurl packages to make the calls.

I have the following situation that looks crazy and cannot figure out how to
solve it:

#SOAP Request
library(RCurl)

h = basicTextGatherer()

body='
http://schemas.xmlsoap.org/soap/encoding/";
xmlns:SOAP-ENV="http://schemas.xmlsoap.org/soap/envelope/";
xmlns:xsd=" http://www.w3.org/2001/XMLSchema";
xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance";
xmlns:SOAP-ENC=" http://schemas.xmlsoap.org/soap/encoding/";>
http://tempuri.org";>
MY ID
"Cancer"
en-US
Off


Web
0
10
Title

'

curlPerform(url="http://soap.search.msn.com/webservices.asmx?wsdl";,
httpheader=c(Accept="text/xml", Accept="multipart/*",
 'Content-Type' = "text/xml; charset=utf-8"),
postfields=body,
writefunction = h$update,
verbose = FALSE
)

##Put the responce in body.
body = h$value()
h$reset()

An now the crazy stuff. This code works well for some keywords but not for
all.

So for the Query:"Cancer" the code works fine, but for the Query: "Tiger" I
get the following error message: 
  Error in curlPerform(url =
"http://soap.search.msn.com/webservices.asmx?wsdl";,  : 
  Failed writing body

I cannot figure that out. Please help

Thanks 
Konstantinos Vougas

-- 
View this message in context: 
http://www.nabble.com/Microsoft-SOAP---Help%21%21-tf4792961.html#a13711324
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] Resample a matrix (table)

[Luca Penasa]

I have a matrix of values... for exemple the matrix could be 100 rows
and 50 columns...
I would like to resample the matrix to obtain a new matrix -> a 10x5 matrix
So i need to calculate a new value for every cell, i would like to use a
weightened mean to do this. I explain this thing better:
In the exemple the original matrix would be divided into 50
partial-matrix, for every partial-matrix should be extracted a value
obtained from a weightened mean: the cell at center of the
partial-matrix should have the max importance into the mean, the fhurter
cell should have the lower importance (something like a gaussian shape,
where i can decide the sigma used) In this way i'll obtain a
resampled matrix with the values calculated using a weightened mean of
the single partial-matrix...

Is there anybody who can suggest me a way to do this??? is there a
function that can do this?

PS: sorry for my english... and... im a newbie, please try to explain in
a simple way.. thank you

Luca Penasa,
geology student,
University of Padua, Italy


 
 
 --
 Email.it, the professional e-mail, gratis per te: http://www.email.it/f
 
 Sponsor:
 In REGALO 'All the Good Thing' di NELLY FURTADO
 Clicca qui: http://adv.email.it/cgi-bin/foclick.cgi?mid=6617&d=12-11

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

[billycorg]

I am a very beginner user of R. I have yet sample weights that I would like
to apply to all my data set.
Which function can I apply? The survey package begins from the survey
design.. I have yet, in my dataset, the weights column.

thanks a lot,
Vincenzo
-- 
View this message in context: 
http://www.nabble.com/a-very-little-help-with-sample-weights-tf4792687.html#a13710423
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] Package fEcofin fails to download

[Piotr Chmielowski]

Download of the fEcofin fails from:

http://cran.r-project.org/bin/windows/contrib/r-release/fEcofin_260.72.z
ip

results in IE complaining

Internet Explorer cannot display the webpage 
   
   Most likely causes:
You are not connected to the Internet. 
The website is encountering problems. 
There might be a typing error in the address. 
 
I set up the site to be trusted, and switched off all anti-virus
software.

Installation from RGui fails too. 
However, other packages such as fOptions, fCalendar etc
download/download-install just fine.


Any help would be appreciated

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

[R] help in long loops

[Mahmudul Haque]

hi,

please help me out in the following case. seems like it stuck in some 
where(already 7 hrs passed). what I want is to combine 4 matrix in to one 
matrix of desired length.

final_matrix<-function(ob_feat,content_feat,link_feat,link_feat_transformed){
complete_feat<-matrix(rep(-1),nrow=11402,ncol=278)


for(i in 1:8944)
{q<-c(0)
for(j in 1:11402)
{
if(as.matrix(content_feat[i,2])==as.matrix(ob_feat[j,2]))
{complete_feat[i,1]=as.matrix(ob_feat[j,2])
complete_feat[i,2:97]=as.matrix(content_feat[i,3:98])
complete_feat[i,98:99]=as.matrix(ob_feat[j,3:4])
complete_feat[i,100:140]=as.matrix(link_feat[j,3:43])
complete_feat[i,141:278]=as.matrix(link_feat_transformed[j,3:140])
   q<-1}
if(q==1)
break;
}


}
for (i in 8945:11402){
complete_feat[i,1]=as.matrix(ob_feat[i,2])
complete_feat[i,98:99]=as.matrix(ob_feat[i,3:4])
complete_feat[i,100:140]=as.matrix(link_feat[i,3:43])
complete_feat[i,141:278]=as.matrix(link_feat_transformed[i,3:140])
}

list(complete_feat=complete_feat)
}

kind regards,

mahmudul haque

 __



[[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] producing output as *.spo (spss output format)

[Ivan Uemlianin]

Muenchen, Robert A (Bob) wrote:
> You probably don't want to spend time figuring out the .spo format. From
> SPSS 16 on, that format is obsolete and replaced by the Unicode
> XML-based .spv file. SPSS 16 users need a separate Legacy Viewer to read
> .spo files. -Bob

Great.  I've had a quick Google around but have been unable to find out 
anything about this .spv format.  Is it going to be as closed as .spo? 
Can users of 'legacy' versions of spss get an spv viewer?  I suppose 
these are all questions for spss.

Thanks for all comments.  I'm ready to switch to R - with html or odf 
output where necessary.

Best wishes

Ivan

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

[Roland Kaiser]

Hi all!

I encountered precision problems using pdf().

So far, I found out, that pdf() sets the Mediabox
to even values in dtp points (1/72 inch).

This can be seen in the following example.

page.width <- 13.1/2.54 # in centimeters
page.height <- 19/2.54

pdf(file = "foo.pdf", width = page.width, height = page.height)
plot(0,0)
dev.off()

# The size of the Mediabox in foo.pdf is
# /MediaBox [0 0 371 538]
# but ist should be

dtp.pt <- 1/72

page.width/dtp.pt
# [1] 371.3386
page.height/dtp.pt
# [1] 538.5827

Any ideas to cope with this problem.

Thanks

Roland

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

Re: [R] Largest N Values Efficiently?

[David Katz]

x is a 1XN sparse matrix of numerics. I am using the Matrix package to
represent as a sparse matrix; the representation has a numeric vector
representing the positions within the matrix. My goal is find the columns
with the n largest values, here positive  correlations. Part of my strategy
is to only sort the nonzeros which are available as a numeric vector. 

Thanks for your interest and input.



Prof Brian Ripley wrote:
> 
> What is 'x' here?  What type?  Does it contain NAs?  Are there ties?  R's 
> ordering functions are rather general, and you can gain efficiency by 
> ruling some of these out.
> 
> See ?sort, look at the 'partial' argument, including the comments in the 
> Details.  And also look at ?sort.list.
> 
> sort.int(x) is more efficient than x[order(x)], and x[order(x)[1:n]] is 
> more efficient than x[order(x)][1:n] for most types.
> 
> Finally, does efficiency matter?  As the examples in ?sort show, R can 
> sort a vector of length 2000 is well under 1ms, and 1e7 random normals in 
> less time than they take to generate.  There are not many tasks where 
> gaining efficiency over x[order(x)][1:n] will be important.  E.g.
> 
>> system.time(x <- rnorm(1e6))
> user  system elapsed
> 0.440.000.44
>> system.time(x[order(x)][1:4])
> user  system elapsed
> 1.720.001.72
>> system.time(x2 <- sort.int(x, method = "quick")[1:4])
> user  system elapsed
> 0.310.000.32
>> system.time(min(x))
> user  system elapsed
> 0.020.000.02
>> system.time(x2 <- sort.int(x, partial=1)[1])
> user  system elapsed
> 0.070.000.07
> 
> and do savings of tenths of a second matter?  (There is also 
> quantreg::kselect, if you work out how to use it, which apparently is 
> a bit faster at partial sorting on MacOS X but not elsewhere.)
> 
> 
> On Sun, 11 Nov 2007, David Katz wrote:
> 
>>
>> What is the most efficient alternative to x[order(x)][1:n] where
>> length(x)>>n?
> 
> That is the smallest n values, pace your subject line.
> 
>> I also need the positions of the mins/maxs perhaps by preserving names.
>>
>> Thanks for any suggestions.
>>
> 
> -- 
> 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.
> 
> 

-- 
View this message in context: 
http://www.nabble.com/Largest-N-Values-Efficiently--tf4788033.html#a13708965
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] Rcmdr Mac OS X

[blemcke]

Thanks for your answer.

I downloaded the recent version today (1.3-8)

B

Am 12.11.2007 um 17:11 schrieb John Fox:

> Dear Birgit,
>
> I'm not a Mac user, so can only offer limited help. I've just  
> verified,
> however, that the most recent version of the Rcmdr package on CRAN
> (version 1.3-8) installs and runs fine under R 2.6.0 on my Ubuntu  
> Linux
> system.
>
> What is the version of the Rcmdr package that you installed?
>
> I'm copying this response to Rob Goedman, who may be able to help
> further.
>
> Regards,
>  John
>
> On Mon, 12 Nov 2007 13:31:30 +0100
>  blemcke <[EMAIL PROTECTED]> wrote:
>> Sorry forgot subject.
>>
>> B
>>
>>
>>
>>> Von: Birgit Lemcke <[EMAIL PROTECTED]>
>>> Datum: 12. November 2007 13:26:55 GMT+01:00
>>> An: R Hilfe <[EMAIL PROTECTED]>
>>>
>>> Hello members of the mailinglist!
>>>
>>> I am running R 2.6.0 GUI 1.21 (4815) (4815) on a PPC with Mac OS X
>>
>>> 10.4.
>>>
>>> I just tried to use the R commander.
>>> I followed these instructions:
>>>
>>> 4) Start R (e.g. double-click the R icon in /Applications).
>>> 5) Start X11 (e.g. click the X11 button on the R Console or double-
>>
>>> click the X11 icon in /Applications/Utilities).
>>> 6) Close the X11 terminal (optional).
>>> 7) Type library(Rcmdr) in the R Console and press Return.
>>>
>>> X11 started and also the Rcmdr window opened, but I got this
>>> messsages:
>>>
>>> Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht
>>> gefunden #not found
>>> Error : .onAttach in 'attachNamespace' fehlgeschlagen #failed
>>> Fehler: Laden von Paket/Namensraum für 'Rcmdr' fehlgeschlagen
>>> #Error: Loading of Paket/Namensraum for Rcmdr failed
>>> Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht
>>> gefunden #not found
>>> Error in get(x, envir = RcmdrEnv(), mode = mode, inherits = FALSE)
>> :
>>>   Variable "messagesWindow" nicht gefunden
>>>   #variable "messagesWindow" not found
>>>
>>> When I tried to load the Prestige dataset from the car package I
>>> got this:
>>>
>>> Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht
>>> gefunden
>>>
>>> What am I doing wrong?
>>>
>>> Thanks in advance for your help.
>>>
>>> B
>>>
>>> Birgit Lemcke
>>> Institut für Systematische Botanik
>>> Zollikerstrasse 107
>>> CH-8008 Zürich
>>> Switzerland
>>> Ph: +41 (0)44 634 8351
>>> [EMAIL PROTECTED]
>>>
>>>
>>>
>>>
>>>
>>
>> Birgit Lemcke
>> Altwiesenstrasse 89
>> CH-8051 Zürich
>> Switzerland
>> phone +0041 43 540 16 35
>> [EMAIL PROTECTED]
>>
>>
>> Hier noch ein kleiner Homepagetip: http://www.marlem-software.de/
>> kontakt.html
>>
>> __
>> R-help@r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/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 Fox, Professor
> Department of Sociology
> McMaster University
> Hamilton, Ontario, Canada
> http://socserv.mcmaster.ca/jfox/

Birgit Lemcke
Altwiesenstrasse 89
CH-8051 Zürich
Switzerland
phone +0041 43 540 16 35
[EMAIL PROTECTED]


Hier noch ein kleiner Homepagetip: http://www.marlem-software.de/ 
kontakt.html

__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
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 on drawing a tree with "ape"?

[Hua Li]

Dear Simon and everyone,

Thanks for response. 

Specifically, I have a tree with species "A, B, C" and
can be written as
"((A:45.15,C:45.15):46.19,B:91.34):0.0;" . If I use
the R command


> mytreeABC

[1] "((A:45.15,C:45.15):46.19,B:91.34):0.0;"

> plot(read.tree(text = mytreeABC))

I would get a tree as shown in "Rresult.ps" file.
However, I would like the tree to show like in
"treeABC.jpg" file.

In case the r-hlep list do not allow attachment, here
is what Rresult.ps look like:

|-B
|
|
|  |--C
|--|
   |--A

for the tree I want (in treeABC.jpg), I would mark the
edge AC as "b", and the edge AB as "a". 


Thank you very much for help!

Hua


--- Simon Blomberg <[EMAIL PROTECTED]> wrote:

> ?edgelabels
> 
> On Sun, 2007-11-11 at 16:51 -0800, Hua Li wrote:
> > Dear all,
> > 
> > I'm using the "ape" package in R and want to draw
> a
> > phylogenetic tree with not only the tip labels but
> > also some labels for the edges. e.g. Mark the edge
> AB
> > as "m" in the tree ABC.  
> > 
> > Couldn't find a way to do that. Can someone help?
> > 
> > Thanks,
> > 
> > Hua
> > 
> > __
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained,
> reproducible code.
> -- 
> Simon Blomberg, BSc (Hons), PhD, MAppStat. 
> Lecturer and Consultant Statistician 
> Faculty of Biological and Chemical Sciences 
> The University of Queensland 
> St. Lucia Queensland 4072 
> Australia
> Room 320 Goddard Building (8)
> T: +61 7 3365 2506 
> email: S.Blomberg1_at_uq.edu.au
> 
> Policies:
> 1.  I will NOT analyse your data for you.
> 2.  Your deadline is your problem.
> 
> The combination of some data and an aching desire
> for 
> an answer does not ensure that a reasonable answer
> can 
> be extracted from a given body of data. - John
> Tukey.
> 
> 


__




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

Re: [R] Rcmdr Mac OS X

[John Fox]

Dear Birgit,

I'm not a Mac user, so can only offer limited help. I've just verified,
however, that the most recent version of the Rcmdr package on CRAN
(version 1.3-8) installs and runs fine under R 2.6.0 on my Ubuntu Linux
system.

What is the version of the Rcmdr package that you installed?

I'm copying this response to Rob Goedman, who may be able to help
further.

Regards,
 John

On Mon, 12 Nov 2007 13:31:30 +0100
 blemcke <[EMAIL PROTECTED]> wrote:
> Sorry forgot subject.
> 
> B
> 
> 
> 
> > Von: Birgit Lemcke <[EMAIL PROTECTED]>
> > Datum: 12. November 2007 13:26:55 GMT+01:00
> > An: R Hilfe <[EMAIL PROTECTED]>
> >
> > Hello members of the mailinglist!
> >
> > I am running R 2.6.0 GUI 1.21 (4815) (4815) on a PPC with Mac OS X
>  
> > 10.4.
> >
> > I just tried to use the R commander.
> > I followed these instructions:
> >
> > 4) Start R (e.g. double-click the R icon in /Applications).
> > 5) Start X11 (e.g. click the X11 button on the R Console or double-
> 
> > click the X11 icon in /Applications/Utilities).
> > 6) Close the X11 terminal (optional).
> > 7) Type library(Rcmdr) in the R Console and press Return.
> >
> > X11 started and also the Rcmdr window opened, but I got this  
> > messsages:
> >
> > Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht  
> > gefunden #not found
> > Error : .onAttach in 'attachNamespace' fehlgeschlagen #failed
> > Fehler: Laden von Paket/Namensraum für 'Rcmdr' fehlgeschlagen
> > #Error: Loading of Paket/Namensraum for Rcmdr failed
> > Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht  
> > gefunden #not found
> > Error in get(x, envir = RcmdrEnv(), mode = mode, inherits = FALSE)
> :
> >   Variable "messagesWindow" nicht gefunden
> >   #variable "messagesWindow" not found
> >
> > When I tried to load the Prestige dataset from the car package I  
> > got this:
> >
> > Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht  
> > gefunden
> >
> > What am I doing wrong?
> >
> > Thanks in advance for your help.
> >
> > B
> >
> > Birgit Lemcke
> > Institut für Systematische Botanik
> > Zollikerstrasse 107
> > CH-8008 Zürich
> > Switzerland
> > Ph: +41 (0)44 634 8351
> > [EMAIL PROTECTED]
> >
> >
> >
> >
> >
> 
> Birgit Lemcke
> Altwiesenstrasse 89
> CH-8051 Zürich
> Switzerland
> phone +0041 43 540 16 35
> [EMAIL PROTECTED]
> 
> 
> Hier noch ein kleiner Homepagetip: http://www.marlem-software.de/ 
> kontakt.html
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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 Fox, Professor
Department of Sociology
McMaster University
Hamilton, Ontario, Canada
http://socserv.mcmaster.ca/jfox/

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

Re: [R] strange `nls' behaviour

[Duncan Murdoch]

On 11/12/2007 9:14 AM, Joerg van den Hoff wrote:
> On Mon, Nov 12, 2007 at 07:36:34AM -0500, Duncan Murdoch wrote:
>> On 11/12/2007 6:51 AM, Joerg van den Hoff wrote:
>> >I initially thought, this should better be posted to r-devel
>> >but alas! no response. 
>> 
>> I think the reason there was no response is that your example is too 
>> complicated.  You're doing a lot of strange things (fitfunc as a result 
>> of deriv, using as.name, as.call, as.formula, etc.)  You should simplify 
> 
> thanks for the feedback.
> 
> concerning  "lot  of  strange  things":  OK.  I  thought the
> context might be important ("why, for heaven's sake  do  you
> want  to  do  this!?"), but, then, maybe not. so the easiest
> way to trigger a similar (not the  identical)  behaviour  is
> something like
> 
> f <- function (n) {
>#-
>#define n data points for a (hardcoded) model:
>#---
>x <- seq(0, 5, length  = n)
>y <- 2 * exp(-1*x) + 2; 
>y <- rnorm(y,y, 0.01*y)
> 
>#the model (same as the above hardcoded one):
>model <- y ~ a * exp (-b*x) + c
> 
>#call nls as usual:
>res1 <- try(nls(model, start = c(a=2, b=1, c=2)))
> 
>#call it a bit differently:
>res2 <- nls(y ~ eval(model[[3]]), start = c(a=2, b=1, c=2))
> 
>list(res1 = res1, res2 = res2)
>#-
> }

I'd say the problem is relying on the default for the envir parameter to 
eval.  It's reasonable to expect nls to set things up so that terms in 
the model formula are taken from the right place, but when your model 
formula is playing with evaluation, you should be very careful to make 
sure the evaluation takes place in the right context.

The default for envir is "parent.frame()", and that can't be right: that 
will see local variables in whatever function called it, so if one of 
them has a variable called "model", you'll subscript it.

If I were trying to do what you're doing, I would construct the formula 
before the call to nls, and pass that.  I.e. something like the 
following silly code:

model2 <- model
model2[[3]] <- model[[3]] # The eval() is implicit
res2 <- nls(model2, start = c(a=2, b=1, c=2))

If you really want to put eval() in a formula, then I can't see how you 
could avoid an explicit specification of the envir parameter.  So I'd 
call this a bug in your code.

Duncan Murdoch


> this is without all the overhead of my first example and now
> (since not quite the same) the problem  arises  at  n  ==  3
> where  the  fit  can't  really  procede  (there  are  also 3
> parameters -- the first example was more  relevant  in  this
> respect)  but  anyway  the  second nls-call does not procede
> beyond the parsing phase of `model.frame'.
> 
> so,  I  can't  see  room for a real bug in my code, but very
> well a chance that I misuse `nls'  (i.e.  not  understanding
> what really is tolerable for the `model' argument of `nls').
> but  if the latter is not the case, I would think there is a
> bug in `nls'  (or,  actually,  rather  in  `model.frame'  or
> whatever)  when  parsing  the  nls call.


> 
> 
>> it down to isolate the bug.  Thats a lot of work, but you're the one in 
>> the best position to do it.  I'd say there's at least an even chance 
>> that the bug is in your code rather than in nls().
>> 
>> And 2.5.0 *is* ancient; please confirm the bug exists in R-patched if it 
>> turns out to be an R bug.
> 
> if  need  be,  I'll  do  that  (if  I  get it compiled under
> macosX). but assuming  that  you  have  R-patched  installed
> anyway, I would appreciate if you would copy the new example
> above or the old one  below  to  your  R-  prompt  and  see,
> whether  it  crashes  with  the same error message if called
> with the special number of data points (3 for the new, 5 for
> the  old  example)?  and  if  so,  maybe  bring  this to the
> attention of d. bates?
> 
> 
> j. van den hoff
>> 
>> Duncan Murdoch
>> 
> 
>> 
>> 
>> 
>> so  I  try  it  here.  sory  for  the
>> >lengthy explanation but it seems unavoidable. to quickly see
>> >the problem simply copy the litte example below and execute
>> >
>> >f(n=5)
>> >
>> >which  crashes. called with n !=  5 (and of course n>3 since
>> >there are 3 parameters in the model...) everything is as  it
>> >should be.
>> >
>> >in detail:
>> >I  stumbled over the follwing _very_ strange behaviour/error
>> >when using `nls' which  I'm  tempted  (despite  the  implied
>> >"dangers") to call a bug:
>> >
>> >I've  written a driver for `nls' which allows specifying the
>> >model and the data vectors using arbitrary  symbols.   these
>> >are  internally  mapped  to  consistent names, which poses a
>> >slight complication when using `deriv' to  provide  analytic
>> >derivatives. the following fragment gives the idea:
>> >
>> >#-
>> >f <- function(n = 4) {
>> >
>> >   x <- seq(0, 5, length  = n)
>> >
>> >   y <- 2 * exp(-1

Re: [R] plot

[Benilton Carvalho]

replace the second plot() by

lines(d7, (a*d7)+b)

b

On Nov 12, 2007, at 10:42 AM, elyakhlifi mustapha wrote:

> Hi,
> I have been searched how to do to display in the same chart this  
> following functions:
>
> plot(d7,d6,type="b",col="dark red")
> plot(d7,(a*d7)+b,type="l")
>
> in using par it display in 2 diffenrent charts.
>
>
>   
> _
>
> l
>   [[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

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

[R] plot

[elyakhlifi mustapha]

Hi,
I have been searched how to do to display in the same chart this following 
functions:

plot(d7,d6,type="b",col="dark red")
plot(d7,(a*d7)+b,type="l")

in using par it display in 2 diffenrent charts.


  
_ 

l 
[[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] strange `nls' behaviour

[Katharine Mullen]

On about line 547 in nls.R there is
  mf$formula <-  # replace by one-sided linear model formula
  as.formula(paste("~", paste(varNames[varIndex], collapse = "+")),
   env = environment(formula))

If this is replaced with

mf$formula <-  # replace by one-sided linear model formula
  as.formula(paste("~", paste(varNames[1], collapse = "+")),
   env = environment(formula))

then the behviour seems to be OK at least for Joerg's examples.

It remains to be determined whether this is a _good_ solution (that is, if
it's a general solution to get the desired behavior).

On Mon, 12 Nov 2007, Martin Maechler wrote:

> > "DM" == Duncan Murdoch <[EMAIL PROTECTED]>
> > on Mon, 12 Nov 2007 07:36:34 -0500 writes:
>
> DM> On 11/12/2007 6:51 AM, Joerg van den Hoff wrote:
> >> I initially thought, this should better be posted to r-devel
> >> but alas! no response.
>
> DM> I think the reason there was no response is that your example is too
> DM> complicated.  You're doing a lot of strange things (fitfunc as a 
> result
> DM> of deriv, using as.name, as.call, as.formula, etc.)  You should 
> simplify
> DM> it down to isolate the bug.  Thats a lot of work, but you're the one 
> in
> DM> the best position to do it.  I'd say there's at least an even chance
> DM> that the bug is in your code rather than in nls().
>
> yes.. and.. no :
> - His code is quite peculiar, but I think only slightly too complicated
>
> - one could argue that the bug is in Joerg's thinking that
>   something like
>   nls(y ~ eval(fitfunc), )
>
>   should be working at all.
>   But then he had found by experiment that it (accidentally I   d'say)
>   does work in many cases.
>
> DM> And 2.5.0 *is* ancient; please confirm the bug exists in R-patched if 
> it
> DM> turns out to be an R bug.
>
> You are right, but indeed (as has Kate just said) it *does*
> exist in current R versions.
>
> I agree that the behavior of nls() here is sub-optimal.
> It *should* be consistent, i.e. work the same for n=4,5,6,..
>
> I had spent about an hour after Joerg's R-devel posting,
> and found to be too busy with more urgent matters --
> unfortunately forgetting to give *some* feedback about my findings.
>
> It may well be that we find that nls() should give an
> (intelligible) error message in such eval() cases - rather than
> only in one case...
>
> Martin Maechler
>
>
> DM> Duncan Murdoch
>
>
>
>
> DM> so  I  try  it  here.  sory  for  the
> >> lengthy explanation but it seems unavoidable. to quickly see
> >> the problem simply copy the litte example below and execute
> >>
> >> f(n=5)
> >>
> >> which  crashes. called with n !=  5 (and of course n>3 since
> >> there are 3 parameters in the model...) everything is as  it
> >> should be.
> >>
> >> in detail:
> >> I  stumbled over the follwing _very_ strange behaviour/error
> >> when using `nls' which  I'm  tempted  (despite  the  implied
> >> "dangers") to call a bug:
> >>
> >> I've  written a driver for `nls' which allows specifying the
> >> model and the data vectors using arbitrary  symbols.   these
> >> are  internally  mapped  to  consistent names, which poses a
> >> slight complication when using `deriv' to  provide  analytic
> >> derivatives. the following fragment gives the idea:
> >>
> >> #-
> >> f <- function(n = 4) {
> >>
> >> x <- seq(0, 5, length  = n)
> >>
> >> y <- 2 * exp(-1*x) + 2;
> >> y <- rnorm(y,y, 0.01*y)
> >>
> >> model <- y ~ a * exp (-b*x) + c
> >>
> >> fitfunc <- deriv(model[[3]], c("a", "b", "c"), c("a", "b", "c", "x"))
> >>
> >> #"standard" call of nls:
> >> res1 <- nls(y ~ fitfunc(a, b, c, x), start = c(a=1, b=1, c=1))
> >>
> >> call.fitfunc <-
> >> c(list(fitfunc), as.name("a"), as.name("b"), as.name("c"), 
> as.name("x"))
> >> call.fitfunc <- as.call(call.fitfunc)
> >> frml <- as.formula("y ~ eval(call.fitfunc)")
> >>
> >> #"computed" call of nls:
> >> res2 <- nls(frml, start = c(a=1, b=1, c=1))
> >>
> >> list(res1 = res1, res2 = res2)
> >> }
> >> #-
> >>
> >> the  argument  `n'   defines  the number of (simulated) data
> >> points x/y which are going to be fitted by some model ( here
> >> y ~ a*exp(-b*x)+c )
> >>
> >> the first call to `nls' is the standard way of calling `nls'
> >> when knowing all the variable and parameter names.
> >>
> >> the second call (yielding `res2') uses a constructed formula
> >> in `frml' (which in this example is of course not necessary,
> >> but  in  the general case 'a,b,c,x,y' are not a priori known
> >> names).
> >>
> >> now, here is the problem: the call
> >>
> >> f(4)
> >>
> >> runs fine/consistently, as doe

Re: [R] bootstrap

[Prof Brian Ripley]

On Mon, 12 Nov 2007, Stefano Ghirlanda wrote:

> i am using the boot package for some bootstrap calculations in place
> of anovas. one reason is that my dependent variable is distributed
> bimodally, but i would also like to learn about bootstrapping in
> general (i have ordered books but they have not yet arrived).
>
> i get the general idea of bootstrapping but sometimes i do not know
> how to define suitable statistics to test specific hypotheses.

That's a basic issue in statistics.  Bootstrapping is only another way to 
assess the variability of a pre-determined statistic (and incidentally is 
not much used for testing, more often for confidence intervals).

> two examples follow.
>
> 1) comparing the means of more than two groups. a suitable statistics
>   could be the sum of squared deviations of group means from the
>   grand mean. does this sound reasonable?

No.  That means nothing by itself, but needs to be compared to the 
residual variation (e.g. by an F statistic).

> 2) testing for interactions. e.g., i want to see whether an
>   independent variable has the same effect in two different
>   samples. in an anova this would be expressed as the significance,
>   or lack thereof, of the interaction between a "sample" factor and
>   another factor for the independent variable. how would i do this
>   with a bootstrap calculation?
>
> my problem with 2) is that when one fits a linear model to the data,
> from which sums of squares for the anova are calculated, the
> interaction between the two factors corresponds to many regression
> coefficients in the linear model (e.g., i actually have three samples
> and an independent variable with four levels). i do not know how to
> summarize these in a single statistics.

Any good book on statistics with R (e.g. MASS) would point you at the 
anova() function to compare models.

> i have seen somewhere that some people calculate F ratios
> nevertheless, but then test them against a bootstrapped distribution
> rather than against the F distribution. is this a sensible approach?
> could one also use sums of squares directly as the bootstrapped
> statistics?

It is not totally off the wall.  The problems are

- How you bootstrap, given that you don't have a single homogeneous group.
   You seem to want to test, so you need to emulate the null-hypothesis
   distribution.  The most common way to do that is to fit a model, find
   some sort of residuals, bootstrap those and use them as errors to
   simulate from the null hypothesis.  At that point you will have to work
   hard to convince many statisticians that you have improved over the
   standard theory or a simulation from a long-tailed error distribution.

- That if you don't believe in the normal distribution of your errors
   (and not response), you probably should not be using least-squares
   based statistical methodology.  And remember that the classical ANOVA
   tests are supported by permutation arguments which are very similar to
   the bootstrap (just sampling without replacement instead of with).

These points are discussed with examples in the linear models chapter of 
MASS (the book) and also in the Davison-Hinkley book which the 'boot' 
package supports.

[Shame about the broken shift key, although it seems to work with F: 
keyboards are really cheap to replace these days.]

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

[R] R 2.6.0 & RMySQL

[Joachim Claudet]

Dear Prof. Ripley (and other readers) --

Unfortunately, I already did this (using the packages menu to install 
RMySQL and copying libMySQL.dll in the RMySQL/libs directory). But after 
the installation, when I load the package, I have the same error message.
Regards,

- Joachim Claudet

-- 
<º)))><  <º)))><  <º)))><  <º)))><  <º)))><  <º)))><  <º)))><  <º)))><
Joachim Claudet - PhD
http://joachim.claudet.free.fr/
Skype: joachim.claudet
Laboratory of Zoology and Marine Biology
Department of Biological and Environmental Science and Technologies
University of Salento - via Monteroni 73100 Lecce, Italy
Phone: +39 0832 298853
Fax: +39 0832 298626
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] strange `nls' behaviour

[Joerg van den Hoff]

On Mon, Nov 12, 2007 at 03:25:38PM +0100, Martin Maechler wrote:
> > "DM" == Duncan Murdoch <[EMAIL PROTECTED]>
> > on Mon, 12 Nov 2007 07:36:34 -0500 writes:
> 
> DM> On 11/12/2007 6:51 AM, Joerg van den Hoff wrote:
> >> I initially thought, this should better be posted to r-devel
> >> but alas! no response. 
> 
> DM> I think the reason there was no response is that your example is too 
> DM> complicated.  You're doing a lot of strange things (fitfunc as a 
> result 
> DM> of deriv, using as.name, as.call, as.formula, etc.)  You should 
> simplify 
> DM> it down to isolate the bug.  Thats a lot of work, but you're the one 
> in 
> DM> the best position to do it.  I'd say there's at least an even chance 
> DM> that the bug is in your code rather than in nls().
> 
> yes.. and.. no : 
> - His code is quite peculiar, but I think only slightly too complicated

thank's  for  this  response. I just tried to give an easier
example in response to duncan's mail.

> 
> - one could argue that the bug is in Joerg's thinking that
>   something like
>   nls(y ~ eval(fitfunc), )
> 
>   should be working at all.
>   But then he had found by experiment that it (accidentally I   d'say)
>   does work in many cases.
> 
> DM> And 2.5.0 *is* ancient; please confirm the bug exists in R-patched if 
> it 
> DM> turns out to be an R bug.
> 
> You are right, but indeed (as has Kate just said) it *does*
> exist in current R versions.
> 
> I agree that the behavior of nls() here is sub-optimal.
> It *should* be consistent, i.e. work the same for n=4,5,6,..
> 
> I had spent about an hour after Joerg's R-devel posting,
> and found to be too busy with more urgent matters --
> unfortunately forgetting to give *some* feedback about my findings.
> 
> It may well be that we find that nls() should give an
> (intelligible) error message in such eval() cases - rather than
> only in one case...

well,  if at all possible, it would _really_ be nice to have
the opportunity  to  pass  `eval'  constructs  throught  the
`model'  argument  to `nls'. I'm far from understanding what
actually is going on in `model.frame', but if it parses  the
thing  correctly  ("accidentally"  or  not)  except  in  one
special case, is'nt  there  a  clean  way  to  "repair"  the
special  case?  in my application I found the ability to use
`eval' constructs here really helpful. at least  I  did  not
find another way to be able to use `deriv' information _and_
allow the user to use arbitrary symbols  in  specifying  the
model  formula. at least I'd prefer the present situation of
"works nearly always" to a consistent "works never".

joerg

> 
> Martin Maechler
> 
> 
> DM> Duncan Murdoch
> 
> 
> 
> 
> DM> so  I  try  it  here.  sory  for  the
> >> lengthy explanation but it seems unavoidable. to quickly see
> >> the problem simply copy the litte example below and execute
> >> 
> >> f(n=5)
> >> 
> >> which  crashes. called with n !=  5 (and of course n>3 since
> >> there are 3 parameters in the model...) everything is as  it
> >> should be.
> >> 
> >> in detail:
> >> I  stumbled over the follwing _very_ strange behaviour/error
> >> when using `nls' which  I'm  tempted  (despite  the  implied
> >> "dangers") to call a bug:
> >> 
> >> I've  written a driver for `nls' which allows specifying the
> >> model and the data vectors using arbitrary  symbols.   these
> >> are  internally  mapped  to  consistent names, which poses a
> >> slight complication when using `deriv' to  provide  analytic
> >> derivatives. the following fragment gives the idea:
> >> 
> >> #-
> >> f <- function(n = 4) {
> >> 
> >> x <- seq(0, 5, length  = n)
> >> 
> >> y <- 2 * exp(-1*x) + 2; 
> >> y <- rnorm(y,y, 0.01*y)
> >> 
> >> model <- y ~ a * exp (-b*x) + c
> >> 
> >> fitfunc <- deriv(model[[3]], c("a", "b", "c"), c("a", "b", "c", "x"))
> >> 
> >> #"standard" call of nls:
> >> res1 <- nls(y ~ fitfunc(a, b, c, x), start = c(a=1, b=1, c=1))
> >> 
> >> call.fitfunc <- 
> >> c(list(fitfunc), as.name("a"), as.name("b"), as.name("c"), 
> as.name("x"))
> >> call.fitfunc <- as.call(call.fitfunc)
> >> frml <- as.formula("y ~ eval(call.fitfunc)")
> >> 
> >> #"computed" call of nls:
> >> res2 <- nls(frml, start = c(a=1, b=1, c=1))
> >> 
> >> list(res1 = res1, res2 = res2)
> >> }
> >> #-
> >> 
> >> the  argument  `n'   defines  the number of (simulated) data
> >> points x/y which are going to be fitted by some model ( here
> >> y ~ a*exp(-b*x)+c )
> >> 
> >> the first call to `nls' is the standard way of calling `nls'
> >> when knowing all the variable and parameter names.
> >> 
> >> the second call (yielding `res2') uses a constructed formula
> 

Re: [R] R 2.6.0 & RMySQL

[Prof Brian Ripley]

I don't understand why you are trying to install a 2.5.x binary on R 
2.6.0.  Perhaps if you tried not to make things so hard for yourself and 
just used the Packages menu to install RMySQL this might work: it does for 
me and quite a few others.


Otherwise, you do need libMySQL.dll on the path or copied to the same 
directory as RMySQL.dll to be sure that is not the issue.  I suspect it is 
rather an R version mismatch.


On Mon, 12 Nov 2007, Joachim Claudet wrote:


Dear list members --

I have problems in loading the pre-compiled library RMySQL_0.5-7 taken
here (http://stat.bell-labs.com/RS-DBI/download/index.html) or
RMySQL_0.6-0 taken here
(http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.5/) on R 2.6.0
under Windows XP. I am using the MySQL features of EasyPHP1-8 (but the
problem is not (for the moment) in connecting with a database but just
in loading the package)



I have seen all the previous posts explaining to copy the libMySQL.dll &
RMySQL.dll of the RMySQL library in the Windows/System32 folder or to
change the path variable in the Environment Variables of the system but
nothing helps.


I do wish people would not write things like that.  Give URLs of what you 
have read, since


- very often it is not 'all'.
- sometimes it is misread.
- some of the self-appointed experts are founts of misinformation (as
  apparently here), so we need to know the credentials.


I still have the following message (I am sorry but, unfortunately, it is
in French; with my translation):


You can always start R in English (see the rw-FAQ) to quote English 
messages.  I would have expected there to be a message box as well (at 
least under XP) saying what the missing procedure was (and that would be 
in French, no doubt).




Error in dyn.load(file, ...) :
impossible de charger la bibliothèque partagée (="impossible to charge
the shared library")
'C:/PROGRA~1/R/R-26~1.0/library/RMySQL/libs/RMySQL.dll':
LoadLibrary failure:  La procédure spécifiée est introuvable (="the
specified procedure cannot be found").
Erreur : le chargement du package / espace de noms a échoué pour
'RMySQL' (="loading package/space of names failed for")


I have just changed my computer and was previously on R 2.4, where I did
not have any problems in installing and loading pre-compiled versions
RMySQL.

Any help will be welcome.

- Joachim Claudet




--
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] strange `nls' behaviour

[Martin Maechler]

> "DM" == Duncan Murdoch <[EMAIL PROTECTED]>
> on Mon, 12 Nov 2007 07:36:34 -0500 writes:

DM> On 11/12/2007 6:51 AM, Joerg van den Hoff wrote:
>> I initially thought, this should better be posted to r-devel
>> but alas! no response. 

DM> I think the reason there was no response is that your example is too 
DM> complicated.  You're doing a lot of strange things (fitfunc as a result 
DM> of deriv, using as.name, as.call, as.formula, etc.)  You should 
simplify 
DM> it down to isolate the bug.  Thats a lot of work, but you're the one in 
DM> the best position to do it.  I'd say there's at least an even chance 
DM> that the bug is in your code rather than in nls().

yes.. and.. no : 
- His code is quite peculiar, but I think only slightly too complicated

- one could argue that the bug is in Joerg's thinking that
  something like
nls(y ~ eval(fitfunc), )

  should be working at all.
  But then he had found by experiment that it (accidentally I   d'say)
  does work in many cases.

DM> And 2.5.0 *is* ancient; please confirm the bug exists in R-patched if 
it 
DM> turns out to be an R bug.

You are right, but indeed (as has Kate just said) it *does*
exist in current R versions.

I agree that the behavior of nls() here is sub-optimal.
It *should* be consistent, i.e. work the same for n=4,5,6,..

I had spent about an hour after Joerg's R-devel posting,
and found to be too busy with more urgent matters --
unfortunately forgetting to give *some* feedback about my findings.

It may well be that we find that nls() should give an
(intelligible) error message in such eval() cases - rather than
only in one case...

Martin Maechler


DM> Duncan Murdoch




DM> so  I  try  it  here.  sory  for  the
>> lengthy explanation but it seems unavoidable. to quickly see
>> the problem simply copy the litte example below and execute
>> 
>> f(n=5)
>> 
>> which  crashes. called with n !=  5 (and of course n>3 since
>> there are 3 parameters in the model...) everything is as  it
>> should be.
>> 
>> in detail:
>> I  stumbled over the follwing _very_ strange behaviour/error
>> when using `nls' which  I'm  tempted  (despite  the  implied
>> "dangers") to call a bug:
>> 
>> I've  written a driver for `nls' which allows specifying the
>> model and the data vectors using arbitrary  symbols.   these
>> are  internally  mapped  to  consistent names, which poses a
>> slight complication when using `deriv' to  provide  analytic
>> derivatives. the following fragment gives the idea:
>> 
>> #-
>> f <- function(n = 4) {
>> 
>> x <- seq(0, 5, length  = n)
>> 
>> y <- 2 * exp(-1*x) + 2; 
>> y <- rnorm(y,y, 0.01*y)
>> 
>> model <- y ~ a * exp (-b*x) + c
>> 
>> fitfunc <- deriv(model[[3]], c("a", "b", "c"), c("a", "b", "c", "x"))
>> 
>> #"standard" call of nls:
>> res1 <- nls(y ~ fitfunc(a, b, c, x), start = c(a=1, b=1, c=1))
>> 
>> call.fitfunc <- 
>> c(list(fitfunc), as.name("a"), as.name("b"), as.name("c"), as.name("x"))
>> call.fitfunc <- as.call(call.fitfunc)
>> frml <- as.formula("y ~ eval(call.fitfunc)")
>> 
>> #"computed" call of nls:
>> res2 <- nls(frml, start = c(a=1, b=1, c=1))
>> 
>> list(res1 = res1, res2 = res2)
>> }
>> #-
>> 
>> the  argument  `n'   defines  the number of (simulated) data
>> points x/y which are going to be fitted by some model ( here
>> y ~ a*exp(-b*x)+c )
>> 
>> the first call to `nls' is the standard way of calling `nls'
>> when knowing all the variable and parameter names.
>> 
>> the second call (yielding `res2') uses a constructed formula
>> in `frml' (which in this example is of course not necessary,
>> but  in  the general case 'a,b,c,x,y' are not a priori known
>> names).
>> 
>> now, here is the problem: the call
>> 
>> f(4)
>> 
>> runs fine/consistently, as does every call with n > 5.
>> 
>> BUT: for n = 5 (i.e. issuing f(5))
>> 
>> the second fit leads to the error message:
>> 
>> "Error in model.frame(formula, rownames, variables, varnames, extras, 
extranames,  : 
>> invalid type (language) for variable 'call.fitfunc'"
>> 
>> I  cornered  this  to a spot in `nls' where a model frame is
>> constructed in variable `mf'.  the parsing/constructing here
>> seems  simply  to  be messed up for n = 5: `call.fitfunc' is
>> interpreted as variable.
>> 
>> I,  moreover, empirically noted that the problem occurs when
>> the total number of  parameters  plus  dependent/independent
>> variables  equals  the number of data points (in the present
>> example a,b,c,x,y).
>> 
>> so it is not the 'magic' number of 5 but rather the identity

Re: [R] strange `nls' behaviour

[Joerg van den Hoff]

On Mon, Nov 12, 2007 at 07:36:34AM -0500, Duncan Murdoch wrote:
> On 11/12/2007 6:51 AM, Joerg van den Hoff wrote:
> >I initially thought, this should better be posted to r-devel
> >but alas! no response. 
> 
> I think the reason there was no response is that your example is too 
> complicated.  You're doing a lot of strange things (fitfunc as a result 
> of deriv, using as.name, as.call, as.formula, etc.)  You should simplify 

thanks for the feedback.

concerning  "lot  of  strange  things":  OK.  I  thought the
context might be important ("why, for heaven's sake  do  you
want  to  do  this!?"), but, then, maybe not. so the easiest
way to trigger a similar (not the  identical)  behaviour  is
something like

f <- function (n) {
   #-
   #define n data points for a (hardcoded) model:
   #---
   x <- seq(0, 5, length  = n)
   y <- 2 * exp(-1*x) + 2; 
   y <- rnorm(y,y, 0.01*y)

   #the model (same as the above hardcoded one):
   model <- y ~ a * exp (-b*x) + c

   #call nls as usual:
   res1 <- try(nls(model, start = c(a=2, b=1, c=2)))

   #call it a bit differently:
   res2 <- nls(y ~ eval(model[[3]]), start = c(a=2, b=1, c=2))

   list(res1 = res1, res2 = res2)
   #-
}

this is without all the overhead of my first example and now
(since not quite the same) the problem  arises  at  n  ==  3
where  the  fit  can't  really  procede  (there  are  also 3
parameters -- the first example was more  relevant  in  this
respect)  but  anyway  the  second nls-call does not procede
beyond the parsing phase of `model.frame'.

so,  I  can't  see  room for a real bug in my code, but very
well a chance that I misuse `nls'  (i.e.  not  understanding
what really is tolerable for the `model' argument of `nls').
but  if the latter is not the case, I would think there is a
bug in `nls'  (or,  actually,  rather  in  `model.frame'  or
whatever)  when  parsing  the  nls call.


> it down to isolate the bug.  Thats a lot of work, but you're the one in 
> the best position to do it.  I'd say there's at least an even chance 
> that the bug is in your code rather than in nls().
> 
> And 2.5.0 *is* ancient; please confirm the bug exists in R-patched if it 
> turns out to be an R bug.

if  need  be,  I'll  do  that  (if  I  get it compiled under
macosX). but assuming  that  you  have  R-patched  installed
anyway, I would appreciate if you would copy the new example
above or the old one  below  to  your  R-  prompt  and  see,
whether  it  crashes  with  the same error message if called
with the special number of data points (3 for the new, 5 for
the  old  example)?  and  if  so,  maybe  bring  this to the
attention of d. bates?


j. van den hoff
> 
> Duncan Murdoch
> 

> 
> 
> 
> so  I  try  it  here.  sory  for  the
> >lengthy explanation but it seems unavoidable. to quickly see
> >the problem simply copy the litte example below and execute
> >
> >f(n=5)
> >
> >which  crashes. called with n !=  5 (and of course n>3 since
> >there are 3 parameters in the model...) everything is as  it
> >should be.
> >
> >in detail:
> >I  stumbled over the follwing _very_ strange behaviour/error
> >when using `nls' which  I'm  tempted  (despite  the  implied
> >"dangers") to call a bug:
> >
> >I've  written a driver for `nls' which allows specifying the
> >model and the data vectors using arbitrary  symbols.   these
> >are  internally  mapped  to  consistent names, which poses a
> >slight complication when using `deriv' to  provide  analytic
> >derivatives. the following fragment gives the idea:
> >
> >#-
> >f <- function(n = 4) {
> >
> >   x <- seq(0, 5, length  = n)
> >
> >   y <- 2 * exp(-1*x) + 2; 
> >   y <- rnorm(y,y, 0.01*y)
> >
> >   model <- y ~ a * exp (-b*x) + c
> >
> >   fitfunc <- deriv(model[[3]], c("a", "b", "c"), c("a", "b", "c", "x"))
> >
> >   #"standard" call of nls:
> >   res1 <- nls(y ~ fitfunc(a, b, c, x), start = c(a=1, b=1, c=1))
> >
> >   call.fitfunc <- 
> >   c(list(fitfunc), as.name("a"), as.name("b"), as.name("c"), as.name("x"))
> >   call.fitfunc <- as.call(call.fitfunc)
> >   frml <- as.formula("y ~ eval(call.fitfunc)")
> >
> >   #"computed" call of nls:
> >   res2 <- nls(frml, start = c(a=1, b=1, c=1))
> >
> >   list(res1 = res1, res2 = res2)
> >}
> >#-
> >
> >the  argument  `n'   defines  the number of (simulated) data
> >points x/y which are going to be fitted by some model ( here
> >y ~ a*exp(-b*x)+c )
> >
> >the first call to `nls' is the standard way of calling `nls'
> >when knowing all the variable and parameter names.
> >
> >the second call (yielding `res2') uses a constructed formula
> >in `frml' (which in this example is of course not necessary,
> >but  in  the general case 'a,b,c,x,y' are not a priori known
> >names).
> >
> >now, here is the problem: the call
> >
> >f(4)
> >
> >runs fine/consistently, as does every 

[R] R 2.6.0 & RMySQL

[Joachim Claudet]

Dear list members --

I have problems in loading the pre-compiled library RMySQL_0.5-7 taken 
here (http://stat.bell-labs.com/RS-DBI/download/index.html) or 
RMySQL_0.6-0 taken here 
(http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/2.5/) on R 2.6.0 
under Windows XP. I am using the MySQL features of EasyPHP1-8 (but the 
problem is not (for the moment) in connecting with a database but just 
in loading the package)
I have seen all the previous posts explaining to copy the libMySQL.dll & 
RMySQL.dll of the RMySQL library in the Windows/System32 folder or to 
change the path variable in the Environment Variables of the system but 
nothing helps.
I still have the following message (I am sorry but, unfortunately, it is 
in French; with my translation):

Error in dyn.load(file, ...) :
impossible de charger la bibliothèque partagée (="impossible to charge 
the shared library") 
'C:/PROGRA~1/R/R-26~1.0/library/RMySQL/libs/RMySQL.dll':
LoadLibrary failure:  La procédure spécifiée est introuvable (="the 
specified procedure cannot be found").
Erreur : le chargement du package / espace de noms a échoué pour 
'RMySQL' (="loading package/space of names failed for")


I have just changed my computer and was previously on R 2.4, where I did 
not have any problems in installing and loading pre-compiled versions 
RMySQL.

Any help will be welcome.

- Joachim Claudet

-- 
<º)))><  <º)))><  <º)))><  <º)))><  <º)))><  <º)))><  <º)))><  <º)))><
Joachim Claudet - PhD
http://joachim.claudet.free.fr/
Skype: joachim.claudet
Laboratory of Zoology and Marine Biology
Department of Biological and Environmental Science and Technologies
University of Salento - via Monteroni 73100 Lecce, Italy
Phone: +39 0832 298853
Fax: +39 0832 298626
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.

[R] Resid() and estimable() functions with lmer

[Adriana Puentes]

Hi all,

Two questions:

1. Is there a way to evaluate models from lmer() with a poisson  
distribution? I get the following error message:

library(lme4)
lmer(tot.fruit~infl.treat+def.treat+(1|initial.size),family=poisson)->model
plot(fitted(model),resid(model))
Error: 'resid' is not implemented yet

Are there any other options?

2. Why doesn't the function estimable() in gmodels work with lmer()  
using a poisson distribution? I know that my coefficients are right  
because it works fine if I run the model without poisson distribution.

#CODE (with latest R version and package updates)#
#Without Poisson distribution#

library(lme4)
library(gmodels)
lmer(tot.fruit~infl.treat+def.treat+(1|initial.size))->model1
summary(model1)

Linear mixed-effects model fit by REML
Formula: tot.fruit ~ infl.treat + def.treat + (1 | initial.size)
   AIC  BIC logLik MLdeviance REMLdeviance
  2449 2471  -1218   2447 2437
Random effects:
  Groups   NameVariance Std.Dev.
  initial.size (Intercept) 71.585   8.4608
  Residual 49.856   7.0608
number of obs: 323, groups: initial.size, 292

Fixed effects:
 Estimate Std. Error t value
(Intercept)  12.8846 1.3028   9.890
infl.treat1  -0.4738 1.1819  -0.401
def.treat2   -3.5522 1.6022  -2.217
def.treat3   -2.1757 1.6461  -1.322
def.treat4   -2.1613 1.7003  -1.271

Correlation of Fixed Effects:
 (Intr) infl.1 df.tr2 df.tr3
infl.treat1 -0.413
def.treat2  -0.616  0.013
def.treat3  -0.641  0.002  0.493
def.treat4  -0.638  0.028  0.469  0.524

#Coefficients#

mean.infl.treat0=c(1,0,1/4,1/4,1/4)
mean.infl.treat1=c(1,1,1/4,1/4,1/4)
mean.def.treat1=c(1,1/2,0,0,0)
mean.def.treat2=c(1,1/2,1,0,0)
mean.def.treat3=c(1,1/2,0,1,0)
mean.def.treat4=c(1,1/2,0,0,1)
means=rbind(mean.infl.treat0,mean.infl.treat1,mean.def.treat1,mean.def.treat2,mean.def.treat3,mean.def.treat4)
estimable(model1,means,sim.lmer=TRUE,n.sim=1000)

   Estimate Std. Error p value
(1 0 0.25 0.25 0.25) 10.897825  0.8356260   0
(1 1 0.25 0.25 0.25) 10.421633  0.9364839   0
(1 0.5 0 0 0)12.577408  1.2055979   0
(1 0.5 1 0 0) 9.097063  1.1769760   0
(1 0.5 0 1 0)10.523485  1.2238635   0
(1 0.5 0 0 1)10.440959  1.2031459   0

#With Poisson distribution#

lmer(tot.fruit~infl.treat+def.treat+(1|initial.size),family=poisson)->model2
summary(model2)

Generalized linear mixed model fit using Laplace
Formula: tot.fruit ~ infl.treat + def.treat + (1 | initial.size)
  Family: poisson(log link)
   AIC  BIC logLik deviance
  1073 1096 -530.5 1061
Random effects:
  Groups   NameVariance Std.Dev.
  initial.size (Intercept) 0.84201  0.91761
number of obs: 323, groups: initial.size, 292

Estimated scale (compare to  1 )  1.130529

Fixed effects:
 Estimate Std. Error z value Pr(>|z|)
(Intercept)  2.088650.10478  19.934  < 2e-16 ***
infl.treat1  0.106150.09412   1.128  0.25941
def.treat2  -0.373540.11398  -3.277  0.00105 **
def.treat3  -0.185660.12679  -1.464  0.14310
def.treat4  -0.106240.13451  -0.790  0.42962
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1

Correlation of Fixed Effects:
 (Intr) infl.1 df.tr2 df.tr3
infl.treat1 -0.418
def.treat2  -0.554  0.016
def.treat3  -0.599 -0.001  0.467
def.treat4  -0.623  0.055  0.460  0.548

estimable(model2,lsmeans,sim.lmer=TRUE,n.sim=1000)
Error in FUN(newX[, i], ...) :
   `param' has no names and does not match number of coefficients of  
model. Unable to construct coefficient vector

#End of CODE#

Any ideas as to why this is? It does say in the help of ?estimable  
that it should work for lmer objects.

Thanks in advance for any help,

Adriana Puentes

-- 
Adriana Puentes (MSc.)
PhD student
Plant Ecology
Evolutionary Biology Centre (EBC)
Uppsala University
Villavägen 14
75236, Uppsala,
Sweden

Webpage: http://www.vaxtbio.uu.se/resfold/puentes.htm
Work: (+46)18-471-2882
Mobile: (+46)(0)76-853-5861
Fax: (+46)18-553-419

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

Re: [R] producing output as *.spo (spss output format)

[Muenchen, Robert A (Bob)]

You probably don't want to spend time figuring out the .spo format. From
SPSS 16 on, that format is obsolete and replaced by the Unicode
XML-based .spv file. SPSS 16 users need a separate Legacy Viewer to read
.spo files. -Bob

=
Bob Muenchen (pronounced Min'-chen), 
Manager, Statistical Consulting Center 
U of TN Office of Information Technology
Stokely Management Center, Suite 200
916 Volunteer Blvd., Knoxville, TN 37996-0520
Voice: (865) 974-5230
FAX: (865) 974-4810
Email: [EMAIL PROTECTED]
Web: http://oit.utk.edu/scc
Map: http://www.utk.edu/maps 
News: http://listserv.utk.edu/archives/statnews.html
=

> -Original Message-
> From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED]
> project.org] On Behalf Of Ivan Uemlianin
> Sent: Saturday, November 10, 2007 6:04 AM
> To: r-help@r-project.org
> Subject: [R] producing output as *.spo (spss output format)
> 
> Dear All
> 
> I am considering moving from SPSS to R as my stats environment of
> choice.  I have read around and everything looks favourable.  There is
> just one issue on which I have been unable to find information.
> 
> Many clients ask me to send them output (tables, graphs, etc) as an
> spss
> output file (ie .spo).  I haven't asked them why, I've just said yes.
> I
> know R can produce graphics as nice as SPSS, and presumably they can
be
> output in some portable format for pasting into a word-processor
> document.  I need to find out why the client wants spo.  In the
> meantime
> let's assume they have a good reason.
> 
> Can R write .spo files?  Failing that does any one know of any spo
> writers that I could wire up to R (eg with some python gluecode)?
> Failing that any suggestions for overcoming the output hurdle would be
> welcome, as R looks very attractive (platform independent, easy to use
> and to automate, fast).
> 
> Best wishes
> 
> Ivan
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] White's test again

[Achim Zeileis]

On Fri, 9 Nov 2007, David Kaplan wrote:

> Hi all,
>
> It seems that I can get White's (HC3) test using MASS.   The syntax I
> used for the particular problem is
>
> anova(scireg3, white.adjust="hc3")

I don't think this is true. My guess is that you are using Anova() from
package "car".

> where scireg3 is an object from the lm function.  But, the anova summary
> table is all I get.  I don't get the new estimates or standard errors
> correcting for heteroskedasticity.  Is there a way to get that information?

Look at
  vignette("sandwich", package = "sandwich")
which discusses various kinds of sandwich covariances (including HC3). You
can plug these into:
  - coeftest() from package "lmtest": partial Wald tests for each
coefficients (as in the usual summary()),
  - waldtest() from package "lmtest": Wald tests for nested models
(similar to anova())
  - linear.hypothesis() from package "car": Wald tests for linear
hypotheses.

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.

[R] Estimate a logit of shares

[Trevor Marcel]

Michael,

I look at it.

Thanks for this first path.

Pierre-Olivier




  
_ 

l 
[[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] removing nodes from a dendrogram

[Jesse Lecy]

Hello,

I am trying to remove nodes from a dendrogram - specifically nodes that have
less than 3 members.  I am trying to create a function to do so, but I am
not sure how to reference the node inside of the function.  Here is the
code:

  rem.nodes <<- function(n) {
  if(!is.leaf(n) & attr(n, "members")<3) {

  REMOVE N ?? }}

reduced.dend <- dendrapply(dend1, rem.nodes)

Any ideas are appreciated.

Jesse

[[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] Rcmdr Mac OS X

[blemcke]

Sorry forgot subject.

B



> Von: Birgit Lemcke <[EMAIL PROTECTED]>
> Datum: 12. November 2007 13:26:55 GMT+01:00
> An: R Hilfe <[EMAIL PROTECTED]>
>
> Hello members of the mailinglist!
>
> I am running R 2.6.0 GUI 1.21 (4815) (4815) on a PPC with Mac OS X  
> 10.4.
>
> I just tried to use the R commander.
> I followed these instructions:
>
> 4) Start R (e.g. double-click the R icon in /Applications).
> 5) Start X11 (e.g. click the X11 button on the R Console or double- 
> click the X11 icon in /Applications/Utilities).
> 6) Close the X11 terminal (optional).
> 7) Type library(Rcmdr) in the R Console and press Return.
>
> X11 started and also the Rcmdr window opened, but I got this  
> messsages:
>
> Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht  
> gefunden #not found
> Error : .onAttach in 'attachNamespace' fehlgeschlagen #failed
> Fehler: Laden von Paket/Namensraum für 'Rcmdr' fehlgeschlagen
> #Error: Loading of Paket/Namensraum for Rcmdr failed
> Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht  
> gefunden #not found
> Error in get(x, envir = RcmdrEnv(), mode = mode, inherits = FALSE) :
>   Variable "messagesWindow" nicht gefunden
>   #variable "messagesWindow" not found
>
> When I tried to load the Prestige dataset from the car package I  
> got this:
>
> Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht  
> gefunden
>
> What am I doing wrong?
>
> Thanks in advance for your help.
>
> B
>
> Birgit Lemcke
> Institut für Systematische Botanik
> Zollikerstrasse 107
> CH-8008 Zürich
> Switzerland
> Ph: +41 (0)44 634 8351
> [EMAIL PROTECTED]
>
>
>
>
>

Birgit Lemcke
Altwiesenstrasse 89
CH-8051 Zürich
Switzerland
phone +0041 43 540 16 35
[EMAIL PROTECTED]


Hier noch ein kleiner Homepagetip: http://www.marlem-software.de/ 
kontakt.html

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

Re: [R] strange `nls' behaviour

[Duncan Murdoch]

On 11/12/2007 6:51 AM, Joerg van den Hoff wrote:
> I initially thought, this should better be posted to r-devel
> but alas! no response. 

I think the reason there was no response is that your example is too 
complicated.  You're doing a lot of strange things (fitfunc as a result 
of deriv, using as.name, as.call, as.formula, etc.)  You should simplify 
it down to isolate the bug.  Thats a lot of work, but you're the one in 
the best position to do it.  I'd say there's at least an even chance 
that the bug is in your code rather than in nls().

And 2.5.0 *is* ancient; please confirm the bug exists in R-patched if it 
turns out to be an R bug.

Duncan Murdoch




so  I  try  it  here.  sory  for  the
> lengthy explanation but it seems unavoidable. to quickly see
> the problem simply copy the litte example below and execute
> 
> f(n=5)
> 
> which  crashes. called with n !=  5 (and of course n>3 since
> there are 3 parameters in the model...) everything is as  it
> should be.
> 
> in detail:
> I  stumbled over the follwing _very_ strange behaviour/error
> when using `nls' which  I'm  tempted  (despite  the  implied
> "dangers") to call a bug:
> 
> I've  written a driver for `nls' which allows specifying the
> model and the data vectors using arbitrary  symbols.   these
> are  internally  mapped  to  consistent names, which poses a
> slight complication when using `deriv' to  provide  analytic
> derivatives. the following fragment gives the idea:
> 
> #-
> f <- function(n = 4) {
> 
>x <- seq(0, 5, length  = n)
> 
>y <- 2 * exp(-1*x) + 2; 
>y <- rnorm(y,y, 0.01*y)
> 
>model <- y ~ a * exp (-b*x) + c
> 
>fitfunc <- deriv(model[[3]], c("a", "b", "c"), c("a", "b", "c", "x"))
> 
>#"standard" call of nls:
>res1 <- nls(y ~ fitfunc(a, b, c, x), start = c(a=1, b=1, c=1))
> 
>call.fitfunc <- 
>c(list(fitfunc), as.name("a"), as.name("b"), as.name("c"), as.name("x"))
>call.fitfunc <- as.call(call.fitfunc)
>frml <- as.formula("y ~ eval(call.fitfunc)")
> 
>#"computed" call of nls:
>res2 <- nls(frml, start = c(a=1, b=1, c=1))
> 
>list(res1 = res1, res2 = res2)
> }
> #-
> 
> the  argument  `n'   defines  the number of (simulated) data
> points x/y which are going to be fitted by some model ( here
> y ~ a*exp(-b*x)+c )
> 
> the first call to `nls' is the standard way of calling `nls'
> when knowing all the variable and parameter names.
> 
> the second call (yielding `res2') uses a constructed formula
> in `frml' (which in this example is of course not necessary,
> but  in  the general case 'a,b,c,x,y' are not a priori known
> names).
> 
> now, here is the problem: the call
> 
> f(4)
> 
> runs fine/consistently, as does every call with n > 5.
> 
> BUT: for n = 5 (i.e. issuing f(5))
> 
> the second fit leads to the error message:
> 
> "Error in model.frame(formula, rownames, variables, varnames, extras, 
> extranames,  : 
>   invalid type (language) for variable 'call.fitfunc'"
> 
> I  cornered  this  to a spot in `nls' where a model frame is
> constructed in variable `mf'.  the parsing/constructing here
> seems  simply  to  be messed up for n = 5: `call.fitfunc' is
> interpreted as variable.
> 
> I,  moreover, empirically noted that the problem occurs when
> the total number of  parameters  plus  dependent/independent
> variables  equals  the number of data points (in the present
> example a,b,c,x,y).
> 
> so it is not the 'magic' number of 5 but rather the identity
> of data vector length and number of parameters+variables  in
> the model which leads to the problem.
> 
> this  is  with  2.5.0  (which  hopefully  is  not considered
> ancient) and MacOSX 10.4.10.
> 
> any ideas?
> 
> thanks
> 
> joerg
> 
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

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

[R] (no subject)

[Birgit Lemcke]

Hello members of the mailinglist!

I am running R 2.6.0 GUI 1.21 (4815) (4815) on a PPC with Mac OS X 10.4.

I just tried to use the R commander.
I followed these instructions:

4) Start R (e.g. double-click the R icon in /Applications).
5) Start X11 (e.g. click the X11 button on the R Console or double- 
click the X11 icon in /Applications/Utilities).
6) Close the X11 terminal (optional).
7) Type library(Rcmdr) in the R Console and press Return.

X11 started and also the Rcmdr window opened, but I got this messsages:

Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht  
gefunden #not found
Error : .onAttach in 'attachNamespace' fehlgeschlagen #failed
Fehler: Laden von Paket/Namensraum für 'Rcmdr' fehlgeschlagen
#Error: Loading of Paket/Namensraum for Rcmdr failed
Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht  
gefunden #not found
Error in get(x, envir = RcmdrEnv(), mode = mode, inherits = FALSE) :
   Variable "messagesWindow" nicht gefunden
   #variable "messagesWindow" not found

When I tried to load the Prestige dataset from the car package I got  
this:

Error in eval(expr, envir, enclos) : objekt ".__M__coerce" nicht  
gefunden

What am I doing wrong?

Thanks in advance for your help.

B

Birgit Lemcke
Institut für Systematische Botanik
Zollikerstrasse 107
CH-8008 Zürich
Switzerland
Ph: +41 (0)44 634 8351
[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] removing nodes from a dendrogram

[Jesse D Lecy]

>
> Hello,
>


 I am trying to remove nodes from a dendrogram - specifically nodes that
> have less than 3 members.  I am trying to create a function to do so, but I
> am not sure how to reference the node inside of the function.  Here is the
> code:
>
>   rem.nodes <<- function(n) {
>   if(!is.leaf(n) & attr(n, "members")<3) {
>
>   REMOVE N ?? }}
>
> reduced.dend <- dendrapply(dend1, rem.nodes)
>
> Any ideas are appreciated.
>
> Jesse
>
>
>

[[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] decompose () vs. stl ()

[Nuno Prista]

Hi,

I have been having trouble finding references that compare the classical
decomposition methods [decompose()] with the STL method of
Cleveland et al. (1990) [stl()].

Can anybody indicate some references that make comparisons between the two
methods and the advantage/disadvantages of each one in different situations?

Thanks in advance,

Nuno

__ 
Instituto de Oceanografia - FCUL, Portugal
Center for Quantitative Fisheries Ecology - ODU, USA

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

Re: [R] how to extract the original data of a glm object

[vito muggeo]

See and *read* the help file ?glm
the object returned by glm() includes the `data' component

hence:

aa<-glm(..)
aa$data

or also eval(aa$call$data)



leffgh ha scritto:
>   my function is 
>  glm(a~log(b)+c+d+e,family=binomial,data=f)->aa
>  
>   
>   I want to extract the original data set of aa. How to do it ?
>   
> You may suggest the model.frame() function. In fact ,i have tried it.
> model.frame returns a data frame of containing a,log(b) NOT b,c,d,e 
> I want to extract a data frame containing a,b,c,d,e,which is exactly the
> same as "f"
> How can I achieve this result?
> 
> 
> I want to do this because I need to extract the formular and act on another
> data set ,whose predict variables are the same as those of "f", but the
> response variable is different .

-- 

Vito M.R. Muggeo
Dip.to Sc Statist e Matem `Vianelli'
Università di Palermo
viale delle Scienze, edificio 13
90128 Palermo - ITALY
tel: 091 6626240
fax: 091 485726/485612

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

Re: [R] strange `nls' behaviour

[Katharine Mullen]

I can confirm this behavior on R-2.6.0 but don't have time to look into it
further at the moment.

On Mon, 12 Nov 2007, Joerg van den Hoff wrote:

>
> I initially thought, this should better be posted to r-devel
> but alas! no response. so  I  try  it  here.  sory  for  the
> lengthy explanation but it seems unavoidable. to quickly see
> the problem simply copy the litte example below and execute
>
> f(n=5)
>
> which  crashes. called with n !=  5 (and of course n>3 since
> there are 3 parameters in the model...) everything is as  it
> should be.
>
> in detail:
> I  stumbled over the follwing _very_ strange behaviour/error
> when using `nls' which  I'm  tempted  (despite  the  implied
> "dangers") to call a bug:
>
> I've  written a driver for `nls' which allows specifying the
> model and the data vectors using arbitrary  symbols.   these
> are  internally  mapped  to  consistent names, which poses a
> slight complication when using `deriv' to  provide  analytic
> derivatives. the following fragment gives the idea:
>
> #-
> f <- function(n = 4) {
>
>x <- seq(0, 5, length  = n)
>
>y <- 2 * exp(-1*x) + 2;
>y <- rnorm(y,y, 0.01*y)
>
>model <- y ~ a * exp (-b*x) + c
>
>fitfunc <- deriv(model[[3]], c("a", "b", "c"), c("a", "b", "c", "x"))
>
>#"standard" call of nls:
>res1 <- nls(y ~ fitfunc(a, b, c, x), start = c(a=1, b=1, c=1))
>
>call.fitfunc <-
>c(list(fitfunc), as.name("a"), as.name("b"), as.name("c"), as.name("x"))
>call.fitfunc <- as.call(call.fitfunc)
>frml <- as.formula("y ~ eval(call.fitfunc)")
>
>#"computed" call of nls:
>res2 <- nls(frml, start = c(a=1, b=1, c=1))
>
>list(res1 = res1, res2 = res2)
> }
> #-
>
> the  argument  `n'   defines  the number of (simulated) data
> points x/y which are going to be fitted by some model ( here
> y ~ a*exp(-b*x)+c )
>
> the first call to `nls' is the standard way of calling `nls'
> when knowing all the variable and parameter names.
>
> the second call (yielding `res2') uses a constructed formula
> in `frml' (which in this example is of course not necessary,
> but  in  the general case 'a,b,c,x,y' are not a priori known
> names).
>
> now, here is the problem: the call
>
> f(4)
>
> runs fine/consistently, as does every call with n > 5.
>
> BUT: for n = 5 (i.e. issuing f(5))
>
> the second fit leads to the error message:
>
> "Error in model.frame(formula, rownames, variables, varnames, extras, 
> extranames,  :
>   invalid type (language) for variable 'call.fitfunc'"
>
> I  cornered  this  to a spot in `nls' where a model frame is
> constructed in variable `mf'.  the parsing/constructing here
> seems  simply  to  be messed up for n = 5: `call.fitfunc' is
> interpreted as variable.
>
> I,  moreover, empirically noted that the problem occurs when
> the total number of  parameters  plus  dependent/independent
> variables  equals  the number of data points (in the present
> example a,b,c,x,y).
>
> so it is not the 'magic' number of 5 but rather the identity
> of data vector length and number of parameters+variables  in
> the model which leads to the problem.
>
> this  is  with  2.5.0  (which  hopefully  is  not considered
> ancient) and MacOSX 10.4.10.
>
> any ideas?
>
> thanks
>
> joerg
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] how to extract the original data of a glm object

[Henrique Dallazuanna]

Try:

aa$data

On 12/11/2007, leffgh <[EMAIL PROTECTED]> wrote:
>
>  my function is
>  glm(a~log(b)+c+d+e,family=binomial,data=f)->aa
>
>
>  I want to extract the original data set of aa. How to do it ?
>
> You may suggest the model.frame() function. In fact ,i have tried it.
> model.frame returns a data frame of containing a,log(b) NOT b,c,d,e
> I want to extract a data frame containing a,b,c,d,e,which is exactly the
> same as "f"
> How can I achieve this result?
>
>
> I want to do this because I need to extract the formular and act on another
> data set ,whose predict variables are the same as those of "f", but the
> response variable is different .
> --
> View this message in context: 
> http://www.nabble.com/how-to-extract-the-original-data-of-a-glm-object-tf4790525.html#a13704122
> 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.
>


-- 
Henrique Dallazuanna
Curitiba-Paraná-Brasil
25° 25' 40" S 49° 16' 22" O

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