[R] glmmADMB error
I ran an example using glmmADMB, and I got this error: data(bacteria,package=MASS) bacteria$present - as.numeric(bacteria$y)-1 (bfit - glmmadmb(present ~ trt + I(week 2), random = ~ 1 | ID, + family = binomial, data = bacteria)) Error in paste0(symbol1, paste0(paste0(var, collapse = symbol2))) : argument symbol1 is missing, with no default Any ideas? R version 3.0.2 (2013-09-25) Platform: x86_64-apple-darwin10.8.0 (64-bit) locale: [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8 attached base packages: [1] splines stats graphics grDevices utils datasets methods base other attached packages: [1] multcomp_1.3-1TH.data_1.0-3 survival_2.37-7 mvtnorm_0.9-9997 data.table_1.8.10 glmmADMB_0.7.7R2admb_0.7.10 [8] MASS_7.3-29 foreign_0.8-59 loaded via a namespace (and not attached): [1] grid_3.0.2 lattice_0.20-24 Matrix_1.1-1.1 nlme_3.1-113 sandwich_2.3-0 tools_3.0.2 zoo_1.7-10 -- Sebastián Daza [[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] CTM and survival analysis with heterogeneity
Hello R experts, I wonder if there is any package to estimate this kind of models in R: Multi-state Multi-spell Survival Models with Heterogeneity One of the most powerful programs for survival models is CTM, the Continuous Time Model, developed for NIH at NORC under the direction of James J. Heckman. Generalizing the competing risks model, CTM allows for transitions between any number of states, repeated spells within a state, time varying covariates, person and state specific heterogeneity, and arbitrary duration dependence specified in a highly flexible hazard model. Although in all cases the baseline hazard is fully parametric, it can be specified as a sufficiently rich function of time to capture a very wide set of duration dependencies. Any references? Thanks! -- Sebastián Daza __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] survey package question
Hello, I have got a cluster sample using an election dataset where I already had the final results of a county-specific election. I am trying to figure out what would be the best sampling design for my data. The structure of the dataset is: 1) polling station (in general schools where people vote, for a county, for example, there are 15 polling stations) 2) inside each polling station, there are voting units, where people actually vote (on average there are about 40 voting units for polling station) 3) for each voting unit I have the total votes by candidate (e.g., candidate 1 =322, candidate 2=122, candidate 3= 89) The initial sampling design is: 1) selection of 5 polling stations PPS (based on number of voters) 2) selection of 10 voting units (SRS) I am interested in estimating the proportion of votes by candidate (let's assume we have 3 candidates). My naive estimate would be: votes for candidate 1 / all valid votes = proportion e.g. candidate 1= 2132 / 10874= .1906 candidate 2= 5323 / 10874= .4895 candidate 3= 3419 / 10874= .3144 In this case, the unit of analysis is voters (or votes). If I specify the sampling design using the survey package in this way... design -svydesign(id=~station + unit fpc=~probstation +probunit, data=sample, pps=brewer) svyciprop(~I(candidate1/totalVotes), design) ... I am assuming that the unit of analysis is the voting unit, right? and I am estimating an average among voting units? I should expand my database at individual level (voters), or I just have to include a unit weight according to the number of voters for voting unit? In other words, is there a way to estimate, for instance, votes for candidate 1 / all valid votes = proportion, directly from the survey package or I have to expand the database at people level (voters), and then estimate the proportion using svymean and the respective design. I would appreciate any advice or help. Sebastian __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] survey package question
Hello Thomas, I use both svymean (with the expanded sample = people), and svyratio (voting unit level), using the same design: design -svydesign(id=~station + unit, fpc=~probstation+probunits, data=sample, pps=brewer) I got different results using the same sample: svyratio (voting unit) Ratio 2.5% 97.5% Result Cand1 0.05252871 0.04537301 0.05968441 0.05181146 Cand20.47226973 0.45215097 0.49238849 0.49041590 Cand3 0.47520156 0.45460831 0.49579482 0.45777264 svymean (expanded sample, individuals or votes) Mean SE 2.5 % 97.5 %Results Cand1 0.0528433 0.004562755 0.04390047 0.06178614 0.05181146 Cand2 0.4717504 0.010201398 0.45175605 0.49174480 0.49041590 Cand30.4754063 0.010429222 0.45496538 0.49584718 0.45777264 Point estimators are different, and confidence intervals are more narrow using svyratio. Could you give me any clue about what is going on? Thank you in advance. Sebastian On Thu, Oct 11, 2012 at 7:50 PM, Sebastián Daza sebastian.d...@gmail.com wrote: Hello Thomas, I use both svymean (with the expanded sample = people), and svyratio (voting unit level), using the same design: design -svydesign(id=~station + unit, fpc=~probstation+probunits, data=sample, pps=brewer) I got different results using the same sample: svyratio (voting unit) Ratio 2.5% 97.5% Result Cand1 0.05252871 0.04537301 0.05968441 0.05181146 Cand20.47226973 0.45215097 0.49238849 0.49041590 Cand3 0.47520156 0.45460831 0.49579482 0.45777264 svymean (expanded sample, individuals or votes) Mean SE 2.5 % 97.5 %Results Cand1 0.0528433 0.004562755 0.04390047 0.06178614 0.05181146 Cand2 0.4717504 0.010201398 0.45175605 0.49174480 0.49041590 Cand30.4754063 0.010429222 0.45496538 0.49584718 0.45777264 Point estimators are different, and confidence intervals are more narrow using svyratio. Could you give me any clue about what is going on? Thank you in advance. Sebastian On Thu, Oct 11, 2012 at 3:56 PM, Sebastián Daza sebastian.d...@gmail.com wrote: Thank you Thomas! On Thu, Oct 11, 2012 at 2:33 PM, Thomas Lumley tlum...@uw.edu wrote: On Fri, Oct 12, 2012 at 6:56 AM, Sebastián Daza sebastian.d...@gmail.com wrote: Hello, I have got a cluster sample using an election dataset where I already had the final results of a county-specific election. I am trying to figure out what would be the best sampling design for my data. The structure of the dataset is: 1) polling station (in general schools where people vote, for a county, for example, there are 15 polling stations) 2) inside each polling station, there are voting units, where people actually vote (on average there are about 40 voting units for polling station) 3) for each voting unit I have the total votes by candidate (e.g., candidate 1 =322, candidate 2=122, candidate 3= 89) The initial sampling design is: 1) selection of 5 polling stations PPS (based on number of voters) 2) selection of 10 voting units (SRS) I am interested in estimating the proportion of votes by candidate (let's assume we have 3 candidates). My naive estimate would be: votes for candidate 1 / all valid votes = proportion e.g. candidate 1= 2132 / 10874= .1906 candidate 2= 5323 / 10874= .4895 candidate 3= 3419 / 10874= .3144 In this case, the unit of analysis is voters (or votes). If I specify the sampling design using the survey package in this way... design -svydesign(id=~station + unit fpc=~probstation +probunit, data=sample, pps=brewer) svyciprop(~I(candidate1/totalVotes), design) ... I am assuming that the unit of analysis is the voting unit, right? and I am estimating an average among voting units? You want a ratio estimator svyratio(~candidate1, ~totalVotes, design) -thomas -- Thomas Lumley Professor of Biostatistics University of Auckland -- Sebastián Daza -- Sebastián Daza -- Sebastián Daza __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] Loop question
Hello, I have a dataframe (Lx) with 5 Lb, and 5 Lw variables. I want to create several variables according to the loop presented below (data attached). Lx - read.csv(Lx.csv, header=T, sep=,) for (i in 1:20) { Lx$sb1[i] - Lx$Lb1[i+1]/Lx$Lb1[i] Lx$sb2[i] - Lx$Lb2[i+1]/Lx$Lb2[i] Lx$sb3[i] - Lx$Lb3[i+1]/Lx$Lb3[i] Lx$sb4[i] - Lx$Lb4[i+1]/Lx$Lb4[i] Lx$sb5[i] - Lx$Lb5[i+1]/Lx$Lb5[i] Lx$sw1[i] - Lx$Lw1[i+1]/Lx$Lw1[i] Lx$sw2[i] - Lx$Lw2[i+1]/Lx$Lw2[i] Lx$sw3[i] - Lx$Lw3[i+1]/Lx$Lw3[i] Lx$sw4[i] - Lx$Lw4[i+1]/Lx$Lw4[i] Lx$sw5[i] - Lx$Lw5[i+1]/Lx$Lw5[i] } How I could also create the variable names (letters b and w, and numbers from 1 to 5 in s and L variables) using a loop in R? In Stata I can use this command: foreach r in b w {; foreach s of numlist 0(5)40 {; foreach ed of numlist 1/5 {; local f = `s'+5; scalar define rL`r'`ed'_`f'over`s' = L`r'`ed'_`f' / L`r'`ed'_`s'; }; }; }; but I do not how to do it in R. Thank you in advance! -- Sebastián Daza __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] trellis plot
Thank you Peter. The problem with you solution is that it doesn't represent the actual values of disruption. Look this example: person - rep(1:2, each=4) income - c(100, 120, 150, 200, 90, 100,120, 150) disruption - c(0,0,0,1,0,1,1,0) time - rep(c(1:4),2) dat - as.data.frame(cbind(person,time, income, disruption)) mycols - c(2, 5) library(lattice) xyplot(income~time|as.factor(person),data=dat, type=c(p,g,o), col.line=black, col.symbol = mycols[dat$disruption+1], xlab=Time, ylab=Familiar Income) I am looking for a way to get the colors of dots according to actual disruption values, and at the same time to draw the lines between all dots for each person. Thank you! Hi everyone, I am just trying to figure out how to do a xyplot where in addition to dots and lines I can change dots' colors according to an individual variable (e.g., marital disruption across time, a dummy 0/1). When I use groups specification (see below), I get two different lines for each individual based on groups, and what I want is to get one line connecting dots, and different dots' colors according to marital disruption. Any ideas about how to do that? person- rep(1:2, each=4) income- c(100, 120, 150, 200, 90, 100,120, 150) disruption- rep(c(0,1), 4) time- rep(c(1:4),2) dat- as.data.frame(cbind(person,time, income, disruption)) library(lattice) xyplot(income~time|as.factor(person),data=dat, type=c(p,g,o), col.line=black, xlab=Time, ylab=Familiar Income) # I just want to change dots' colors according to disruption, not to get two different lines: xyplot(income~time|as.factor(person),data=dat, type=c(p,g,o), col.line=black, groups=disruption, xlab=Time, ylab=Familiar Income) If I understand correctly what you're after, try adding the 'col.symbol' argument instead of the 'groups' argument. col.symbol = disruption + 1 or, better, for arbitrary colours: mycols - c(2, 5) xyplot(, col.symbol = mycols[disruption + 1], ) Peter Ehlers Thank you in advance! -- Sebastián Daza -- Sebastián Daza __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] trellis plot
Hi everyone, I am just trying to figure out how to do a xyplot where in addition to dots and lines I can change dots' colors according to an individual variable (e.g., marital disruption across time, a dummy 0/1). When I use groups specification (see below), I get two different lines for each individual based on groups, and what I want is to get one line connecting dots, and different dots' colors according to marital disruption. Any ideas about how to do that? person - rep(1:2, each=4) income - c(100, 120, 150, 200, 90, 100,120, 150) disruption - rep(c(0,1), 4) time - rep(c(1:4),2) dat - as.data.frame(cbind(person,time, income, disruption)) library(lattice) xyplot(income~time|as.factor(person),data=dat, type=c(p,g,o), col.line=black, xlab=Time, ylab=Familiar Income) # I just want to change dots' colors according to disruption, not to get two different lines: xyplot(income~time|as.factor(person),data=dat, type=c(p,g,o), col.line=black, groups=disruption, xlab=Time, ylab=Familiar Income) Thank you in advance! -- Sebastián Daza __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] calculations combining values from different rows
Hi everyone, I looking for functions or systematic ways to do calculations between different columns and rows in R, as one can do easily in Excel. For example, I have two variables, a and b, where a1 represents an a value in row 1, and b2 represents a b value in row 2, etc. a - c(4,3,5,5,6,7,3,2,1,4) b - c(2,4,1,2,5,3,1,8,7,5) data - cbind(a,b) I have to calculate something like this: x1 = NA x2 = -b1 /24 * a1 + b2 /2 * a2 + b3 /24 * a3 x3 = -b2 /24 *a2 + b3 /2 * a3 + b4 /24 * a4 x4 = -b3 /24 *a3 + b4 /2 * a4 + b5 /24 * a5 ... x9 = -b8 /24* a8 + b9 /2 * a9 + b10 /24 * a10 x10= NA For example, x2 would be equal to: -2/24*4 +4/2*3 + 1/24 *5 Any ideas? Thank you in advance. -- Sebastián Daza __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] data management question
Hi everyone, I would like to have a function to compute some values in a dataset. First, I have to define a value for the lx variable in row 1 (e.g., 100,000), npx is a given proportion. lx of row 2 is equal to lx of row 1 times npx of row 1. I can do this row by row... data[1,lx] - 10 data[2,lx] - data[1,lx]*data[1,npx] data[3,lx] - data[2,lx]*data[2,npx] data[4,lx] - data[3,lx]*data[3,npx] ... data[19,lx] - data[18,lx]*data[18,npx] Any ideas about how to define this in a function or in a more systematic way? Thank you in advance. -- Sebastián Daza __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] raking weighting
Hi everyone, Does anyone know if there is a package to compute raking weights using R? What I need is to create a variable with weights base in some demographic variables (e.g. sexo, age group, area) using the raking procedure. Thank you in advance! -- Sebastián Daza __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] raking weighting
I did. However, I am not sure that I can create a weight variable in my database with the rake function you mention. That is why I am asking. El día 9 de noviembre de 2011 13:44, David Winsemius dwinsem...@comcast.net escribió: On Nov 9, 2011, at 2:27 PM, Sebastián Daza wrote: Hi everyone, Does anyone know if there is a package to compute raking weights using R? What I need is to create a variable with weights base in some demographic variables (e.g. sexo, age group, area) using the raking procedure. There is a `rake` function in package survey. In the future you ought to do this at the R console _before_ sending a question to R-help: RSiteSearch(raking) Thank you in advance! -- Sebastián Daza __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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, MD West Hartford, CT -- Sebastián Daza __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] identifying groups
Hi everyone, Does anyone know if there is a package to perform a Moody's Crowds routine to identify groups using R, or other algorithms designed to search groups by maximizing modularity scores? Otherwise, other packages to identify groups, or other programs to do it. I have about 80 networks clustered in schools. So, I have to identify groups for each school. Thank you in advance! -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] data network format and grouping analysis
Hi everyone, I have a dataset of friendship with this format: ego alter 47461 2 97421 3 14738 1NA 47472NA 974323 14739 21 4748313 97443 5 14740 314 47494NA 97454NA 14741 4NA 47505NA 9746513 14742 510 47516 12 97476 7 ... NA means that individuals don't select any friend. Does anyone know how to format this dataset to use sna or igraph packages? I don't know how to convert it into a matrix or a edgelist in R without losing isolated individuals . Next question, anyone knows if there is a package to perform a Moody's Crowds routine to identify groups using R, or other algorithms designed to search groups by maximizing modularity scores? Thank you in advance! -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] setCoefTemplate
Hi everyone, I am trying to build a table putting standard errors horizontally. I haven't been able to do it. library(memisc) berkeley - aggregate(Table(Admit,Freq)~.,data=UCBAdmissions) berk0 - glm(cbind(Admitted,Rejected)~1,data=berkeley,family=binomial) berk1 - glm(cbind(Admitted,Rejected)~Gender,data=berkeley,family=binomial) berk2 - glm(cbind(Admitted,Rejected)~Gender+Dept,data=berkeley,family=binomial) setCoefTemplate(est.se=c(est = ($est:#)($se:#))) mtable(berk0,berk1,berk2, + coef.style=est.se, + summary.stats=c(Deviance,AIC,N)) Error in dim(ans) - newdims : dims [product 1] do not match the length of object [2] Thank you in advance. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] setCoefTemplate
Hi everyone, I am trying to build a table putting standard errors horizontally. I haven't been able to do it. library(memisc) berkeley - aggregate(Table(Admit,Freq)~.,data=UCBAdmissions) berk0 - glm(cbind(Admitted,Rejected)~1,data=berkeley,family=binomial) berk1 - glm(cbind(Admitted,Rejected)~Gender,data=berkeley,family=binomial) berk2 - glm(cbind(Admitted,Rejected)~Gender+Dept,data=berkeley,family=binomial) setCoefTemplate(est.se=c(est = ($est:#)($se:#))) mtable(berk0,berk1,berk2, + coef.style=est.se, + summary.stats=c(Deviance,AIC,N)) Error in dim(ans) - newdims : dims [product 1] do not match the length of object [2] Thank you in advance. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] coefficients style
Hi everyone, I am trying to build a table putting standard errors horizontally. I haven't been able to do it. library(memisc) berkeley - aggregate(Table(Admit,Freq)~.,data=UCBAdmissions) berk0 - glm(cbind(Admitted,Rejected)~1,data=berkeley,family=binomial) berk1 - glm(cbind(Admitted,Rejected)~Gender,data=berkeley,family=binomial) berk2 - glm(cbind(Admitted,Rejected)~Gender+Dept,data=berkeley,family=binomial) setCoefTemplate(est.se=c(est = ($est:#)($se:#))) mtable(berk0,berk1,berk2, + coef.style=est.se, + summary.stats=c(Deviance,AIC,N)) Error in dim(ans) - newdims : dims [product 1] do not match the length of object [2] Thank you in advance. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] style question
Hi everyone, I am trying to build a table putting standard errors horizontally. I haven't been able to do it. library(memisc) berkeley - aggregate(Table(Admit,Freq)~.,data=UCBAdmissions) berk0 - glm(cbind(Admitted,Rejected)~1,data=berkeley,family=binomial) berk1 - glm(cbind(Admitted,Rejected)~Gender,data=berkeley,family=binomial) berk2 - glm(cbind(Admitted,Rejected)~Gender+Dept,data=berkeley,family=binomial) setCoefTemplate(est.se=c(est = ($est:#)($se:#))) mtable(berk0,berk1,berk2, + coef.style=est.se, + summary.stats=c(Deviance,AIC,N)) Error in dim(ans) - newdims : dims [product 1] do not match the length of object [2] Thank you in advance. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] simple if question
Hi everyone, I have just got different samples from a dataframe (independent and exclusive, there aren't common elements among them). I want to create a variable that indicate the sampling selection of the elements in the original dataframe (for example, 0 = no selected, 1= sample 1, 2=sample 2, etc.). I have tried to do it with ifelse command, but the problem is that the second line replaces the values of the first line, and I haven't been able to do it with the if command (I got this error: In if (data$school %in% sample1) { : the condition has length 1 and only the first element will be used) data$selection - ifelse(data$school %in% sample1, 1, 0) data$selection - ifelse(data$school %in% sample2, 2, 0) Any ideas? Thank you in advance. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] pattern of panel waves using R
Finally, I solved my problem using the following procedure and a database called ej ej$a - 1 head(ej) ano nunico a 1 2008 1 1 2 2009 1 1 3 2008 2 1 4 2009 2 1 5 2008 3 1 6 2009 3 1 library(reshape) dej - cast(ej, nunico ~ ano, sum, margins = FALSE) head(dej) nunico 2008 2009 2010 1 1110 2 2110 3 3110 4 4111 5 5101 6 6111 dej[dej==0] - NA library(Hmisc) na.pattern(dej[,c(2:4)]) pattern 000 001 010 011 110 3385 1073 203 338 573 That way I can review the pattern of my panel data. Sebastian __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] pattern of panel waves using R
Hi everyone, Is there any command to identify the pattern of responses of a database with this format: year id 20081 20091 20082 20092 20083 20093 20084 20094 20104 I just need the frequency of the patterns grouped by id: 2008 2009 2010 = 80 2009 2010 = 30 2008 2009 = 10 and so on Thank you in advance! -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] lavaan diagram
Hi everyone, I got the following error when I tried to diagram a confirmatory factor analysis using lavaan: lavaan.diagram(conf) Error in strwidth(xvars) : plot.new has not been called yet Does anyone have any idea about what is happening? Thank you in advance. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] lavaan diagram
Hi everyone, I get the following message when I try to get a diagram of a CFA: Error in lavaan.diagram(fit) : no slot of name GLIST for this object of class Fit Does anyone know what the problem is? Regards. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] aggregate in R
plyr is very useful to aggregate data... I strongly recommend it. On 2/22/2011 5:59 PM, Jeff Newmiller wrote: Use ?plyr::ddply --- Jeff Newmiller The . . Go Live... DCN:jdnew...@dcn.davis.ca.us Basics: ##.#. ##.#. Live Go... Live: OO#.. Dead: OO#.. Playing Research Engineer (Solar/Batteries O.O#. #.O#. with /Software/Embedded Controllers) .OO#. .OO#. rocks...1k --- Sent from my phone. Please excuse my brevity. Hongwei Dongpdxd...@gmail.com wrote: Hi, R users, I'm wondering how I can aggregate data in R with different functions for different columns. For example: x-rep(1:5,3) y-cbind(x,a=1:15,b=21:35) y-data.frame(y) I want to aggregate a and b in y by x. With a, I want to use function mean; with b, I want to use function sum. I tried: aggregate(y,x,mean(y$a),sum(y$b)) But I got the error: Error in match.fun(FUN) : 'mean(y$a)' is not a function, character or symbol Anyone can tell me how to fix this problem? Thanks. Gary [[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. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] id number by group and correlative
Hello everyone, I am new in R and I am trying to create a id number (a correlative sequence of numbers) by group, and a correlative sequence of numbers inside each group (my idea is to get statistics by group without having to aggregate the database). Here an example: group id_groupcorrelative_group A 1 1 A 1 2 A 1 3 A 1 4 B 2 1 B 2 2 B 2 3 C 3 1 C 3 2 C 3 3 C 3 4 C 3 5 Unfortunately, I have been able to find an explicit lag function to get id_group (I know I can get it using aggregate and merge but I'm just wondering if there is another way to do it). With regard to the correlative_group, I don't have any clue about how to do it. PD: Last question. Is there any way to save in a variable the internal that R uses (those numbers that appear on left side of a dataframe)? Thank you, Sebastian. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] very basic HLM question
Thank you for your reply and sorry for my ambiguity. I computed: summary(anova1 - aov(math ~ as.factor(schoolid), data=nels88)) ICC1(anova1) ICC1 comes from the multilevel package. I have found an article where it is pointed out that with this formula: 34.011/(34.011+72.256) [1] 0.3200523 I should get practically an identical result than using ICC1. I have used the database of that article and I have got an identical result. But with the dataset I am using (schools) it doesn't work. This is the example of the article mentioned above: library(multilevel) base(bh1996) summary(lmer(WBEING ~ 1 + (1|GRP), data=bh1996)) 0.035801/(0.035801+0.789497) summary(test - aov(WBEING~as.factor(GRP),data=bh1996)) ICC1(test) I have attached the small database where this equality doesn't exist. Thank you in advance. On 2/5/2011 11:07 PM, Paul Johnson wrote: 2011/2/5 Sebastián Dazasebastian.d...@gmail.com: Hi everyone, I need to get a between-component variance (e.g. random effects Anova), but using lmer I don't get the same results (variance component) than using random effects Anova. I am using a database of students, clustered on schools (there is not the same number of students by school). According to the ICC1 command, the interclass correlation is .44 ICC1(anova1) [1] 0.4414491 If you don't tell us exactly what model you are calculating in anova1, how would we guess if there is something wrong? Similarly, I get this ICC1 Error: object 'ICC1' not found so it must mean you've loaded a package or written a function, which you've not shown us. I googled my way to a package called multilevel that has ICC1, and its code for ICC1 shows a formula that does not match the one you used to calculate ICC from lmer. function (object) { MOD- summary(object) MSB- MOD[[1]][1, 3] MSW- MOD[[1]][2, 3] GSIZE- (MOD[[1]][2, 1] + (MOD[[1]][1, 1] + 1))/(MOD[[1]][1, 1] + 1) OUT- (MSB - MSW)/(MSB + ((GSIZE - 1) * MSW)) return(OUT) } I'm not saying that's right or wrong, just not obviously identical to the formula you proposed. However, I cannot get the same ICC from the lmer output: anova2- lmer(math ~ 1 + (1|schoolid), data=nels88) summary(anova2- lmer(math ~ 1 + (1|schoolid), data=nels88)) Instead, do this (same thing, fits model only once): anova2- lmer(math ~ 1 + (1|schoolid), data=nels88) summary(anova2) Note that lmer is going to estimate a normally distributed random effect for each school, as well as an individual observation random effect (usual error term) that is assumed independent of the school-level effect. What is anova1 estimating? -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] very basic HLM question
Hi everyone, I need to get a between-component variance (e.g. random effects Anova), but using lmer I don't get the same results (variance component) than using random effects Anova. I am using a database of students, clustered on schools (there is not the same number of students by school). According to the ICC1 command, the interclass correlation is .44 ICC1(anova1) [1] 0.4414491 However, I cannot get the same ICC from the lmer output: anova2 - lmer(math ~ 1 + (1|schoolid), data=nels88) summary(anova2 - lmer(math ~ 1 + (1|schoolid), data=nels88)) Linear mixed model fit by REML Formula: math ~ 1 + (1 | schoolid) Data: nels88 AIC BIC logLik deviance REMLdev 1878 1888 -935.8 18751872 Random effects: Groups NameVariance Std.Dev. schoolid (Intercept) 34.011 5.8319 Residual 72.256 8.5003 Number of obs: 260, groups: schoolid, 10 Fixed effects: Estimate Std. Error t value (Intercept) 48.861 1.927 25.36 The intercept random effect is 34.011. If I compute the ICC manually I get: 34.011/(34.011+72.256) [1] 0.3200523 According to my Anova analysis, the between-component variance is 59.004. Does anyone know what the problem is? How can I get the 59.004 figure using R? -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] grey scale graphs
Hi everyone, Does anyone know how to get black and white theme (grey scale,, I would say) graphs using lattice or ggplot2, as it is shown in this webpage: http://lmdvr.r-forge.r-project.org/figures/figures.html? I am using Sweave, and I cannot get that color configuration. I have added the following option: trellis.device(color=FALSE) but I got a pdf file with color graphs. Thank in advance. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Latent Class Logit Models in discrete choice experiments
See: https://www.msu.edu/~chunghw/downloads.html Maybe you can find something useful there! Regards On 1/31/2011 12:35 PM, Daniel Vecchiato wrote: Dear R users, I would like to perform Latent Class Logit Models for the analysis of choice experiments in environmental valuation. This kind of analysis is usually performed with NLogit Software (http://www.limdep.com). I attach the results I usually obtain using NLogit and NLogit model specifications. For Random parameter models and Logit Models I usually perform my analysis with the package mlogit ( http://cran.r-project.org/web/packages/mlogit/index.html ). The models I would like to run are presented in this ppt: http://pages.stern.nyu.edu/~wgreene/.../Lectures/Part13-LatentClassModels.ppt and an overview is given in this paper by Hoyos: http://dx.doi.org/10.1016/j.ecolecon.2010.04.011 Anybody has any package to suggest for this kind of analysis? (poLCA does not provide me the same estimates) Thanks in advance Daniel __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] rj packages and eclipse
Hi everyone! Does anyone know how to install the rj package in Windows (7)? I am trying to set up eclipse, but I got [INFO] The R package 'rj' is not available, R-StatET tools cannot be initialized. Thank you in advance! __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] LTA
Hi everyone, Does anyone know if there is a package to run Latent Transitional Analysis using R? Regards! -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Using summaryBy with weighted data
Hi everyone, I am trying to run Sweave.bat (batchfiles_0.6-1) from the command line on Windows, but I get this error: C:\batchfiles_0.6-1Sweave.bat Sweave-test-1 Error: rterm.exe not found I don't know how to set up the path if this one were the problem... I ran rcmd.bat and I got this... so I don't know if it is a path problem. C:\batchfiles_0.6-1Rcmd,bat R_ARCH=/x64 R_ARCH0=x64 R_ARCH0=x64 cmdpath=C:\R\R-2.12.1\bin\x64\Rcmd.exe args=,bat 'bat' is not recognized as an internal or external command, operable program or batch file. the path of rterm.exe in my computer is: C:\R\R-2.12.1\bin\x64 thank you in advance! -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] sweave.bat
Hi everyone, I am trying to run Sweave.bat (batchfiles_0.6-1) from the command line on Windows, but I get this error: C:\batchfiles_0.6-1Sweave.bat Sweave-test-1 Error: rterm.exe not found I don't know how to set up the path if this one were the problem... I ran rcmd.bat and I got this... so I don't know if it is a path problem. C:\batchfiles_0.6-1Rcmd,bat R_ARCH=/x64 R_ARCH0=x64 R_ARCH0=x64 cmdpath=C:\R\R-2.12.1\bin\x64\Rcmd.exe args=,bat 'bat' is not recognized as an internal or external command, operable program or batch file. the path of rterm.exe in my computer is: C:\R\R-2.12.1\bin\x64 thank you in advance! -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] easy loop question
Hi everyone, I am new in R and programming. I have tried to remove the values out of range in some variables using a loop: 1) var - names(est8vo[, 77:83]) # I got the variable names var [1] p16.1 p16.2 p16.3 p16.4 p16.5 p16.6 p16.7 for (i in 1:7) { var.i - var[i] est8vo$var.i[ est8vo$var.i==3] - 99 } I got this error: Error in `$-.data.frame`(`*tmp*`, var.i, value = numeric(0)) : replacement has 0 rows, data has 215700 2) The second step would be to define the factors, but I got the same error: for (i in 1:7) { var.i - var[i] est8vo$var.i- factor(est8vo$var.i, levels=c(0, 1, 2, 99), labels=c(vacío, sí, no, doble marca) ) } I don't know how to do it. Thank you in advance! Sebastian __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/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] integration Sweave and TexMakerX
Hi, Does anyone know how to integrate texmakerx and sweave on Windows? I mean, to run .rnw files directly from texmakerx and get a pdf or dvi file. Thank you in advance, -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] two-part growth analysis
Hi everyone! Does anyone know if there is a package to do two-part growth analysis with R? Regards, Sebastian -- Sebastián Daza sebastian.d...@gmail.com __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.