[R] Error of Cross Validation

[zhu yao]

Dear R users:

Recently, I tried to write a program to calculate cross-validated predicted
value.
My sources are as follows. However, the R reported an error.
Could you please check the sources? Thanks.

set.seed(100)
x<-rnorm(100)
y<-sample(rep(0:1,50),replace=T)
dat<-data.frame(x,y)

library(rms)

fito<-lrm(y~x)
preo<-predict(fito)

pre<-matrix(NA,nrow=100,ncol=200)

for (i in 1:200)
{
sam<-sample(1:nrow(dat))
sam<-split(sam,1:10)
 for (j in 1:10)
{
fit<-lrm(y~x,data=dat[-sam[[j]],])
pre[sam[[j]],i]<-predict(fit,data=dat[sam[[j]],])
}
}








*Yao Zhu*
*Department of Urology
Fudan University Shanghai Cancer Center
Shanghai, China*

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

[Ungku Akashah]

HELLO, anybody... could you help me to check the below coding for volcano.
what is the mistake?
what the plot could not display? 







#volcano_plot.r
#
#Author:Amsha Nahid, Jairus Bowne, Gerard Murray
#Purpose:Produces a volcano plot
#
#Input:Data matrix as specified in Data-matrix-format.pdf
#Output:Plots log2(fold change) vs log10(t-test P-value)
#
#Notes:Group value for control must be alphanumerically first
#  Script will return an error if there are more than 2 groups

#
#Load the data matrix
#
# Read in the .csv file
data<-read.csv("input5.csv", sep=",", row.names=1, header=TRUE)
# Get groups information
groups<-data[,1]
# Get levels for groups
grp_levs<-levels(groups)
if (length(levels(groups)) > 2)
print("Number of groups is greater than 2!") else {

#
#Split the matrix by group
#
new_mats<-c()
for (ii in 1:length(grp_levs))
new_mats[ii]<-list(data[which(groups==levels(groups)[ii]),])

#
#Calculate the means
#
# For each matrix, calculate the averages per column
submeans<-c()
# Preallocate a matrix for the means
means<-matrix(
nrow = 2,
ncol = length(colnames(data[,-1])),
dimnames = list(grp_levs,colnames(data[,-1]))
)
# Calculate the means for each variable per sample
for (ii in 1:length(new_mats))
{submeans[ii]<-list(apply(new_mats[[ii]][,-1],2,mean,na.rm=TRUE))
means[ii,]<-submeans[[ii]]}

#
#Calculate the fold change
#
folds<-matrix(
nrow=length(means[,1]),
ncol=length(means[1,]),
dimnames=list(rownames(means),colnames(means))
)
for (ii in 1:length(means[,1]))
for (jj in 1:length(means[1,]))
folds[ii,jj]<-means[ii,jj]/means[1,jj]

#
#t-test P value data
#

pvals<-matrix(nrow=ncol(data[,-1]),ncol=1,dimnames=list(colnames(data[-1]),"P-Value"))


#
#Perform the t-Test
#
for(ii in 1:nrow(pvals)) {
pvals[ii,1]<-t.test(new_mats[[1]][,ii+1],new_mats[[2]][,ii+1])$p.value
}

m<-length(pvals)
x_range<-range(c(
min(
range(log2(folds[2,])),
range(c(-1.5,1.5))
),
max(
range(log2(folds[2,])),
range(c(-1.5,1.5))
)
))
y_range<-range(c(
min(range(-log10(pvals)),
range(c(0,2))
),
max(range(-log10(pvals)),
range(c(0,2))
)
))

#
#Plot data
#
# Define a function, since it's rather involved
volcano_plot<-function(fold, pval)
{plot(x_range, # x-dim 
y_range,   # y-dim
type="n",  # empty plot
xlab="log2 Fold Change",   # x-axis title
ylab="-log10 t-Test P-value",  # y-axis title
main="Volcano Plot",   # plot title
)
abline(h=-log10(0.05),col="green",lty="44")# horizontal line at 
P=0.05
abline(v=c(-1,1),col="violet",lty="1343")  # vertical lines at 
2-fold
# Plot points based on their values:
for (ii in 1:m)
# If it's below 0.05, we're not overly interested: purple.
if (-log10(pvals[ii])>(-log10(0.05))) {
# Otherwise, more checks;
# if it's greater than 2-fold decrease: blue
if (log2(folds[2,][ii])>(-1)) {
# If it's significant but didn't change much: orange
if (log2(folds[2,][ii])<1) {
points(
log2(folds[2,][ii]),
-log10(pvals[ii]),
col="orange",
pch=20
)
# Otherwise, greater than 2-fold increase: red
} else {
points(
log2(folds[2,][ii]), 
-log10(pvals[ii]),
col="red",
pch=20
)
}
} else {
points(
log2(folds[2,][ii]), 
-log10(pvals[ii]),
col="blue",
pch=20
)
}
} else {
points(
log2(folds[2,][ii]),
-log10(pvals[ii]),
col="purple",

Re: [R] Trouble with Optimization in "Alabama" Package

[mcgrete]

Hello Erik and Ravi,

I am curious if Erick or Ravi have a solution to Erick's question.

I have a similar problem.  I wish to use auglag to solve a nonlinear
optimization problem.  The cost function, grad function, and constraints are
all function of some parameters that are dependent upon the system I wish to
analyze.  For any given system, these parameters are fixed, e.g. x,y,z are
parameters that are problem or system dependent.  The variable f is a vector
that I wish to find such that my cost function is maximized or minimized.

I have used auglag without any issues as follows.  At first, I called auglag
without trying to provide arguments for x,y,z, but rather only the initial
parameter fo.  My functions (fn, gr, hin, heq, etc) were all defined in the
format of:
foo<-function(f){
# code...
}
The parameters x,z,z were obtained from what I understand to be the global
environment.  I was NOT calling 'auglag' from a function, but simply from
the command line in 'Rkward'. In other words, parameters x,y,z were
calculated outside of any functions but were utilized in the function
fn,gr,hin,heq,etc.

Now, I have constructed a function, foobar, that utilizes inputs a,b,c that
are used to calculate parameters x,y,z.  And within foobar, I attempt to
call 'auglag' to find an optimum solution.  For example:
foobar<-function(a,b,c){
# code to calculate x,y,z as function of a,b,c
ans<-auglag(par=fo,fn=foo1,gr=foo2,hin=foo3,heq=foo4)
return(ans)
}

However, auglag appears to utilize x,y,z from the global environment rather
than my local variables x,y,z inside foobar.

I have attempted to construct functions fn,gr,hin,heq,... by several ways
without success.
For example
foo<-function(f,x,y,z){code...},
ans<-auglag(par=fo,fn=foo1,gr=foo2,hin=foo3,heq=foo4), or
ans<-auglag(par=c(fo,x,y,z),fn=foo1(,gr=foo2,hin=foo3,heq=foo4)

and I tried
foo<-function(f,...){code...},
ans<-auglag(par=fo,fn=foo1,gr=foo2,hin=foo3,heq=foo4), or

foo<-function(f,...){code...},
ans<-auglag(par=c(fo,x,y,z),fn=foo1(,gr=foo2,hin=foo3,heq=foo4)

All failed.

Is it simply not feasible to use 'auglag' in this manner?  Or, have I missed
how to pass arguments to auglag, which in turn can pass these arguments to
fn,gr,hin,heq,...

Any help would be greatly appreciated.

Regards,
Tim



--
View this message in context: 
http://r.789695.n4.nabble.com/Trouble-with-Optimization-in-Alabama-Package-tp2548801p3610402.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] Skyline plots from several trees in newick format

[zale]

No, unfortunately not. It works perfectly fine with other kinds of files.
Everything, but these ones...The program that generates them does not have
any helpful specifications regarding this format ..*.trees

On Mon, Jun 20, 2011 at 5:55 AM, Bhoom [via R] <
ml-node+3610421-1897267921-239...@n4.nabble.com> wrote:

> Does read.nexus work for you?
>
> Bhoom
>
> On Sun, Jun 19, 2011 at 5:04 PM, Andra Tolbus <[hidden 
> email]>wrote:
>
>
> > Dear all,
> >
> >  I am trying to create a consensus skyline plot using the "ape"
> > package(newbie).
> >  I have a nexus file that looks like the one from above containing
> many
> > trees.  (my_file.trees)
> >
> >  #NEXUS
> > begin trees;  [Treefile generated by sim_coal.exe (Laurent Excoffier)]
> >
> > tree true_tree_1 = [&U] (((13.1:6, 43.1:5):12, 28.1:14):284, (21.1:0,
>
> > (20.1:1, 4.1:0):0):1, ((37.1:0, (23.1:0, 29.1:0):0):0, 15.1:0):3):7,
> > (36.1:3, (40.1:0, 26.1:1):2):1):2, ((16.1:0, 33.1:0):3, (10.1:0,
> > 18.1:0):2):8):46, 6.1:1, (11.1:0, 34.1:0):0):2, 2.1:1):19, (((38.1:4,
>
> > (35.1:4, ((8.1:0, 42.1:0):1, 22.1:0):3):2):0, ((5.1:2, (27.1:0,
> > 39.1:0):0):0, (17.1:0, 30.1:0):1):4):8, ((31.1:0, ((41.1:0, (7.1:0,
> > 25.1:1):0):0, 3.1:0):3):0, (12.1:0, (14.1:0, 1.1:0):0):0):13):10):14,
> > (((24.1:0, 9.1:0):0, 32.1:0):1, 19.1:0):25):20):206);
> > tree true_tree_2 = [&U] 37.1:4, ((36.1:0, 25.1:0):0, 39.1:0):2):0,
> > 33.1:1):2, (1.1:0, 43.1:0):1):81, (21.1:1, (38.1:1, 14.1:1):0):3,
> > (((29.1:0, 35.1:0):1, 4.1:1):1, ((41.1:1, (24.1:0, ((8.1:0, 32.1:0):0,
> > 7.1:0):0):1):0, (10.1:2, (6.1:1, (26.1:1, (15.1:0,
> > 16.1:0):0):0):0):2):0):4):1, (28.1:1, (3.1:1, 13.1:0):0):7):12,
> (((22.1:5,
> > (((2.1:1, 34.1:1):0, 19.1:0):3, (18.1:0, 42.1:0):3):0):0, 17.1:0):3,
> > (27.1:14, (((30.1:0, 23.1:0):1, 31.1:0):1, 5.1:2):5):3):5):9, ((11.1:0,
> > 9.1:1):0, (12.1:0, (20.1:0, 40.1:1):1):3):23):67);
> >
> > I haven't managed to read this file content using the *read.tree*
> function
> > ..
> >
> >  Error in if (tp[3] != "") obj$node.label <- tp[3] :
> >  missing value where TRUE/FALSE needed
> >
> > Does anyone know first of all, why I cannot read this file (not even with
>
> > read.nexus function)? And secondly how can one generate an "average"
> > skyline
> > plot from several trees?
> >
> > Thank you in advance!
> > Cheers!
> >
> >[[alternative HTML version deleted]]
> >
> > __
> > [hidden email] 
> > mailing list
> > https://stat.ethz.ch/mailman/listinfo/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]]
>
> __
> [hidden email] mailing 
> list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>
> --
>  If you reply to this email, your message will be added to the discussion
> below:
>
> http://r.789695.n4.nabble.com/Skyline-plots-from-several-trees-in-newick-format-tp3609929p3610421.html
>  To unsubscribe from Skyline plots from several trees in newick format, click
> here.
>
>


--
View this message in context: 
http://r.789695.n4.nabble.com/Skyline-plots-from-several-trees-in-newick-format-tp3609929p3610452.html
Sent from the R help mailing list archive at Nabble.com.
[[alternative HTML version deleted]]

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

[R] Model Formulae Evaluation

[albeam]

Hello All,

I have searched but haven't been able to find an answer to this question.
I'm writing a function that needs to be able evaluate the right-hand side of
a general non-linear formula for different parameter values, similar to what
nls() does. I looked through the nls() code but it isn't clear what is going
on. Could someone point me to a reference or a way I might do this?

Thanks

--
View this message in context: 
http://r.789695.n4.nabble.com/Model-Formulae-Evaluation-tp3610328p3610328.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] (no subject)

[Bhoom Suktitipat]

Hi Jim,

Just want to make sure that I understand correctly.What do you mean by SNP =
Single Nucleotide Polymorphism?

Bhoom

On Sun, Jun 19, 2011 at 5:57 PM, Jim Silverton wrote:

> Hello, Does anyone has software to simulate SNP data? Specifically I would
> like any package that produces  2 x 3 SNP data.
>
> --
> Thanks,
> Jim.
>
>[[alternative HTML version deleted]]
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

[[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] Skyline plots from several trees in newick format

[Bhoom Suktitipat]

Does read.nexus work for you?

Bhoom

On Sun, Jun 19, 2011 at 5:04 PM, Andra Tolbus wrote:

> Dear all,
>
>  I am trying to create a consensus skyline plot using the "ape"
> package(newbie).
>  I have a nexus file that looks like the one from above containing many
> trees.  (my_file.trees)
>
>  #NEXUS
> begin trees;  [Treefile generated by sim_coal.exe (Laurent Excoffier)]
>
> tree true_tree_1 = [&U] (((13.1:6, 43.1:5):12, 28.1:14):284, (21.1:0,
> (20.1:1, 4.1:0):0):1, ((37.1:0, (23.1:0, 29.1:0):0):0, 15.1:0):3):7,
> (36.1:3, (40.1:0, 26.1:1):2):1):2, ((16.1:0, 33.1:0):3, (10.1:0,
> 18.1:0):2):8):46, 6.1:1, (11.1:0, 34.1:0):0):2, 2.1:1):19, (((38.1:4,
> (35.1:4, ((8.1:0, 42.1:0):1, 22.1:0):3):2):0, ((5.1:2, (27.1:0,
> 39.1:0):0):0, (17.1:0, 30.1:0):1):4):8, ((31.1:0, ((41.1:0, (7.1:0,
> 25.1:1):0):0, 3.1:0):3):0, (12.1:0, (14.1:0, 1.1:0):0):0):13):10):14,
> (((24.1:0, 9.1:0):0, 32.1:0):1, 19.1:0):25):20):206);
> tree true_tree_2 = [&U] 37.1:4, ((36.1:0, 25.1:0):0, 39.1:0):2):0,
> 33.1:1):2, (1.1:0, 43.1:0):1):81, (21.1:1, (38.1:1, 14.1:1):0):3,
> (((29.1:0, 35.1:0):1, 4.1:1):1, ((41.1:1, (24.1:0, ((8.1:0, 32.1:0):0,
> 7.1:0):0):1):0, (10.1:2, (6.1:1, (26.1:1, (15.1:0,
> 16.1:0):0):0):0):2):0):4):1, (28.1:1, (3.1:1, 13.1:0):0):7):12, (((22.1:5,
> (((2.1:1, 34.1:1):0, 19.1:0):3, (18.1:0, 42.1:0):3):0):0, 17.1:0):3,
> (27.1:14, (((30.1:0, 23.1:0):1, 31.1:0):1, 5.1:2):5):3):5):9, ((11.1:0,
> 9.1:1):0, (12.1:0, (20.1:0, 40.1:1):1):3):23):67);
>
> I haven't managed to read this file content using the *read.tree* function
> ..
>
>  Error in if (tp[3] != "") obj$node.label <- tp[3] :
>  missing value where TRUE/FALSE needed
>
> Does anyone know first of all, why I cannot read this file (not even with
> read.nexus function)? And secondly how can one generate an "average"
> skyline
> plot from several trees?
>
> Thank you in advance!
> Cheers!
>
>[[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] analysis with the data from mysql database

[amrita gs]

Hi everyone,

 I have certain values retrieved from mysql database.How can i do
certain analysis like histogram, correlational analysis tc using this data,I
tried it. But when i tried to plot a histogram it actually showed error the
data is not numeric evenif it was stored as integer data type in mysql db.

Please anyone help. Can anyone specify the steps for doing data analysis
after retrieving data
from the database.

[[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] For loop by factor.

[Bill.Venables]

If I understand you correctly, you are trying to find the cumulative maximum 
from the end within each level of the factor.  If this is what you are trying 
to do, then here is one way you might like to do it.  First, define the 
function:

> cumMax <- function(x) Reduce(max, x, right = TRUE, accumulate = TRUE)

Here is how you might use it:

> test
  A B
1 a 3
2 a 2
3 a 1
4 b 3
5 b 2
6 c 2
7 c 3
8 c 1
9 c 1
> (test <- within(test, B <- ave(B, A, FUN = cumMax))
  A B
1 a 3
2 a 2
3 a 1
4 b 3
5 b 2
6 c 3
7 c 3
8 c 1
9 c 1 

-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On 
Behalf Of Christopher Peters
Sent: Monday, 20 June 2011 6:21 AM
To: r-help@r-project.org
Subject: [R] For loop by factor.

I have a data.frame as follows:

a  3
a  2
a  1
b  3
b  2
c  2
c  3
c  1
c  1

Each factor (a, b, c) should be monotonically decreasing, notice that factor
'c' is not.

I could use some help to figure out how to form a logical structure (mostly
just syntax), that will check each 'next value' for each factor to see if it
is less than the previous value.  If it is less than the previous value, do
nothing, else subtract 'next value' from 'current value', add that amount to
the starting value and each previous value to the 'next value' is greater
than 'previous value'.

So basically the data.frame would look like:

a 3
a 2
a 1
b 3
b 2
c 3
c 3
c 1
c 1

Thanks for your help!
__
Christopher P. Peters
Research Associate
Center for Energy Studies
http://www.enrg.lsu.edu/staff/peters
Energy, Coast & Environment Building
Louisiana State University
Baton Rouge, LA 70803
Telephone: 225-578-4400
Fax: 225-578-4541

[[alternative HTML version deleted]]

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

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

Re: [R] For loop by factor.

[jim holtman]

try this:

>  test <- data.frame(A=c("a", "a", "a", "b", "b", "c", "c", "c", "c"), 
> B=c(3,2,1,3,2,2,3,1,1))
>  test
  A B
1 a 3
2 a 2
3 a 1
4 b 3
5 b 2
6 c 2
7 c 3
8 c 1
9 c 1
>  # determine which group is not decreasing
>  tapply(test$B, test$A, function(x) any(diff(x) > 0))
a b c
FALSE FALSE  TRUE
>


On Sun, Jun 19, 2011 at 4:20 PM, Christopher Peters  wrote:
> I have a data.frame as follows:
>
> a  3
> a  2
> a  1
> b  3
> b  2
> c  2
> c  3
> c  1
> c  1
>
> Each factor (a, b, c) should be monotonically decreasing, notice that factor
> 'c' is not.
>
> I could use some help to figure out how to form a logical structure (mostly
> just syntax), that will check each 'next value' for each factor to see if it
> is less than the previous value.  If it is less than the previous value, do
> nothing, else subtract 'next value' from 'current value', add that amount to
> the starting value and each previous value to the 'next value' is greater
> than 'previous value'.
>
> So basically the data.frame would look like:
>
> a 3
> a 2
> a 1
> b 3
> b 2
> c 3
> c 3
> c 1
> c 1
>
> Thanks for your help!
> __
> Christopher P. Peters
> Research Associate
> Center for Energy Studies
> http://www.enrg.lsu.edu/staff/peters
> Energy, Coast & Environment Building
> Louisiana State University
> Baton Rouge, LA 70803
> Telephone: 225-578-4400
> Fax: 225-578-4541
>
>        [[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
Data Munger Guru

What is the problem that 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] extract data from a column

[Nanami]

I figured it out:


x<-sub("^.*:([[:digit:]]+)..([[:digit:]]+).*", "\\1 \\2", CTSS$V4)

--
View this message in context: 
http://r.789695.n4.nabble.com/extract-data-from-a-column-tp3609890p3610147.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] extract data from a column

[Nanami]

So if I am given some data that look like this:

> head(CTSS)
V1 V2 V3V4   V5 V6  V7   
1 chr1 564563 564598 chr1:564588..564589,+ 1336  
2 chr1 564620 564649 chr1:564644..564645,+   94 
3 chr1 565369 565404 chr1:565371..565372,+  217  
4 chr1 565463 565541 chr1:565480..565481,+ 1214 
5 chr1 565653 565697 chr1:565662..565663,+ 1031  
6 chr1 565861 565922 chr1:565883..565884,+  316 

what I would like to do is to obtain two columns from column V4 like:
start_coord  stop_coord
564588564589
564644564645

I will try and see if it works what you suggested. 
Thank you,
Best,
Nanami




--
View this message in context: 
http://r.789695.n4.nabble.com/extract-data-from-a-column-tp3609890p3610049.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] extract data from a column

[Sam Albers]

If I understand what you want (which I may very well not) you could use
something like this:

If this is an example of your type of data:
564589,+

substr(x, 1, 6)
as.numeric(x)

Please try to post something more thorough if you would like a better
answer.

Sam

--
View this message in context: 
http://r.789695.n4.nabble.com/extract-data-from-a-column-tp3609890p3610030.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] (no subject)

[Jim Silverton]

Hello, Does anyone has software to simulate SNP data? Specifically I would
like any package that produces  2 x 3 SNP data.

-- 
Thanks,
Jim.

[[alternative HTML version deleted]]

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

Re: [R] hello about proxy configuration

[Duncan Murdoch]

On 19/06/2011 7:29 PM, p_conno...@slingshot.co.nz wrote:

Quoting 王海生:


dear everyone
system:windows XP
R2.13.0

I download the windows binary from website and successfully install
it ,because I use a proxy server ,when I follow the instrction as
follows:
I set a system property  R_HOME=C:\Program Files\R\R-2.13.0
" R_HOME\bin\i386\Rgui.exe  http_proxy=http://211.83.105.140:808/";


It's more likely to be 8080, not 808.

If that doesn't work, there is an easier way.  When you install the Windows
binary, there is an option called 'Internet access'.  Change from the default
to internet2, and R will then use your browser settings.  Assuming your
browser
works, R will have the correct proxy setting.


You can make this change any time, not only at install time. 
Temporarily in a session using the setInternet2 function, permanently by 
editing the shortcut to include the --internet2 option on the command line.


Duncan Murdoch



HTH




error:: '\i' is an unrecognized escape in character string starting
"R_HOME\bin\i"
what is wrong with my command? I just want to know is there any other
way such as mannually edit a configuration file to set the proxy?
kind regards
   edwin_uestc
Center of E-Health technology and Engineering in UESTC


[[alternative HTML version deleted]]

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


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


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

Re: [R] hello about proxy configuration

[p_connolly]

Quoting 王海生 :


dear everyone
system:windows XP
R2.13.0

I download the windows binary from website and successfully install 
it ,because I use a proxy server ,when I follow the instrction as 
follows:

I set a system property  R_HOME=C:\Program Files\R\R-2.13.0
" R_HOME\bin\i386\Rgui.exe  http_proxy=http://211.83.105.140:808/";


It's more likely to be 8080, not 808.

If that doesn't work, there is an easier way.  When you install the Windows
binary, there is an option called 'Internet access'.  Change from the default
to internet2, and R will then use your browser settings.  Assuming your 
browser

works, R will have the correct proxy setting.

HTH



error:: '\i' is an unrecognized escape in character string starting 
"R_HOME\bin\i"
what is wrong with my command? I just want to know is there any other 
way such as mannually edit a configuration file to set the proxy?

kind regards
  edwin_uestc
Center of E-Health technology and Engineering in UESTC


[[alternative HTML version deleted]]

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


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

Re: [R] For loop by factor.

[Sarah Goslee]

This works, but I'm still hunting for a more elegant final step:
> test <- data.frame(A=c("a", "a", "a", "b", "b", "c", "c", "c", "c"), 
> B=c(3,2,1,3,2,2,3,1,1))
> test2 <- lapply(split(test$B, test$A), sort, dec=TRUE)
> test3 <- data.frame(A=rep(names(test2), times=lapply(test2, length)), 
> B=unlist(test2))

It will also group data by factor, which was already done in your example.

Sarah

On Sun, Jun 19, 2011 at 4:20 PM, Christopher Peters  wrote:
> I have a data.frame as follows:
>
> a  3
> a  2
> a  1
> b  3
> b  2
> c  2
> c  3
> c  1
> c  1
>
> Each factor (a, b, c) should be monotonically decreasing, notice that factor
> 'c' is not.
>
> I could use some help to figure out how to form a logical structure (mostly
> just syntax), that will check each 'next value' for each factor to see if it
> is less than the previous value.  If it is less than the previous value, do
> nothing, else subtract 'next value' from 'current value', add that amount to
> the starting value and each previous value to the 'next value' is greater
> than 'previous value'.
>
> So basically the data.frame would look like:
>
> a 3
> a 2
> a 1
> b 3
> b 2
> c 3
> c 3
> c 1
> c 1
>
-- 
Sarah Goslee
http://www.functionaldiversity.org

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

[R] roundoff error in recursion

[Zhongyi Yuan]

Dear R users,

I run into the following problem and hope someone can help me out. Thank you
in advance for your time.

I have a function defined recursively as follows:
recFun = function(k,x){ #i is less than 10. x can be any real number. I need
the value of recFun(k=1,2,3,4, x=0.12).
 if(k==1)
  integrand = function(z) (1 + pi*exp(z+0.9*x+0.012))^(-1.5) * dnorm(z)
 else
  integrand = function(z) (1 + pi*exp(z+0.9*x+0.012))^(-1.5) *
recFun(k-1,z+0.9*x+0.012) * dnorm(z)
 integrate(integrand, -Inf, Inf)
}

The goal of the function is to evaluate a quantity satisfying the recursive
relation that
$$I(1,x)=\int_{-\infty }^{\infty }\left( 1+\pi
e^{z+0.9x+0.012}\right)^{-1.5}\Phi (dz),$$ and
$$I(k,x)=\int_{-\infty }^{\infty }\left( 1+\pi e^{z+0.9x+0.012}\right)
^{-1.5} I(k-1,z+0.9x+0.012)\Phi (dz),$$
where $\Phi$ is the distribution function of the standard normal
distribution.

But for k>=2, a roundoff error occurs:
> recFun(2,0.12)
Error in integrate(integrand, -Inf, Inf) : roundoff error was detected

Can someone suggest a solution to this problem? Thanks a lot.

Best,
Zhongyi

[[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] Skyline plots from several trees in newick format

[Andra Tolbus]

Dear all,

  I am trying to create a consensus skyline plot using the "ape"
package(newbie).
  I have a nexus file that looks like the one from above containing many
trees.  (my_file.trees)

 #NEXUS
begin trees;  [Treefile generated by sim_coal.exe (Laurent Excoffier)]

tree true_tree_1 = [&U] (((13.1:6, 43.1:5):12, 28.1:14):284, (21.1:0,
(20.1:1, 4.1:0):0):1, ((37.1:0, (23.1:0, 29.1:0):0):0, 15.1:0):3):7,
(36.1:3, (40.1:0, 26.1:1):2):1):2, ((16.1:0, 33.1:0):3, (10.1:0,
18.1:0):2):8):46, 6.1:1, (11.1:0, 34.1:0):0):2, 2.1:1):19, (((38.1:4,
(35.1:4, ((8.1:0, 42.1:0):1, 22.1:0):3):2):0, ((5.1:2, (27.1:0,
39.1:0):0):0, (17.1:0, 30.1:0):1):4):8, ((31.1:0, ((41.1:0, (7.1:0,
25.1:1):0):0, 3.1:0):3):0, (12.1:0, (14.1:0, 1.1:0):0):0):13):10):14,
(((24.1:0, 9.1:0):0, 32.1:0):1, 19.1:0):25):20):206);
tree true_tree_2 = [&U] 37.1:4, ((36.1:0, 25.1:0):0, 39.1:0):2):0,
33.1:1):2, (1.1:0, 43.1:0):1):81, (21.1:1, (38.1:1, 14.1:1):0):3,
(((29.1:0, 35.1:0):1, 4.1:1):1, ((41.1:1, (24.1:0, ((8.1:0, 32.1:0):0,
7.1:0):0):1):0, (10.1:2, (6.1:1, (26.1:1, (15.1:0,
16.1:0):0):0):0):2):0):4):1, (28.1:1, (3.1:1, 13.1:0):0):7):12, (((22.1:5,
(((2.1:1, 34.1:1):0, 19.1:0):3, (18.1:0, 42.1:0):3):0):0, 17.1:0):3,
(27.1:14, (((30.1:0, 23.1:0):1, 31.1:0):1, 5.1:2):5):3):5):9, ((11.1:0,
9.1:1):0, (12.1:0, (20.1:0, 40.1:1):1):3):23):67);

I haven't managed to read this file content using the *read.tree* function
..

  Error in if (tp[3] != "") obj$node.label <- tp[3] :
  missing value where TRUE/FALSE needed

Does anyone know first of all, why I cannot read this file (not even with
read.nexus function)? And secondly how can one generate an "average" skyline
plot from several trees?

Thank you in advance!
Cheers!

[[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] For loop by factor.

[Christopher Peters]

I have a data.frame as follows:

a  3
a  2
a  1
b  3
b  2
c  2
c  3
c  1
c  1

Each factor (a, b, c) should be monotonically decreasing, notice that factor
'c' is not.

I could use some help to figure out how to form a logical structure (mostly
just syntax), that will check each 'next value' for each factor to see if it
is less than the previous value.  If it is less than the previous value, do
nothing, else subtract 'next value' from 'current value', add that amount to
the starting value and each previous value to the 'next value' is greater
than 'previous value'.

So basically the data.frame would look like:

a 3
a 2
a 1
b 3
b 2
c 3
c 3
c 1
c 1

Thanks for your help!
__
Christopher P. Peters
Research Associate
Center for Energy Studies
http://www.enrg.lsu.edu/staff/peters
Energy, Coast & Environment Building
Louisiana State University
Baton Rouge, LA 70803
Telephone: 225-578-4400
Fax: 225-578-4541

[[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] Transfer my package in CRAN

[fadia bougacha]

Hello,


I found a problem to transfer my package on the site
ftp://cran.rproject.org/incoming/ .

Access is not allowed. Please,could you  give me the settings for connecting
to your server to transfer my package using the FileZilla.


Best Regards

Fadia

[[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] extract data from a column

[Nanami 13]

Hi all,
I have a column that has the following format:
 chr1:564588..564589,+  and I want to extract only the coordinates; I have
tried writing a regular expression but I couldn't figure out how I should
write it. Does anyone know?
Thank you,
Best,
Nanami

[[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] hello about proxy configuration

[王海生]

dear everyone
system:windows XP
R2.13.0

I download the windows binary from website and successfully install it ,because 
I use a proxy server ,when I follow the instrction as follows:
I set a system property  R_HOME=C:\Program Files\R\R-2.13.0
" R_HOME\bin\i386\Rgui.exe  http_proxy=http://211.83.105.140:808/";
error:: '\i' is an unrecognized escape in character string starting 
"R_HOME\bin\i"
what is wrong with my command? I just want to know is there any other way such 
as mannually edit a configuration file to set the proxy?
kind regards
  edwin_uestc
Center of E-Health technology and Engineering in UESTC


[[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] Multivariate HPD credible volume -- is it computable in R?

[Mike Marchywka]

> Date: Sun, 19 Jun 2011 19:06:20 +1200
> From: agw1...@gmail.com
> To: r-help@r-project.org
> Subject: [R] Multivariate HPD credible volume -- is it computable in R?
> 
> Hi all,
> 
> I'm new to the list and am hoping to get some advice. I have a set of
> multivariate data and would like to find the densest part of the data cloud
> containing 95% of the data, like a 95% HPD credible volume. Is there any R
> code available to compute that?

It looks like the LaplacesDemon pkg was just updated FWIW,

http://www.google.com/search?hl=en&q=hpd+credible+volume+cran

If youjust  want to find the density of your data in some n-dim space, that 
sounds like
multi dim binning or histogram would work (you can check google but IIRC there 
is
no general n-dim binning pacakge.  I also mentioned a version based
on a density field being associated with each point which could be summed over
all points to get density at arbirtary point but you need to implrement that in 
some
intelligent way for it to be fast ). I guess you could bin using aggregate see 
if this
example would work,

df<-data.frame(z=runif(100), x<-runif(100)*10,y<-runif(100)*5,c<-rep(1,100))
d<-aggregate(c~floor(x)+floor(y),df,FUN="sum")
d
which (d$c==max(d$c))

now at this point you can imagine you may be able to find a surface that 
encloses 
a certain amount of your data



> 
> Thank you very much! Your help and patience are much appreciated.
> 
> G.S.
> 
>   [[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] Colors in a dendrogram

[Ayoub Maatallaoui]

Hello,
i'm trying to draw a dendrogram to show the clustering hierarchy. i'm 
working on 2 classes but each class contains more than 3000 instances.
i'm using the hclust function to do that! but i can't read anything from 
this figure.
i would like to color each class with its own color so at least i can 
extract some information from the clustering.

does any one of you did something like this before?

Regards,
Ayoub.

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

Re: [R] how to get the miminum value in the list

[Dennis Murphy]

Hi:

It's just an extra step:

y <- list(list(c(1,5),c(2,3,4)), list(c(1, 3, 4), c(5, 7)))
lapply(seq_len(length(y)), function(i) lapply(y[[i]], min))
[[1]]
[[1]][[1]]
[1] 1

[[1]][[2]]
[1] 2

[[2]]
[[2]][[1]]
[1] 1

[[2]][[2]]
[1] 5

unlist(lapply(seq_len(length(y)), function(i) lapply(y[[i]], min)))
[1] 1 2 1 5

HTH,
Dennis

On Sun, Jun 19, 2011 at 8:25 AM, jiliguala  wrote:
>
> but in my case, the list is a two-variable list, list[[j]][[i]]
>
> when i use
>
> lapply(list, min)
>
> it appears
> """Error in FUN(X[[1L]], ...) : invalid 'type' (list) of argument"""
>
> thanks
>
> --
> View this message in context: 
> http://r.789695.n4.nabble.com/how-to-get-the-miminum-value-in-the-list-tp3609013p3609433.html
> Sent from the R help mailing list archive at Nabble.com.
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

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

[R] Accessor functions in lattice graphics

[Frank Harrell]

I know about the current.row, current.column, and panel.number functions that
are useful within panel functions written for lattice.  Are there easy ways
to obtain the names of the conditioning variables (those appearing after |)
and their values for the current panel?
Thanks
Frank

-
Frank Harrell
Department of Biostatistics, Vanderbilt University
--
View this message in context: 
http://r.789695.n4.nabble.com/Accessor-functions-in-lattice-graphics-tp3609435p3609435.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] how to get the miminum value in the list

[jiliguala]

but in my case, the list is a two-variable list, list[[j]][[i]]

when i use 

lapply(list, min)

it appears  
"""Error in FUN(X[[1L]], ...) : invalid 'type' (list) of argument"""

thanks

--
View this message in context: 
http://r.789695.n4.nabble.com/how-to-get-the-miminum-value-in-the-list-tp3609013p3609433.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] Using MCMC sampling to estimate p values with a mixed model

[Robert A LaBudde]

If your alternative hypothesis is unequal variances (2-sided), both F 
< 1 and F > 1 are of interest, and rejection of the equal variance 
null can occur on either side.


The usual ANOVA F test is 1-sided, with an alternative the numerator 
variance exceeds the denominator one, so this is perhaps why you are confused.


At 02:40 PM 6/18/2011, Woodcock, Helena wrote:

Hi everyone,

Apologies if this is a silly question but I am a student and this is 
my first time using R so I am still trying to educate myself on 
commands, models e.t.c


I have a mixed model with four dichotomous fixed factors and subject 
as a random factor (as each person completed four vignettes, with 
factors crossed across vignettes).


I have run an lmer model and used the Monte Carlo method to see if 
there are any significant main effects or interactions. However, 
when I looked at the p values some are showing as significant 
although the F value is less than 1. Is it possible to have a 
significant effect with an F value below 1?.


I have a sample size of 150 and have read that the pMCMC values can 
be anti-conservative so wonder if it is because my sample size may 
be too small?.


Thank you for any help

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



Robert A. LaBudde, PhD, PAS, Dpl. ACAFS  e-mail: r...@lcfltd.com
Least Cost Formulations, Ltd.URL: http://lcfltd.com/
824 Timberlake Drive Tel: 757-467-0954
Virginia Beach, VA 23464-3239Fax: 757-467-2947

"Vere scire est per causas scire"

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

Re: [R] save and load in R

[Mary Kindall]

Thanks  Jeff and Duncan
Assign worked for me. I did not check the other methods suggested by you.  I
am repreducing my code below:

##
filenames = list.files(path = ".", pattern = '.txt$', all.files = FALSE,
full.names = TRUE, ignore.case = FALSE)  #all input files
numFiles = length(filenames)

outfilenames <- paste("./file", 1:numFiles, '.Rdata', sep="") # output files


for ( i in 1:numFiles)
{
dataFile = read.table(filenames[i], header=TRUE, sep='\t');
save(dataFile, file = outfilenames[i])   #Saving into the output files
}


newnames <- paste("file", 1:numFiles, sep="")  #output variables to load
into

for ( i in 1:numFiles)
{
load(file = outfilenames[i]);
assign(newnames[i], dataFile)   #assign into corresponding output variables;

}
#


Regards
-
M



On Sun, Jun 19, 2011 at 10:38 AM, Duncan Murdoch
wrote:

> On 11-06-19 10:26 AM, Mary Kindall wrote:
>
>> I have a list of txt files that I want to convert into .rdata R data
>> object.
>>
>> filenames
>> 1. "./file1.txt"
>> 2. "./file2.txt"
>> 3. "./file3.txt"
>> 4. "./file4.txt"
>> 5. "./file5.txt"
>> 6. "./file6.txt"
>> 7. "./file7.txt"
>> 8. "./file8.txt"
>> 9. "./file9.txt"
>> 10. "./file10.txt"
>>
>> I saved these files as
>>
>> for ( i in 1:10)
>> {
>> dataFile = read.table(filenames[i], header=TRUE, sep='\t');
>> save (dataFile, file = outfilenames[i])
>> }
>>
>> The inpt files are saves as:
>> outfilenames
>> 1. "./file1.Rdata"
>> 2. "./file2.Rdata"
>> 3. "./file3.Rdata"
>> 4. "./file4.Rdata"
>> 5. "./file5.Rdata"
>> 6. "./file6.Rdata"
>> 7. "./file7.Rdata"
>> 8. "./file8.Rdata"
>> 9. "./file9.Rdata"
>> 10. "./file10.Rdata"
>>
>>
>> Now I want to load these out files in such a way that the data is loaded
>> into a variable that is same as the file name without extension.
>>
>> file1 = load (file = './file1.Rdata')
>> file2 = load (file = './file2.Rdata')
>> file3 = load (file = './file3.Rdata')
>> file4 = load (file = './file4.Rdata')
>>
>> How can I do that.
>>
>
> When you load() a file, the variables in it are restored with the same
> names that were saved.  So you would need something like
>
> newnames <- paste("file", 1:10, sep="") # file1, file2, etc.
>
>
> for (i in 1:10) {
>  load(file=outfilenames[i]) # assuming that's still around...
>  assign(newnames[i], dataFile)
> }
>
> It would be a little simpler to use saveRDS() and readRDS() to save and
> load your files.  They don't save the object names.
>
> A more R-like version of this would be to create a list of datasets, e.g.
>
> files <- list()
>
> for (i in 1:10) {
>  load(file=outfilesnames[i])
>  files[[i]] <- dataFile
> }
>
> Then you don't end up creating 10 objects, but you can still access them
> separately.
>
> Duncan Murdoch
>
>


-- 
-
Mary Kindall
Yorktown Heights, NY
USA

[[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] save and load in R

[Duncan Murdoch]

On 11-06-19 10:26 AM, Mary Kindall wrote:

I have a list of txt files that I want to convert into .rdata R data
object.

filenames
1. "./file1.txt"
2. "./file2.txt"
3. "./file3.txt"
4. "./file4.txt"
5. "./file5.txt"
6. "./file6.txt"
7. "./file7.txt"
8. "./file8.txt"
9. "./file9.txt"
10. "./file10.txt"

I saved these files as

for ( i in 1:10)
{
dataFile = read.table(filenames[i], header=TRUE, sep='\t');
save (dataFile, file = outfilenames[i])
}

The inpt files are saves as:
outfilenames
1. "./file1.Rdata"
2. "./file2.Rdata"
3. "./file3.Rdata"
4. "./file4.Rdata"
5. "./file5.Rdata"
6. "./file6.Rdata"
7. "./file7.Rdata"
8. "./file8.Rdata"
9. "./file9.Rdata"
10. "./file10.Rdata"


Now I want to load these out files in such a way that the data is loaded
into a variable that is same as the file name without extension.

file1 = load (file = './file1.Rdata')
file2 = load (file = './file2.Rdata')
file3 = load (file = './file3.Rdata')
file4 = load (file = './file4.Rdata')

How can I do that.


When you load() a file, the variables in it are restored with the same 
names that were saved.  So you would need something like


newnames <- paste("file", 1:10, sep="") # file1, file2, etc.

for (i in 1:10) {
  load(file=outfilenames[i]) # assuming that's still around...
  assign(newnames[i], dataFile)
}

It would be a little simpler to use saveRDS() and readRDS() to save and 
load your files.  They don't save the object names.


A more R-like version of this would be to create a list of datasets, e.g.

files <- list()
for (i in 1:10) {
  load(file=outfilesnames[i])
  files[[i]] <- dataFile
}

Then you don't end up creating 10 objects, but you can still access them 
separately.


Duncan Murdoch

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

[R] save and load in R

[Mary Kindall]

I have a list of txt files that I want to convert into .rdata R data
object.

filenames
1. "./file1.txt"
2. "./file2.txt"
3. "./file3.txt"
4. "./file4.txt"
5. "./file5.txt"
6. "./file6.txt"
7. "./file7.txt"
8. "./file8.txt"
9. "./file9.txt"
10. "./file10.txt"

I saved these files as

for ( i in 1:10)
{
dataFile = read.table(filenames[i], header=TRUE, sep='\t');
save (dataFile, file = outfilenames[i])
}

The inpt files are saves as:
outfilenames
1. "./file1.Rdata"
2. "./file2.Rdata"
3. "./file3.Rdata"
4. "./file4.Rdata"
5. "./file5.Rdata"
6. "./file6.Rdata"
7. "./file7.Rdata"
8. "./file8.Rdata"
9. "./file9.Rdata"
10. "./file10.Rdata"


Now I want to load these out files in such a way that the data is loaded
into a variable that is same as the file name without extension.

file1 = load (file = './file1.Rdata')
file2 = load (file = './file2.Rdata')
file3 = load (file = './file3.Rdata')
file4 = load (file = './file4.Rdata')

How can I do that.
Regards
-
M

-- 
-
Mary Kindall
Yorktown Heights, NY
USA

[[alternative HTML version deleted]]

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

Re: [R] how to get the miminum value in the list

[Daniel Malter]

Hi, that depends on whether you want to get the minimum within each list
element or the global minimum across all list elements. The first is
achieved by using lapply(). The second can be achieved by unlisting the list
(which assumes that all list elements are numeric) and looking for its
minimum.

x<-list(c(1,5),c(2,3,4))
lapply(x,min)
min(unlist(x))

HTH,
Daniel


jiliguala wrote:
> 
> hi, R users 
> 
> here i have one problem, if i wanna get the minimum value in the normal
> data, i can do this,
> 
> ## which(data1==min(data1)).
> 
> but if i want get the minimum value of a list which has two variables
> ##list1[[j]][[i]]##,
> i tried the codes like this, but it did not work.
> 
> ##  which(list1==min(list1)).
> 
> thanks for helping
> 

--
View this message in context: 
http://r.789695.n4.nabble.com/how-to-get-the-miminum-value-in-the-list-tp3609013p3609325.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] automatically generate the output name of my "for" loops

[jiliguala]

true, it did work using unlist() funciton, it can be found the minimum
value.

and for further step, if i want identify the minimum value with the list, 

which means the minimum value i have found is belonged to which
list[[j]][[i]].

i used the codes like this

##which(list1==min(list1))##

but it did not work.


thanks for helping


--
View this message in context: 
http://r.789695.n4.nabble.com/automatically-generate-the-output-name-of-my-for-loops-tp3592160p3609298.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] Using MCMC sampling to estimate p values with a mixed model

[Ben Bolker]

Woodcock, Helena  liverpool.ac.uk> writes:

> Apologies if this is a silly question but I am a student 
> and this is my first time using R so I am still trying to
> educate myself on commands, models e.t.c
> 
> I have a mixed model with four dichotomous fixed factors and 
> subject as a random factor (as each person
> completed four vignettes, with factors crossed across vignettes).
> 
> I have run an lmer model and used the Monte Carlo method to 
> see if there are any significant main effects or
> interactions. However, when I looked at the p values 
> some are showing as significant although the F value
> is less than 1. Is it possible to have a significant 
> effect with an F value below 1?.
> 
> I have a sample size of 150 and have read that the 
> pMCMC values can be anti-conservative so wonder if it is
> because my sample size may be too small?.
> 

  It's hard to know exactly without more details/seeing the data;
it does sound suspicious.
  Unless there's something you haven't told us, it sounds like
this model is fairly close to a classical experimental design
(randomized block), so I would guess that the answers should (?) be
fairly close to the classical ones.  Have you tried fitting with
lme (from the nlme package) and seeing what it guesses for denominator
degrees of freedom, or working out appropriate denominator degrees
of freedom yourself?
  
  I would recommend that you send follow-ups to 
r-sig-mixed-mod...@r-project.org ...

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

[R] apply kmeans into every level of the tree..

[jiliguala]

hi, r users 

I'm a new r user. i want to apply kmeans into every level of the tree which
i built.

the codes for clusterring only one level of the tree.

#
d  <- dist(data1, method = "euclidean")

h1 <- hclust(d, method="ward")  

c1 <- cutree(h1, 100)  # 100 groups which i want to cluster

cdg <- aggregate(data1, list(c1),mean)[,2:(nd+1)] ## nd is the
attributes i have
3
k1 <- kmeans(daata1, centers=cdg)
##

now i want to cluster every level of the tree in order to search it.

any thought about that, thanks for helping.

--
View this message in context: 
http://r.789695.n4.nabble.com/apply-kmeans-into-every-level-of-the-tree-tp3609250p3609250.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] please help! what are the different using log-link function and log transformation?

[Rubén Roa]

Funny, I was going to quote your paper with Dichmont in Fisheries Research 
regarding this, then I forgot to do it, but you yourself came out and posted 
the explanation.

This forum is great!

Rubén

-Original Message-
From: r-help-boun...@r-project.org on behalf of bill.venab...@csiro.au
Sent: Sun 6/19/2011 2:36 PM
To: gloryk...@hotmail.com; r-help@r-project.org
Subject: Re: [R] please help! what are the different using log-link function 
and log transformation?
 
The two commands you give below are certain to lead to very different results, 
because they are fitting very different models.

The first is a gaussian model for the response with a log link, and constant 
variance.

The second is a gaussian model for a log-transformed response and identity 
link.  On the original scale this model would imply a constant coefficient of 
variation and hence a variance proportional to the square of the mean, and not 
constant.

Your problem is not particularly an R issue, but a difficulty with 
understanding generalized linear models (and hence generalized additive models, 
which are based on them).

Bill Venables.

From: r-help-boun...@r-project.org [r-help-boun...@r-project.org] On Behalf Of 
pigpigmeow [gloryk...@hotmail.com]
Sent: 19 June 2011 17:39
To: r-help@r-project.org
Subject: [R] please help! what are the different using log-link function and 
log transformation?

I'm new R-programming user, I need to use gam function.

y <- gam(a ~ s(b), family = gaussian(link=log),  data)
y <- gam(log(a) ~ s(b), family = gaussian (link=identity), data)

why [do] these two command [give different] results?

I guess these two command results are same, but actally these two command
results are different, Why?

--
View this message in context: 
http://r.789695.n4.nabble.com/please-help-what-are-the-different-using-log-link-function-and-log-transformation-tp3608931p3608931.html
Sent from the R help mailing list archive at Nabble.com.

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

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


[[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] please help! what are the different using log-link function and log transformation?

[Bill.Venables]

The two commands you give below are certain to lead to very different results, 
because they are fitting very different models.

The first is a gaussian model for the response with a log link, and constant 
variance.

The second is a gaussian model for a log-transformed response and identity 
link.  On the original scale this model would imply a constant coefficient of 
variation and hence a variance proportional to the square of the mean, and not 
constant.

Your problem is not particularly an R issue, but a difficulty with 
understanding generalized linear models (and hence generalized additive models, 
which are based on them).

Bill Venables.

From: r-help-boun...@r-project.org [r-help-boun...@r-project.org] On Behalf Of 
pigpigmeow [gloryk...@hotmail.com]
Sent: 19 June 2011 17:39
To: r-help@r-project.org
Subject: [R] please help! what are the different using log-link function and 
log transformation?

I'm new R-programming user, I need to use gam function.

y <- gam(a ~ s(b), family = gaussian(link=log),  data)
y <- gam(log(a) ~ s(b), family = gaussian (link=identity), data)

why [do] these two command [give different] results?

I guess these two command results are same, but actally these two command
results are different, Why?

--
View this message in context: 
http://r.789695.n4.nabble.com/please-help-what-are-the-different-using-log-link-function-and-log-transformation-tp3608931p3608931.html
Sent from the R help mailing list archive at Nabble.com.

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

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

Re: [R] Required libraries

[Duncan Murdoch]

On 11-06-19 3:19 AM, Jim Lemon wrote:

On 06/18/2011 08:47 PM, Trying To learn again wrote:

Hi all,

I´m trying to "recuperate" old files I wrotte and I´m trying to execute on R
version 2.13.0 (2011-04-13), the thing is I execute my old file and nothing
happens...I suposse I need to install some library but there no appears no
message, no text telling "required grafics library" or something like that.
How can I see wich libraries need to execute a file?


Hi Trying,
This is way out in guess field, but you may be running R scripts in a
way that the error messages are not displayed on your screen. If so, try
opening the script in a text editor (e.g. Notepad if you are using
Windows) and cut and paste the commands into the console window one by
one. When something goes wrong, you should get the error messages that
you want.


One minor comment:

Don't use Notepad, use the built-in script editor if you are using the 
GUI in Windows.  It's not a great editor, but it's better than Notepad. 
 You can use Ctrl-R to cut and paste one line at a time.


Duncan Murdoch




Do you recommend to work with the last avaliable CRAN R version?


Probably a good idea.

Jim

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


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

Re: [R] please help! what are the different using log-link function and log transformation?

[Rubén Roa]

The problem is not that you are new to R. This is a conceptual issue.

Let y be the response variable and let {x_i} be a set of predictors. Your first 
model (identity response and log-link) is saying that 

y=f(x_1)f(x_2)...f(x_n) + e, e~Normal(0,sigma)

i.e. this is an additive observation-error model with constant variance.

Your second model (log-response and identity link) is saying that
 
y=f(x_1)f(x_2)...f(x_n)u, u=exp(e), e~Normal(0,sigma)

i.e. this a multiplicative observation-error model with variance proportional 
to the mean.

Plot the data versus response and visually examine whether you have 
heteroscedasticity. If this is true, use your second model, otherwise use the 
first one.

One key to understand these kind of dichotomies is to realize that statistical 
models are composed of a process part and an observation part. In your models 
the process part is deterministic and multiplicative but after that, you still 
have two choices, make the random observation part additive (your first model) 
or multiplicative (your second model).

Needless to say (but I am saying it anyways) these two models will give 
different results, at the very least because one assumes constant variance 
(your first model) whereas the other assumes a variance proportional to the 
mean.

In my experience with multiplicative process models, the random observation 
part shall usually be multiplicative as well because of heteroscedasticity.

HTH

Rubén

-

Rubén H. Roa-Ureta, Ph. D.
AZTI Tecnalia, Txatxarramendi Ugartea z/g,
Sukarrieta, Bizkaia, SPAIN

-Original Message-
From: r-help-boun...@r-project.org on behalf of pigpigmeow
Sent: Sun 6/19/2011 9:39 AM
To: r-help@r-project.org
Subject: [R] please help! what are the different using log-link function and 
log transformation?

I'm new R-programming user, I need to use gam function.

y<-gam(a~s(b),family=gaussian(link=log),data)
y<-gam(loga~s(b), family =gaussian (link=identity),data)
why these two command results are different?
I guess these two command results are same, but actally these two command
results are different, Why?

--
View this message in context: 
http://r.789695.n4.nabble.com/please-help-what-are-the-different-using-log-link-function-and-log-transformation-tp3608931p3608931.html
Sent from the R help mailing list archive at Nabble.com.

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


[[alternative HTML version deleted]]

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

Re: [R] [rgl] Dynamically change the color of 3d objects

[Duncan Murdoch]

On 11-06-18 3:47 PM, Fernando (Est - UnB) wrote:

Hi,

I'm new here.

This question concerns the package rgl.

Is it possible to change the colour of a set of points drawn with
plot3d /without/ removing them from the scene?


No, you can't modify objects other than by removing and recreating them.



The idea is to create a presentation of a clustering algorithm, step
by step, representing the currently assigned cluster by a colour.
Removing all the points with "rgl.pop" and replotting them at each
step works, but is kinda slow...


If you do it this way:

# plot the points, and save the id in pointsID

par3d(skipRedraw=TRUE)
rgl.pop(pointsID)
pointsID <- points3d(  ...  ) # plot the points again
par3d(skipRedraw=FALSE)

it will go faster, and appear to the viewer as though the colour just 
changed, with no removal of points.


Duncan Murdoch

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

[R] how to get the miminum value in the list???

[jiliguala]

hi, R users 

here i have one problem, if i wanna get the minimum value in the normal
data, i can do this,

## which(data1==min(data1)).

but if i want get the minimum value of a list which has two variables
##list1[[j]][[i]]##,
i tried the codes like this, but it did not work.

##  which(list1==min(list1)).

hope anyone to solve that.

thanks

--
View this message in context: 
http://r.789695.n4.nabble.com/how-to-get-the-miminum-value-in-the-list-tp3609013p3609013.html
Sent from the R help mailing list archive at Nabble.com.

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

Re: [R] Applying function to all elements of a formula

[Scott Fortmann-Roe]

Thanks for the tips. I'll give them a try!

On Sat, Jun 18, 2011 at 12:07 PM, Dennis Murphy  wrote:

> Much better..nice!
>
> Dennis
>
> On Sat, Jun 18, 2011 at 1:53 AM, Dimitris Rizopoulos
>  wrote:
> > maybe another way is by reconstructing the formula using paste(), e.g.,
> >
> > data <- data.frame(y = rnorm(5), x1 = runif(5),
> >z = runif(5), age = runif(5))
> >
> > nameRsp <- "y"
> > nams <- names(data)
> > namsX <- nams[!nams %in% nameRsp]
> > form <- as.formula(paste(nameRsp, "~" ,
> >paste("log(", namsX, ")", sep = "", collapse = "+")))
> >
> > lm(form, data)
> >
> >
> > I hope it helps.
> >
> > Best,
> > Dimitris
> >
> >
> > On 6/18/2011 10:41 AM, Dennis Murphy wrote:
> >>
> >> Yes, it's possible, but if you want to do prediction on future
> >> x-values, you will likely have a problem.
> >>
> >> One way to do it would be something like (assuming y is the first column
> >> of dat)
> >>
> >> reg<- lm(y ~ log(as.matrix(dat[, -1])), dat)
> >>
> >> but the output would be pretty ugly (see summary(reg)). Another would
> >> be to construct the matrix outside the data frame and do something
> >> like
> >>
> >> X<- log(as.matrix(dat[, -1]))
> >> reg<- lm(dat$y ~ X)
> >>
> >> There may be better ways, though...
> >>
> >> Dennis
> >>
> >> On Fri, Jun 17, 2011 at 11:08 PM, Scott Fortmann-Roe
> >>  wrote:
> >>>
> >>> Hi,
> >>>
> >>> I would like to do a regression like:
> >>>
> >>> reg<- lm(y~log(.), data)
> >>>
> >>> where the log function is applied to "." in the form:
> >>>
> >>>log(x1)+ log(x2)+ log(x3)...
> >>>
> >>> instead of in the form
> >>>
> >>>   log(x1+x2+x3+...)
> >>>
> >>>
> >>> Is this possible?
> >>>
> >>> Thank you,
> >>> Scott
> >>>
> >>>[[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.
> >>
> >
> > --
> > Dimitris Rizopoulos
> > Assistant Professor
> > Department of Biostatistics
> > Erasmus University Medical Center
> >
> > Address: PO Box 2040, 3000 CA Rotterdam, the Netherlands
> > Tel: +31/(0)10/7043478
> > Fax: +31/(0)10/7043014
> > Web: http://www.erasmusmc.nl/biostatistiek/
> >
>

[[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] Asymetrical Confidence Interval

[Tal Galili]

Hi Muhuri,
What do you suspect the underlying distribution is of your statistic?



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




On Thu, Jun 16, 2011 at 11:43 PM, Muhuri, Pradip (SAMHSA/CBHSQ) <
pradip.muh...@samhsa.hhs.gov> wrote:

>
> Dear List,
>
> I wanted to calculate the asymmetrical confidence interval based on the
> sample statistic and standard error that available from the published report
> (complex survey-based).
>
> The calculation details can be seen from pages 17-18 of the document at the
> following link:
> http://www.oas.samhsa.gov/nsduh/2k5MRB/2k5statInference.pdf.
>
>
> Could someone tell me whether R has any function included in it "survey" or
> other contributed package of R.
>
> Thank you in advance,
>
> Pradip Muhuri
>
> __
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/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] please help! what are the different using log-link function and log transformation?

[pigpigmeow]

I'm new R-programming user, I need to use gam function.

y<-gam(a~s(b),family=gaussian(link=log),data)
y<-gam(loga~s(b), family =gaussian (link=identity),data)
why these two command results are different?
I guess these two command results are same, but actally these two command
results are different, Why?

--
View this message in context: 
http://r.789695.n4.nabble.com/please-help-what-are-the-different-using-log-link-function-and-log-transformation-tp3608931p3608931.html
Sent from the R help mailing list archive at Nabble.com.

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

[R] Multivariate HPD credible volume -- is it computable in R?

[agw1118]

Hi all,

I'm new to the list and am hoping to get some advice. I have a set of
multivariate data and would like to find the densest part of the data cloud
containing 95% of the data, like a 95% HPD credible volume. Is there any R
code available to compute that?

Thank you very much! Your help and patience are much appreciated.

G.S.

[[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] Required libraries

[Jim Lemon]

On 06/18/2011 08:47 PM, Trying To learn again wrote:

Hi all,

I´m trying to "recuperate" old files I wrotte and I´m trying to execute on R
version 2.13.0 (2011-04-13), the thing is I execute my old file and nothing
happens...I suposse I need to install some library but there no appears no
message, no text telling "required grafics library" or something like that.
How can I see wich libraries need to execute a file?


Hi Trying,
This is way out in guess field, but you may be running R scripts in a 
way that the error messages are not displayed on your screen. If so, try 
opening the script in a text editor (e.g. Notepad if you are using 
Windows) and cut and paste the commands into the console window one by 
one. When something goes wrong, you should get the error messages that 
you want.



Do you recommend to work with the last avaliable CRAN R version?


Probably a good idea.

Jim

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

[R] [rgl] Dynamically change the color of 3d objects

[Fernando (Est - UnB)]

Hi,

I'm new here.

This question concerns the package rgl.

Is it possible to change the colour of a set of points drawn with
plot3d /without/ removing them from the scene?

The idea is to create a presentation of a clustering algorithm, step
by step, representing the currently assigned cluster by a colour.
Removing all the points with "rgl.pop" and replotting them at each
step works, but is kinda slow...

Thanks!
Fernando.

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

Re: [R] prediction intervals

[Dave Evens]

Thank you for your email.

The problem is that i'm not quite sure how to specify the model using arima().

Here's an example of my problem:

I would like to fit an ARMA model, but I'm not sure exactly how to fit it. 
  

Here's an example of the problem. 
  
This is my time variable, hourly data  
t <- seq(as.POSIXct("2011-01-01 00:00:00"), as.POSIXct("2011-12-31 23:00:00"), 
by="hour") 
  
My response 
y <- rnorm(length(t), 1000, 500) 
  
Seasonal factors: 
t.h <- as.POSIXlt(t)$hour # hours of the day
t.d <- as.POSIXlt(t)$day  # days of the week
t.m <- as.POSIXlt(t)$mon  # months of the year
This is my regressor 
x.reg <- rnorm(length(t), 10, 1) 
  
  
and I have the following auto-regressive lags (1 to 10, 24 and 48 hours) 
y.lag1 <- lag(y, 1)
y.lag2 <- lag(y, 2)
y.lag3 <- lag(y, 3)
y.lag4 <- lag(y, 4)
y.lag5 <- lag(y, 5)
y.lag6 <- lag(y, 6)
y.lag7 <- lag(y, 7)
y.lag8 <- lag(y, 8)
y.lag9 <- lag(y, 9)
y.lag10 <- lag(y, 10)
y.lag24 <- lag(y, 24)
y.lag48 <- lag(y, 48)



I want to fit an ARMA (auto-regressive moving average) with my 3 seasonal 
factors, 12 lagged variables and the regressor against my response variable.
Does someone know how such an ARMA model can be fit? 


Regards,
Dave





From: Greg Snow 

roject.org>
Sent: Friday, 17 June 2011, 21:40
Subject: RE: [R] prediction intervals


I am not an expert in time series (that is why I referred you to the task view 
rather than give my own inexpert opinion).  I do remember from a textbook that 
covered the basics of time series that prediction 2 time points ahead was 
different from plugging in the next time estimate and predicting one ahead.  
And basic theory of linear models also supports this.  
 
The standard linear model assumes all x’s to be fixed (or independent of each 
other if random) if you add a new x that is not independent and not fixed, then 
this does not hold.  There is also a whole area of linear models for what to 
do if x is not an exact value but measured with error (or predicted), another 
area that I know exists and is complicated, but which I have not studied beyond 
the basics.  
 
Also the standard formula for prediction intervals has a piece to account for 
variability around the mean and another piece to account for the variability in 
the estimate of the mean due to the coefficients being estimates rather than 
known quantities, it only makes sense that there should be another piece to 
account for uncertainty due to the value of x in the prediction when it is not 
known.
 
The best I can recommend is look through the task view.  Maybe someone else 
has a better refrence.
 
-- 
Gregory (Greg) L. Snow Ph.D.
Statistical Data Center
Intermountain Healthcare
greg.s...@imail.org
801.408.8111
 

Sent: Friday, June 17, 2011 1:48 AM
To: Greg Snow; r-help@r-project.org
Subject: Re: [R] prediction intervals
 
Thank you for your post Greg.
Do you have any useful references regarding this variability (papers etc)?
 
Many thanks.
Dave
 
 
From:Greg Snow 

roject.org>
Sent: Thursday, 16 June 2011, 21:32
Subject: RE: [R] prediction intervals

I don't think that this approach is appropriate here.  Each iteration after 
the 1st the lm/predict combination will assume that the new data is exact when 
in fact it is an estimate with some error involved.  To properly do this you 
need to take into account that variability.  There is a time series task view 
on CRAN that may point you to better tools.

-- 
Gregory (Greg) L. Snow Ph.D.
Statistical Data Center
Intermountain Healthcare
greg.s...@imail.org
801.408.8111


> -Original Message-
> From: r-help-boun...@r-project.org [mailto:r-help-bounces@r-
> project.org] On Behalf Of Dave Evens
> Sent: Thursday, June 16, 2011 11:33 AM
> To: r-help@r-project.org
> Subject: [R] prediction intervals
> 
> 
> 
> Dear members,
> 
> I'm fitting linear model using "lm" which has numerous auto-regressive
> terms as well as other explanatory variables. In order to calculate
> prediction intervals, i've used a for-loop as the auto-regressive
> parameters need to be updated each time so that a new forecast and
> corresponding prediction interval can be calculated.
> 
> I'm fitting a number of these models which have different values for
> the response variable and possibly different explanatory variables. The
> response is temperature in fahrenheit (F), and the different models are
> for cities. So each city has its own fitted linear model for
> temperature. I'm assuming that they're independent models for the time
> being, I want to combine the results across all cities and have overall
> prediction intervals. Because I assuming that they're independent can I
> just add together the degrees of freedom from each model (i.e. total
> degrees of freedom=df1+df2+...) and the variance-covariance matrices
> (i.e. V=V1+V2+...) in order to calcalate the overall prediction
> intervals?
> 
> Any help would be most appreciated.
> 
> Regards,
> Dave
>

[R] help with quantification of high density region

[Christine Buske]

Hello,

I am entirely new to this list, and not very experienced with R yet (but very 
excited about it!!). I am struggling with a dataset of x and y coordinate 
points I am plotting as a density map as follows:

trial <- read.csv(file="F12.csv", sep=",", skip=32, head=TRUE)
trial.sub <- subset(trial, Seconds>=30)
trial.sub <- subset(trial.sub, Seconds<=1830)

d <- trial.sub
require(MASS)
dens <- kde2d(d$x, d$y, h=20, n=100)
color <- function(x)rev(heat.colors(x))
filled.contour(dens, color.palette = color)

Each file contains 53000 observations (29 per second), so 53000 x and y 
coordinates. I have posted a sample of my heat map here: 
https://picasaweb.google.com/neurosci.christine/Public?authkey=Gv1sRgCLO3m6Obpsr_sgE#5619604653154721794

My main question I can't seem to find an answer to is what the numbers in the 
legend mean. The heatmap is created for coordinate points of a location 
occupied by an animal, so the high density means the animal spent most of its 
time there, but what percentage of time? if is possible to get this out of the 
data and/or displayed more clearly in the legend?

I addition, I would like to quantify the area size of the highest density 
region. Alternatively or in addition to this, I wouldn't mind measuring how 
many 'zones' the highest density region occupies. This might be more difficult, 
unless I can overlay a grid overtop of the data and count the number of zones 
covered. I have tried different strategies for this, including creating subsets 
of the coordinate point ranges, but I couldn't get it to work and any input 
will be greatly appreciated. To see a sample of another heatmap where this 
would be very relevant is posted here:

https://picasaweb.google.com/neurosci.christine/Public?authkey=Gv1sRgCLO3m6Obpsr_sgE#5619604756338447826


I really appreciate any help I could get with this! I've been working on it for 
a week and not getting any further ahead. 
Christine


Below is some of the data posted in case it might be helpful:

Trackingarea(x1 


125 7   547 482 


Zone1-1 (x1 


127 3   273 130 


Zone2-2 (x1 


272 5   412 131 


Zone3-3 (x1 


413 7   556 128 


Zone4-4 (x1 


128 128 274 239 


Zone5-5 (x1 


273 128 415 245 


Zone6-6 (x1 


417 129 555 245 


Zone7-7 (x1 


127 239 272 348 


Zone8-8 (x1 
 

Re: [R] Bartlett's Test of Sphericity

[Tal Galili]

Sorry for the confusion, thank you for the correction and explanation.


Tal

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




On Sat, Jun 18, 2011 at 7:38 PM, peter dalgaard  wrote:

>
> On Jun 18, 2011, at 15:48 , Tal Galili wrote:
>
> > Hi Thibault,
> > Not that I think you'll use this after the above responses, but for the
> > record, have a look at:
> > ?bartlett.test
>
> That's the "other" Bartlett's test, namely the one for comparison of
> variances. That test is well known to rely on higher moments of the normal
> distribution.
>
> --
> Peter Dalgaard
> Center for Statistics, Copenhagen Business School
> Solbjerg Plads 3, 2000 Frederiksberg, Denmark
> Phone: (+45)38153501
> Email: pd@cbs.dk  Priv: pda...@gmail.com
>
>

[[alternative HTML version deleted]]

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

[R] Using MCMC sampling to estimate p values with a mixed model

[Woodcock, Helena]

Hi everyone,

Apologies if this is a silly question but I am a student and this is my first 
time using R so I am still trying to educate myself on commands, models e.t.c

I have a mixed model with four dichotomous fixed factors and subject as a 
random factor (as each person completed four vignettes, with factors crossed 
across vignettes).

I have run an lmer model and used the Monte Carlo method to see if there are 
any significant main effects or interactions. However, when I looked at the p 
values some are showing as significant although the F value is less than 1. Is 
it possible to have a significant effect with an F value below 1?.

I have a sample size of 150 and have read that the pMCMC values can be 
anti-conservative so wonder if it is because my sample size may be too small?.

Thank you for any help

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