[R] Bayesian PLS-Based Dimension Reduction and Classification Methods for Time to Event Data in R
Dear all, Please I am currently working on Bayesian Partial Least Squares Dimension Reduction and Classification methods for Modelling survival data in R. Can someone advise on useful R packages and good R books to read. I am also interested in working with experienced colleagues who may be interested. Thank you hugely for your anticipated assistance. Kabir Olorede, Department of Statistics and Mathematical Sciences, College of Pure and Applied Sciences, P.M.B. 1530, Kwara State University, Malete, Nigeria Alt. E-mail: kabir.olor...@kwasu.edu.ng [[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] Construction of Optimal Designs in the Language R
Dear All, Please can someone help advise on procedures for construction of optimal designs and estimation of parameters in the language R? Guide on useful R libraries, codes and books would be appreciated. Kabir Olorede, Department of Statistics and Mathematical Sciences, Kwara State University, Malete, Nigeria. [[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] matafor package: adding additional information under a forest plot
1) There is no function that will do that automatically for you, but it's easy to do with text() and paste(). If you have fitted a random-effects model with rma(), then you can find the Q-statistic in res$QE, the I^2 statistic in res$I2, and probably whatever else you need, and then just paste() together what you want and add it to the plot via text(). 2) If you want to add lots of columns, then indeed you will need a wide plot. You probably want to set up your plotting device to be quite wide. For example, with png(), make sure 'width' is large and you may need to decrease the 'pointsize' to still make it all fit. Best, Wolfgang From: r-help-boun...@r-project.org [r-help-boun...@r-project.org] On Behalf Of Dietlinde Schmidt [schmidt.dietli...@web.de] Sent: Friday, May 02, 2014 4:09 PM To: r-help@r-project.org Subject: [R] matafor package: adding additional information under a forest plot Dear R-Users, i have to questions about plotting a forest plot with the metafor package: 1) Is there a (text()) function to add additional information (e.g. heterogeneity statistics) under the forest plot? 2) I want to add various columns of additional information about each study (e.g. quality, time point, measure). I thought i could do that with the ilab-argument, but i keep getting in the way with the line that restricts the plot and although i tweaked everything around, I do not succeed. Basecally, i want the format of the forest plot to be rather wide than tall, because i have little studies but many columns with additional study information. Does someone know how to solve that (over xlim, ylim???)? Thanks a lot for your help Linde __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Return TRUE only for first match of values between matrix and vector.
Hi,
This should be little more faster.
indx - A==B
indx1 - which(indx, arr.ind=TRUE)
indx[indx1[duplicated(indx1[,1]),]]- FALSE
indx
##Speed comparison
##previous method
fun1 - function(mat, vec) {
stopifnot(dim(mat)[1] == length(vec))
indx - mat == vec
t(apply(indx, 1, function(x) {
x[duplicated(x) !is.na(x)] - FALSE
x
}))
}
##modified
fun2 - function(mat, vec) {
stopifnot(dim(mat)[1] == length(vec))
indx - mat == vec
indx1 - which(indx, arr.ind = TRUE)
indx[indx1[duplicated(indx1[, 1]), ]] - FALSE
indx
}
identical(fun1(A,B), fun2(A,B))
#[1] TRUE
set.seed(498)
A1 - matrix(sample(40,1e5*500,replace=TRUE), ncol=500)
set.seed(345)
B1 - sample(70, 1e5, replace=TRUE)
system.time(res1 - fun1(A1,B1))
# user system elapsed
# 7.840 0.344 8.195
system.time(res2 - fun2(A1,B1))
# user system elapsed
# 0.304 0.080 0.382
identical(res1,res2)
#[1] TRUE
which(rowSums(res1,na.rm=TRUE)1)
#integer(0)
A.K.
On Friday, May 2, 2014 7:51 AM, arun smartpink...@yahoo.com wrote:
Hi,
Try:
indx - A==B
t(apply(indx,1,function(x) {x[duplicated(x) !is.na(x)] - FALSE; x}))
# [,1] [,2] [,3]
#[1,] TRUE FALSE FALSE
#[2,] FALSE NA FALSE
#[3,] NA NA NA
#[4,] TRUE NA FALSE
#[5,] FALSE TRUE FALSE
A.K.
On Friday, May 2, 2014 4:47 AM, nevil amos nevil.a...@gmail.com wrote:
I wish to return True in a matrix for only the first match of a value
per row where the value equals that in a vector with the same number of
values as rosw in the matrix
eg:
A-matrix(c(2,3,2,1,1,2,NA,NA,NA,5,1,0,5,5,5),5,3)
B-c(2,1,NA,1,5)
desired result:
[,1] [,2] [,3]
[1,] TRUE FALSE FALSE
[2,] FALSE NA FALSE
[3,] NA NA NA
[4,] TRUE NA FALSE
[5,] FALSE TRUE FALSE
however A==B returns:
[,1] [,2] [,3]
[1,] TRUE TRUE FALSE
[2,] FALSE NA FALSE
[3,] NA NA NA
[4,] TRUE NA FALSE
[5,] FALSE TRUE TRUE
and
apply(A,1,function(x) match (B,x))
returns
[,1] [,2] [,3] [,4] [,5]
[1,] 1 NA 1 NA NA
[2,] 3 NA NA 1 1
[3,] NA 2 2 2 NA
[4,] 3 NA NA 1 1
[5,] NA NA 3 3 2
thanks
[[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] SQL vs R
Thank you very much, Mr Arkell. el On 2014-05-03, 07:11 , Bert Gunter wrote: By making the effort to learn R? See e.g. the Introduction to R tutorial that ships with R. -- Bert Bert Gunter Genentech Nonclinical Biostatistics (650) 467-7374 Data is not information. Information is not knowledge. And knowledge is certainly not wisdom. H. Gilbert Welch On Fri, May 2, 2014 at 2:23 PM, Dr Eberhard Lisse nos...@lisse.na wrote: Hi, How do I do something like this without using sqldf? a - sqldf(SELECT COUNT(*) FROM b WHERE c = 'd') or e - sqldf(SELECT f, COUNT(*) FROM b GROUP BY f ORDER BY f) greetings, el __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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 me
hi please help me for below warning . may model is in frailtypack with frailtypenal mod.sha.gap - frailtyPenal(Surv(s,time,event) ~ cluster(id)+sex1+age2+bmi1+whr1+tg1+hdl1+chol1+por+pr2+car1+fib2+fat2+terminal(d),formula.terminalEvent = ~sex1+age2+bmi1+whr1+tg1+hdl1+chol1+por+pr2+car1+fib2+fat2,Frailty = TRUE, joint=TRUE,n.knots =4, kappa1 = 1,kappa2=1,data = data,recurrentAG = TRUE) I test many number for kappa1 and kappa2 and n.knote but have warning below Be patient. The program is computing ... The program took 0.52 seconds Warning message: In frailtyPenal(Surv(s, time, event) ~ cluster(id) + sex1 + age2 + : Problem in the loglikelihood computation. The program stopped abnormally. Please verify your dataset. please help me how can get a output? __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] auglag store parameters
Hi, I am minimizing a non-linear function subject to linear and nonlinear equalities and inequalities. I managed to understand how to use auglag for this case. Now, I'm using loops to get solutions for different assumptions, however, I'm facing a difficulty in storing just $par part of the answer to a vector form. Could someone point me to how to do this in auglag? Thank you for all your great work! [[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] auglag store parameters
Hi, I am minimizing a non-linear function subject to linear and nonlinear equalities and inequalities. I managed to understand how to use auglag for this case. Now, I'm using loops to get solutions for different assumptions, however, I'm facing a difficulty in storing just $par part of the answer to a vector form. Could someone point me to how to do this in auglag? Thank you for all your great work! [[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] Logistic regression help!
Hi guys, I have a trouble to solve the specificity and senstitivity for a logistic regression model. I really need your help, would you please help me out? :) Thank you!! This is the model I constructed: model=glm(Status ~ Gender.compl+ X2.4.times.per.month+ Once.a.month+Others+Council.tenant+Living.with.friends.family+Living.with.parents+Owner.occupier+Private.tenant+X1.year.to.2.years+X2.to.4.years+X3.months.to.a.year+Less.than.3.months+Over.4.years+Emp.compl+Reqloanamount+EmpNetMonthlyPay+Customer_Age+RTI+acc.compl+iic.compl+irb.compl+jbc.compl+jic.compl+jq.compl+kic.compl+lbc.compl+mbc.compl+njc.compl+or.compl+pq.compl+pr.compl+qic.compl+teb.compl+tpb.compl+vbc.compl+yzb.compl+zr.compl, data=learning.set1, family=binomial) so how to compute the sensititvity and specitivity? Many thanks! Siqi __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] error in r
dear all members
i have error in the code below Error in Y[i, 9] = 0.9 * XI[i, 2] + eps[9]
: subscript out of bounds is there anyone helps me please.
many thanks in advance
thanoon
llibrary(mvtnorm) #Load mvtnorm package
library(R2WinBUGS) #Load R2WinBUGS package
N=500 #Sample size
BZ=numeric(N) #Fixed covariate in structural equation
XI=matrix(NA, nrow=N, ncol=2) #Explanatory latent variables
Eta=numeric(N)#Outcome latent variables
Y=matrix(NA, nrow=N, ncol=8) #Observed variables
#The covariance matrix of xi
phi=matrix(c(1, 0.3, 0.3, 1), nrow=2)
#Estimates and standard error estimates
Eu=matrix(NA, nrow=100, ncol=10);SEu=matrix(NA, nrow=100, ncol=10)
Elam=matrix(NA, nrow=100, ncol=7); SElam=matrix(NA, nrow=100, ncol=7)
Eb=numeric(100); SEb=numeric(100)
Egam=matrix(NA, nrow=100, ncol=5); SEgam=matrix(NA, nrow=100, ncol=5)
Esgm=matrix(NA, nrow=100, ncol=10); SEsgm=matrix(NA, nrow=100, ncol=10)
Esgd=numeric(100); SEsgd=numeric(100)
Ephx=matrix(NA, nrow=100, ncol=3); SEphx=matrix(NA, nrow=100, ncol=3)
R=matrix(c(1.0, 0.3, 0.3, 1.0), nrow=2)
parameters=c(u, lam, b, gam, sgm, sgd, phx)
init1=list(u=rep(0,10), lam=rep(0,7), b=0, gam=rep(0,5), psi=rep(1,10),
psd=1, phi=matrix(c(1, 0, 0, 1), nrow=2))
init2=list(u=rep(1,10), lam=rep(1,7), b=1, gam=rep(1,5), psi=rep(2,10),
psd=2, phi=matrix(c(2, 0, 0, 2), nrow=2))
inits=list(init1, init2)
eps=numeric(10)
for (t in 1:100) {
#Generate Data
for (i in 1:N) {
BZ[i]=rt(1, 5)
XI[i,]=rmvnorm(1, c(0,0), phi)
delta=rnorm(1, 0, sqrt(0.36))
Eta[i]=0.5*BZ[i]+0.4*XI[i,1]+0.4*XI[i,2]+0.3*XI[i,1]*XI[i,2]
+0.2*XI[i,1]*XI[i,1]+0.5*XI[i,2]*XI[i,2]+delta
eps[1:3]=rnorm(3, 0, sqrt(0.3))
eps[4:7]=rnorm(4, 0, sqrt(0.5))
eps[8:10]=rnorm(3, 0, sqrt(0.4))
Y[i,1]=Eta[i]+eps[1]
Y[i,2]=0.9*Eta[i]+eps[2]
Y[i,3]=0.7*Eta[i]+eps[3]
Y[i,4]=XI[i,1]+eps[4]
Y[i,5]=0.9*XI[i,1]+eps[5]
Y[i,6]=0.7*XI[i,1]+eps[6]
Y[i,7]=0.5*XI[i,1]+eps[7]
Y[i,8]=XI[i,2]+eps[8]
Y[i,9]=0.9*XI[i,2]+eps[9]
Y[i,10]=0.7*XI[i,2]+eps[10]
}
#Run WinBUGS
data=list(N=500, zero=c(0,0), z=BZ, R=R, y=Y)
model-bugs(data,inits,parameters,
model.file=C:/Simulation/model.txt,
n.chains=2,n.iter=1,n.burnin=4000,n.thin=1,
bugs.directory=c:/Program Files/WinBUGS14/,
working.directory=C:/Simulation/)
#Save Estimates
Eu[t,]=model$mean$u; SEu[t,]=model$sd$u
Elam[t,]=model$mean$lam; SElam[t,]=model$sd$lam
Eb[t]=model$mean$b;SEb[t]=model$sd$b
Egam[t,]=model$mean$gam; SEgam[t,]=model$sd$gam
Esgm[t,]=model$mean$sgm; SEsgm[t,]=model$sd$sgm
Esgd[t]=model$mean$sgd;SEsgd[t]=model$sd$sgd
Ephx[t,1]=model$mean$phx[1,1]; SEphx[t,1]=model$sd$phx[1,1]
Ephx[t,2]=model$mean$phx[1,2]; SEphx[t,2]=model$sd$phx[1,2]
Ephx[t,3]=model$mean$phx[2,2]; SEphx[t,3]=model$sd$phx[2,2]
}
[[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] error in r
On 03/05/2014, 11:49 AM, thanoon younis wrote:
dear all members
i have error in the code below Error in Y[i, 9] = 0.9 * XI[i, 2] + eps[9]
: subscript out of bounds is there anyone helps me please.
You created Y with 8 columns, then you refer to column 9.
Duncan Murdoch
many thanks in advance
thanoon
llibrary(mvtnorm) #Load mvtnorm package
library(R2WinBUGS) #Load R2WinBUGS package
N=500 #Sample size
BZ=numeric(N) #Fixed covariate in structural equation
XI=matrix(NA, nrow=N, ncol=2) #Explanatory latent variables
Eta=numeric(N)#Outcome latent variables
Y=matrix(NA, nrow=N, ncol=8) #Observed variables
#The covariance matrix of xi
phi=matrix(c(1, 0.3, 0.3, 1), nrow=2)
#Estimates and standard error estimates
Eu=matrix(NA, nrow=100, ncol=10);SEu=matrix(NA, nrow=100, ncol=10)
Elam=matrix(NA, nrow=100, ncol=7); SElam=matrix(NA, nrow=100, ncol=7)
Eb=numeric(100); SEb=numeric(100)
Egam=matrix(NA, nrow=100, ncol=5); SEgam=matrix(NA, nrow=100, ncol=5)
Esgm=matrix(NA, nrow=100, ncol=10); SEsgm=matrix(NA, nrow=100, ncol=10)
Esgd=numeric(100); SEsgd=numeric(100)
Ephx=matrix(NA, nrow=100, ncol=3); SEphx=matrix(NA, nrow=100, ncol=3)
R=matrix(c(1.0, 0.3, 0.3, 1.0), nrow=2)
parameters=c(u, lam, b, gam, sgm, sgd, phx)
init1=list(u=rep(0,10), lam=rep(0,7), b=0, gam=rep(0,5), psi=rep(1,10),
psd=1, phi=matrix(c(1, 0, 0, 1), nrow=2))
init2=list(u=rep(1,10), lam=rep(1,7), b=1, gam=rep(1,5), psi=rep(2,10),
psd=2, phi=matrix(c(2, 0, 0, 2), nrow=2))
inits=list(init1, init2)
eps=numeric(10)
for (t in 1:100) {
#Generate Data
for (i in 1:N) {
BZ[i]=rt(1, 5)
XI[i,]=rmvnorm(1, c(0,0), phi)
delta=rnorm(1, 0, sqrt(0.36))
Eta[i]=0.5*BZ[i]+0.4*XI[i,1]+0.4*XI[i,2]+0.3*XI[i,1]*XI[i,2]
+0.2*XI[i,1]*XI[i,1]+0.5*XI[i,2]*XI[i,2]+delta
eps[1:3]=rnorm(3, 0, sqrt(0.3))
eps[4:7]=rnorm(4, 0, sqrt(0.5))
eps[8:10]=rnorm(3, 0, sqrt(0.4))
Y[i,1]=Eta[i]+eps[1]
Y[i,2]=0.9*Eta[i]+eps[2]
Y[i,3]=0.7*Eta[i]+eps[3]
Y[i,4]=XI[i,1]+eps[4]
Y[i,5]=0.9*XI[i,1]+eps[5]
Y[i,6]=0.7*XI[i,1]+eps[6]
Y[i,7]=0.5*XI[i,1]+eps[7]
Y[i,8]=XI[i,2]+eps[8]
Y[i,9]=0.9*XI[i,2]+eps[9]
Y[i,10]=0.7*XI[i,2]+eps[10]
}
#Run WinBUGS
data=list(N=500, zero=c(0,0), z=BZ, R=R, y=Y)
model-bugs(data,inits,parameters,
model.file=C:/Simulation/model.txt,
n.chains=2,n.iter=1,n.burnin=4000,n.thin=1,
bugs.directory=c:/Program Files/WinBUGS14/,
working.directory=C:/Simulation/)
#Save Estimates
Eu[t,]=model$mean$u; SEu[t,]=model$sd$u
Elam[t,]=model$mean$lam; SElam[t,]=model$sd$lam
Eb[t]=model$mean$b;SEb[t]=model$sd$b
Egam[t,]=model$mean$gam; SEgam[t,]=model$sd$gam
Esgm[t,]=model$mean$sgm; SEsgm[t,]=model$sd$sgm
Esgd[t]=model$mean$sgd;SEsgd[t]=model$sd$sgd
Ephx[t,1]=model$mean$phx[1,1]; SEphx[t,1]=model$sd$phx[1,1]
Ephx[t,2]=model$mean$phx[1,2]; SEphx[t,2]=model$sd$phx[1,2]
Ephx[t,3]=model$mean$phx[2,2]; SEphx[t,3]=model$sd$phx[2,2]
}
[[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] knitr - Highlight code/output
Hello,
I find situations where some segments of the code are displayed in the
output but not entirely highlighted. I guess there should be a way to fix
this, but I could't find it in the options.
Here's an example with knitr/LaTex.
\documentclass{article}
\begin{document}
\section{Example 1}
echo=TRUE, eval=TRUE=
set.seed(1)
x - matrix(rnorm(120), 12, 10)
x
@
\end{document}
Thanks,
Axel.
[[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] auglag store parameters
On May 3, 2014, at 5:07 AM, Delger Enkhbayar wrote: Hi, I am minimizing a non-linear function subject to linear and nonlinear equalities and inequalities. I managed to understand how to use auglag for this case. Now, I'm using loops to get solutions for different assumptions, however, I'm facing a difficulty in storing just $par part of the answer to a vector form. Could someone point me to how to do this in auglag? Let's imagine you are in London and ask the world's R audience how to get to a particular pub. But you don't tell them where you are and most of them are not Londoners so they don't know what the closest Tube stop might be, even if they looked up the pub's name on Google. There are probably many besides me out here that have no idea what 'auglag' or what package it might reside in, but do know enough R to possibly help if only there were an example offered. Thank you for all your great work! Part of that great work is the Posting Guide. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. David Winsemius Alameda, CA, 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] Seeking well-commented versions of mgcv source code
You can download the source code for the package with download.packages(mgcv, type=source, destdir=/tmp) Bill Dunlap TIBCO Software wdunlap tibco.com On Fri, May 2, 2014 at 7:15 PM, Andrew Crane-Droesch andre...@gmail.com wrote: Dear List, I'm looking for well-commented versions of various functions comprising mgcv, so that I can modify a piece of it for a project I'm working on. In particular I'm looking for * testStat * summary.gam * liu2 * simf Obviously I can find these by typing mgcv:::whatever. But there are a lot of nested if statements, making it difficult to follow. Comments in the code describing exactly what is happening at each step would make my life a lot easier. Where can I find more-detailed versions of the code? Thanks, Andrew [[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] mgcv gam bs=re questions
Dear R-helpers, I am working on a project assessing the prevalence and variance (random effects) of linear and nonlinear trends in a data set that has short time series (each time series identified as PID 1 through 5). I am using mgcv (gam) with the bs=”re” option (more on why not gamm or gamm4 at the end). I want to know whether trend is linear or nonlinear in each time series, and whether trend varies significantly over time series (e.g., PID 1 through 5 below)--e.g., whether they have significantly different edfs. The analysis predicts a count outcome (DVY) from time (SessIDX) and treatment (TX), using a quasipoisson approach to correct standard errors for overdispersion. The dataset in the examples below (SID5DVID1) has 5 time series (PID 1 through 5) with 93, 97, 70, 86, and 103 time points each, respectively. E.g., JID SID PID DVID SessIDX DVY TX 1 1 5 11 1 23 0 2 1 5 11 2 20 0 3 1 5 11 3 30 0 4 1 5 11 4 40 0 … 434 1 5 51 100 5 1 435 1 5 51 101 3 1 436 1 5 51 102 12 1 437 1 5 51 103 7 1 My question concerns how to interpret the random effects that result from variations in the model. A model with smoothed trend and no random effect yields M1 - gam(DVY ~ s(SessIDX) + factor(TX), +data = SID5DVID1, +family = quasipoisson(link=log), method=REML) Approximate significance of smooth terms: edf Ref.df F p-value s(SessIDX) 1.876 2.403 11.52 4.34e-06 *** So trend over all five time series is linear and significant. If I allow trend to interact with PID (be different for each time series), without a random effect, I get M1a - gam(DVY ~ s(SessIDX, by = factor(PID)) + factor(TX), data = SID5DVID1, family = quasipoisson(link=log), method=REML) summary(M1a,all.p=TRUE) Approximate significance of smooth terms: edf Ref.df F p-value s(SessIDX):factor(PID)1 1.001 1.001 31.827 2.98e-08 *** s(SessIDX):factor(PID)2 1.000 1.000 2.658 0.103724 s(SessIDX):factor(PID)3 5.100 5.886 8.674 9.42e-09 *** s(SessIDX):factor(PID)4 1.337 1.602 9.822 0.000366 *** s(SessIDX):factor(PID)5 3.998 4.924 7.967 4.40e-07 *** It appears trend may be quite different for each time series. So now I want to test whether trend varies significantly over time series using bs=”re”. One model is M2 - gam(DVY ~ s(SessIDX, bs = re) + factor(TX), data = SID5DVID1, family = quasipoisson(link=log), method=REML) summary(M2,all.p=TRUE) gam.vcomp(M2) Approximate significance of smooth terms: edf Ref.df F p-value s(SessIDX) 0.9609 1 24.48 6.53e-07 *** … Standard deviations and 0.95 confidence intervals: std.dev lower upper s(SessIDX) 0.01470809 0.00347556 0.06224258 scale 3.20933040 3.00277692 3.43009218 Rank: 2/2 So the random effect std.dev is 0.01470809, and it is significant with p-value = 6.53e-07. I had originally thought that this meant that trend varied significantly over the five time series. However, after having read a recent article (Gary J. McKeown and Ian Sneddon, DOI: 10.1037/a0034282), I am not sure I am correctly interpreting this. In their article, this approach (a) changes the shape of the nonlinear trend only slightly compared to not using random effects [i.e., just using s(SessIDX)] but still produces only one trend line for all time series, and (b) increases the width of the confidence bands around the trend line (which of course makes perfect sense). In other words, it does not obviously allow each time series (PID) to have a quite different edf or trend line. So, parallel with M1a above, I can run: M5a - gam(DVY ~ s(SessIDX, by = factor(PID), bs=re) + factor(TX), data = SID5DVID1, family = quasipoisson(link=log), method=REML) summary(M5a,all.p=TRUE) gam.vcomp(M5a) Approximate significance of smooth terms: edf Ref.df F p-value s(SessIDX):factor(PID)1 0.9239 1 38.76 0.000368 *** s(SessIDX):factor(PID)2 0.9892 1 107.00 2e-16 *** s(SessIDX):factor(PID)3 0.9867 1 94.85 2e-16 *** s(SessIDX):factor(PID)4 0.9830 1 108.33 1.64e-13 *** s(SessIDX):factor(PID)5 0.9658 1 73.41 1.32e-07 *** Standard deviations and 0.95 confidence intervals: std.dev lower upper s(SessIDX):factor(PID)1 0.008853821 0.001960054 0.03999387 s(SessIDX):factor(PID)2 0.065276527 0.016064261 0.26524874 s(SessIDX):factor(PID)3 0.048967099 0.012003711 0.19975296 s(SessIDX):factor(PID)4 0.027830597 0.006777380 0.11428341 s(SessIDX):factor(PID)5 0.012218325 0.002893741 0.05158977 scale 2.559988984 2.394424050 2.73700208 This seems to produce five random effects, all apparently significant. But it is not clear to me how to interpret these random effects. For
Re: [R] SQL vs R
Thanks, will try to figure this out :-)-O el On 2014-05-03, 06:40 , Carlos Ortega wrote: Hi, With the new package dplyr you can create equivalent SQL sintaxt queries like the one you need. You can find examples of how to apply it here: http://martinsbioblogg.wordpress.com/2014/03/26/using-r-quickly-calculating-summary-statistics-with-dplyr/ http://martinsbioblogg.wordpress.com/2014/03/27/more-fun-with-and/ Regards, Carlos. 2014-05-02 23:23 GMT+02:00 Dr Eberhard Lisse nos...@lisse.na mailto:nos...@lisse.na: Hi, How do I do something like this without using sqldf? a - sqldf(SELECT COUNT(*) FROM b WHERE c = 'd') or e - sqldf(SELECT f, COUNT(*) FROM b GROUP BY f ORDER BY f) greetings, el __ R-help@r-project.org mailto:R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Saludos, Carlos Ortega www.qualityexcellence.es http://www.qualityexcellence.es __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Seeking well-commented versions of mgcv source code
It turns out that commented versions of the source code are available on GitHub: https://github.com/cran/mgcv/blob/master/R/mgcv.r __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] SQL vs R
On 04/05/14 00:05, Dr Eberhard W Lisse wrote: Thank you very much, Mr Arkell. I don't get it. Can anyone explain the (joke? allusion?) ? cheers, Rolf Turner On 2014-05-03, 07:11 , Bert Gunter wrote: By making the effort to learn R? See e.g. the Introduction to R tutorial that ships with R. -- Bert Bert Gunter Genentech Nonclinical Biostatistics (650) 467-7374 Data is not information. Information is not knowledge. And knowledge is certainly not wisdom. H. Gilbert Welch On Fri, May 2, 2014 at 2:23 PM, Dr Eberhard Lisse nos...@lisse.na wrote: Hi, How do I do something like this without using sqldf? a - sqldf(SELECT COUNT(*) FROM b WHERE c = 'd') or e - sqldf(SELECT f, COUNT(*) FROM b GROUP BY f ORDER BY f) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] SQL vs R
On Sat, May 3, 2014 at 5:42 PM, Rolf Turner r.tur...@auckland.ac.nz wrote: On 04/05/14 00:05, Dr Eberhard W Lisse wrote: Thank you very much, Mr Arkell. I don't get it. Can anyone explain the (joke? allusion?) ? I believe it's a moderately offensive reply from someone who feels unfairly dismissed, derived from British pop culture. But someone who's actually British could better explain, I'm sure. Personally, I'm not sure how much work someone who appears to have not read the posting guide should really expect the list to do on his behalf. But snarky replies to reasonable requests to read the documentation are easier than doing one's own work. Sarah cheers, Rolf Turner On 2014-05-03, 07:11 , Bert Gunter wrote: By making the effort to learn R? See e.g. the Introduction to R tutorial that ships with R. -- Bert Bert Gunter Genentech Nonclinical Biostatistics (650) 467-7374 Data is not information. Information is not knowledge. And knowledge is certainly not wisdom. H. Gilbert Welch On Fri, May 2, 2014 at 2:23 PM, Dr Eberhard Lisse nos...@lisse.na wrote: Hi, How do I do something like this without using sqldf? a - sqldf(SELECT COUNT(*) FROM b WHERE c = 'd') or e - sqldf(SELECT f, COUNT(*) FROM b GROUP BY f ORDER BY f) __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- 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.
Re: [R] SQL vs R
Google Pressdram :-)-O el On 2014-05-03, 23:42 , Rolf Turner wrote: On 04/05/14 00:05, Dr Eberhard W Lisse wrote: Thank you very much, Mr Arkell. I don't get it. Can anyone explain the (joke? allusion?) ? cheers, Rolf Turner __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Question on Cluster Package, agnes() function
Anna F... on Thu, 1 May 2014 22:09:28 + writes: Hi Martin, I am a statistician at National Jewish Health in Colorado, and I have been working on clustering a dataset using Ward's minimum variance. When plotting the dendrogram, the y-axis is labeled as 'height'. Can you explain to me (or point me in the right direction) on how this distance between merging clusters is calculated for the Ward method? I have found the calculation that SAS uses, and I want to check if it is the same in your method. Here is a summary of the code I am using: Agnes(x,method=ward,diss=TRUE) Well, as R is case sensitive, it must be agnes(x,method=ward,diss=TRUE) Interestingly, the new version of R, R 3.1.0 has now two different versions of Ward in hclust() : -- http://stat.ethz.ch/R-manual/R-patched/library/stats/html/hclust.html where it is stated that previously it was basically not using Ward's method unless the user was calling it in a specific way, but agnes() was and is. *The* reference for all basic routines in the 'cluster' package is Kaufman, L. and Rousseeuw, P.J. (1990). _Finding Groups in Data: An Introduction to Cluster Analysis_. Wiley, New York. Alternatively, the source code of R and all packages is open, and for the cluster package, you can either get it from cluster_*.tar.gz from CRAN, or also you can see the (subversion) development version at http://svn.r-project.org/ Specifically, the C code which computes agnes() is https://svn.r-project.org/R-packages/trunk/cluster/src/twins.c and there, case 4: /* 4: ward's method */ ta = (double) kwan[la]; tb = (double) kwan[lb]; tq = (double) kwan[lq]; fa = (ta + tq) / (ta + tb + tq); fb = (tb + tq) / (ta + tb + tq); fc = -tq / (ta + tb + tq); int nab = ind_2(la, lb); dys[naq] = sqrt(fa * dys[naq] * dys[naq] + fb * dys[nbq] * dys[nbq] + fc * dys[nab] * dys[nab]); break; contains the distance calculation for ward. ... [ in private communication with Anna, she agreed that I reply publicly to R-help such that others can chime in and all will be searchable for people with a similar question. MM ] Best regards, Martin Maechler __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] SQL vs R
On 04/05/14 10:16, Dr Eberhard W Lisse wrote: Google Pressdram :-)-O el On 2014-05-03, 23:42 , Rolf Turner wrote: On 04/05/14 00:05, Dr Eberhard W Lisse wrote: Thank you very much, Mr Arkell. I don't get it. Can anyone explain the (joke? allusion?) ? Thank you. cheers, Rolf Turner __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] SQL vs R
On 04/05/14 09:58, Sarah Goslee wrote: SNIP Personally, I'm not sure how much work someone who appears to have not read the posting guide should really expect the list to do on his behalf. But snarky replies to reasonable requests to read the documentation are easier than doing one's own work. Well put. cheers, Rolf Turner __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] panel data: variable selection AND multicollinearity
Dear all, I am running a regression on panel data (127 observations). The Hausman Test indicats that the random-effects model is superior to the fixed-effects model. However, for a linear regression lm() there is a simple command in R to carry out a variable selection procedure: step(). Is there something similar for a regression on panel data plm()? Respectively, what is the command used for variable selection in case of plm? Is there any? Further, for linear regression lm() you can check easily for multicollinearity applying vif() of the car-package. What is the respective command for regression on panel data? Help is really appreciated. Looking forward to your answers. Thank you very much in advance. Regards, Phil -- View this message in context: http://r.789695.n4.nabble.com/panel-data-variable-selection-AND-multicollinearity-tp4689935.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] SQL vs R
By making the effort to learn R?? very constructive and not condescending at all. We, lesser beings, are indebted to you, sir. -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Bert Gunter Sent: Saturday, May 03, 2014 1:12 AM To: Dr Eberhard Lisse Cc: r Subject: Re: [R] SQL vs R By making the effort to learn R? See e.g. the Introduction to R tutorial that ships with R. -- Bert Bert Gunter Genentech Nonclinical Biostatistics (650) 467-7374 Data is not information. Information is not knowledge. And knowledge is certainly not wisdom. H. Gilbert Welch On Fri, May 2, 2014 at 2:23 PM, Dr Eberhard Lisse nos...@lisse.na wrote: Hi, How do I do something like this without using sqldf? a - sqldf(SELECT COUNT(*) FROM b WHERE c = 'd') or e - sqldf(SELECT f, COUNT(*) FROM b GROUP BY f ORDER BY f) greetings, el __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] SQL vs R
?table ?aggregate Also, packages plyr, data.table, and dplyr. You might consider reading [1], but if your interests are really as simple as your examples then the table function should be sufficient. That function is discussed in the Introduction to R document that you really should have read before posting here. [1] http://www.jstatsoft.org/v40/i01/ --- Jeff NewmillerThe . . Go Live... DCN:jdnew...@dcn.davis.ca.usBasics: ##.#. ##.#. Live Go... Live: OO#.. Dead: OO#.. Playing Research Engineer (Solar/BatteriesO.O#. #.O#. with /Software/Embedded Controllers) .OO#. .OO#. rocks...1k --- Sent from my phone. Please excuse my brevity. On May 2, 2014 2:23:13 PM PDT, Dr Eberhard Lisse nos...@lisse.na wrote: Hi, How do I do something like this without using sqldf? a - sqldf(SELECT COUNT(*) FROM b WHERE c = 'd') or e - sqldf(SELECT f, COUNT(*) FROM b GROUP BY f ORDER BY f) greetings, el __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] legendre quadrature
On 01-05-2014, at 09:43, pari hesabi statistic...@hotmail.com wrote:
Hello everybody
I need to approximate the amount of integral by using
legendre quadrature. I have written a program which doesn't give me a
logical answer; Can anybody help me and send the correct program? For
example the approximated amount of integral of ( x ^2) on (-1,1) based
on legendre quad rule.
integrand-function(x) {x^2}
rules - legendre.quadrature.rules( 50 )
order.rule - rules[[50]]
chebyshev.c.quadrature(integrand, order.rule, lower = -1, upper = 1)
You must be using package gaussquad.
Why are you using Legendre rules but doing Chebyshev quadrature (which does not
seem correct)?
Replace the last line of your given code with
legendre.quadrature(integrand, order.rule, lower = -1, upper = 1)
and the result will make more sense.
Berend
__
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/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] SQL vs R
On May 3, 2014, at 9:10 PM, Satish Anupindi Rao wrote: By making the effort to learn R?? very constructive and not condescending at all. We, lesser beings, are indebted to you, sir. For Pete's sake. The OP didn't even express his original request in natural language or offer a working example. Those of us who are not regular SQL users would have needed to parse out the SQL code in order to figure out what was intended. (My guess is that it would have been quite easy to solve if those were what were offered.) But making the effort to divine the intent didn't seem justified by the level of courtesy offered by the questioner. -- David. -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Bert Gunter Sent: Saturday, May 03, 2014 1:12 AM To: Dr Eberhard Lisse Cc: r Subject: Re: [R] SQL vs R By making the effort to learn R? See e.g. the Introduction to R tutorial that ships with R. -- Bert Bert Gunter Genentech Nonclinical Biostatistics (650) 467-7374 Data is not information. Information is not knowledge. And knowledge is certainly not wisdom. H. Gilbert Welch On Fri, May 2, 2014 at 2:23 PM, Dr Eberhard Lisse nos...@lisse.na wrote: Hi, How do I do something like this without using sqldf? a - sqldf(SELECT COUNT(*) FROM b WHERE c = 'd') or e - sqldf(SELECT f, COUNT(*) FROM b GROUP BY f ORDER BY f) greetings, el __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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. David Winsemius Alameda, CA, 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-es] Duda sobre como unificar datos en una oración
Hola a todos, tengo una duda sobre como unificar en R, datos de una oración Por ejemplo tengo la siguiente oración Siempre que, el atributo A es alto y el atributo B es media y el atributo C es alto y el atributo X es bajo; entonces el atributo Z es alto y quiero que me quede de la siguiente forma Siempre que, el atributo A y el atributo C es alto y el atributo B es medio y el atributo X es bajo; entonces el atributo Z es alto -- Daliana Ramos GarcÃa Universidad de las Ciencias Informáticas III Escuela Internacional de Invierno en la UCI del 17 al 28 de febrero del 2014. Ver www.uci.cu [[alternative HTML version deleted]] ___ R-help-es mailing list R-help-es@r-project.org https://stat.ethz.ch/mailman/listinfo/r-help-es
