On Wed, 19 May 2010, Sabatier, Jennifer F. (CDC/OID/NCHHSTP) wrote:
Hi R-help,
Yes, this is my second request for assistance in a single day....
I am attempting to use svychisq() inside a function I made. The goal
of this function is to produce a table of summary statistics that I can
later output to EXCEL (simple frequencies and sample sizes from regular
crosstabulation on dataset "data" but the chi-square using survey
methods on "audit").
Here's my code (I can't supply data for you as I am not that
sophisticated and the real data is not cleared for public consumption -
I really apologize):
You could have checked to see if the same problem occurs using one of the
built-in data sets, but it turns out that the diagnosis is fairly
straightforward. Also, since you use the function 'crosstab' without saying
which package you got it from, it might be hard to get your code to run anyway.
Your problem is that you are passing a bare name and trying to substitute it
into a formula, and that isn't how R function arguments or formulas work.
Simplifying the problem to the essentials, your function would work only if
myFormulaMakingFunction<-function(svyX){
~svyX+SEX
}
myFormulaMakingFunction(con)
returned ~con+SEX. It actually returns ~svyX+SEX
There are various ways around this. My preferred one would be
myPreferredWay<-function(formula){
update(formula, ~.+SEX)
}
myPreferredWay(~con)
which does return ~con+SEX.
Then, as an example:
data(api)
dclus1<-svydesign(id=~dnum, weights=~pw, data=apiclus1, fpc=~fpc)
myIllustrativeExample<-function(X,formula){
tab<-table(X,apiclus1$stype)
chisq<-svychisq(update(formula, ~.+stype), statistic="adjWald",design=dclus1)
list(tab,chisq)
}
myIllustrativeExample(apiclus1$comp.imp, ~comp.imp)
Now, this can be improved: it appears that the vector argument is supposed to
be the same variable as the formula argument, so we shouldn't make the user
supply both of them.
mySupererogatoryEffort <- function(formula, design){
X<-model.frame(formula, model.frame(design))[[1]]
tab<-table(X, model.frame(design)$stype)
chisq<-svychisq(update(formula, ~.+stype), statistic="adjWald",design=design)
list(table=tab, F=chisq$statistic, df=chisq$parameter, p=chisq$p.value)
}
mySupererogatoryEffort(~comp.imp, dclus1)
As a final note, there is no round= argument to svychisq(). There is a round=
argument to svytable(), which is documented on the same help page, but 4 is not
one of the allowed values.
-thomas
# create my svydesign object
audit <- svydesign(id~id, strata=~field, weights=~wt, data=data,
fpc=~AllocProportion)
# my function to create my table
mkMyCrossTable <- function(X, svyX, T) {
tbl <- crosstab(X, data$SEX, prop.c=TRUE)
tbl <- data.frame(cbind(tbl$t, tbl$prop.col))
tbl$var <- rownames(tbl)
chisq <- svychisq(~svyX + SEX, design=audit, statistic="adjWald",
round=4)
chisq <- data.frame(do.call("cbind", chisq)
chisq <- data.frame(chisq[,3])
Table <- data.frame(tbl$var,
paste(formatC(tbl$X0.1*100, format="f", digits=1),
"%", sep=""),
tbl$X0,
paste(formatC(tbl$X1.1*100, format="f",
digits=1), "%", sep=""),
tbl$X1,
chisq[1])
Table[2: length(Table[,1]), 6] <- NA
Table <- NAToUnkown(Table, unknown = " ")
Colnames(Table) <- c(T, "Male (%)", "Male (n)", "Female (%)", "Female
(n)", "p-value")
Table
}
con3 <- mkMyCrossTable(data$con, con, "Constituency")
Thomas Lumley Assoc. Professor, Biostatistics
tlum...@u.washington.edu University of Washington, Seattle
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
