Hello,
Thanks for the indications for use epi.2by2 in svytable, however when running
the code throws warnings
that reference NANs, could help me with this? Following you show the results:
> int19<-svytable(~Perdida_peso_no_intensionada+APES_DICOT+tabaco,Muestra.comp,Ntotal
> = NULL)
> int19
, , tabaco = Nunca fumador
APES_DICOT
Perdida_peso_no_intensionada Buena Mala
No/No sabe 2108.3 599.8
Si 264.8 253.6
, , tabaco = Activo pasivo
APES_DICOT
Perdida_peso_no_intensionada Buena Mala
No/No sabe 200.8 32.0
Si 25.0 32.0
> epi.2by2(dat = int19, method = "case.control", conf.level = 0.95,
+ units = 100, homogeneity = "breslow.day", verbose = TRUE)$OR.homog
test.statistic df p.value
1 6.566923 1 0.01038914
Hubo 12 avisos (use warnings() para verlos)
> warnings(int19)
Warning messages:
1: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
2: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
3: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
4: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
5: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
6: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
7: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
8: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
9: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
10: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
11: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
12: In qbinom(p, size, prob, lower.tail, log.p) : Se han producido NaNs 2108.3
264.8 599.8 253.6 200.8 25 32 32
Thanks
> Date: Sat, 1 Sep 2012 12:34:11 +1200
> Subject: Re: [R] test Breslow-Day for svytable??
> From: tlum...@uw.edu
> To: dwinsem...@comcast.net
> CC: dianamm...@hotmail.com; r-help@r-project.org; m.steven...@massey.ac.nz
>
> On Sat, Sep 1, 2012 at 4:27 AM, David Winsemius <dwinsem...@comcast.net>
> wrote:
> >
> > On Aug 31, 2012, at 7:20 AM, Diana Marcela Martinez Ruiz wrote:
> >
> >> Hi all,
> >>
> >> I want to know how to perform the test Breslow-Day test for homogeneity of
> >> odds ratios (OR) stratified for svytable. This test is obtained with the
> >> following code:
> >>
> >> epi.2by2 (dat = daty, method = "case.control" conf.level = 0.95,
> >
> > missing comma here ...........................^
> >
> >> units = 100, homogeneity = "breslow.day", verbose = TRUE)
> >>
> >> where "daty" is the object type table svytable consider it, but when I
> >> run the code
> >> does not throw the homogeneity test.
> >
> > You are asked in the Posting guide to copy all errors and warnings when
> > asking about unexpected behavior. When I run epi.2y2 on the output of a
> > syvtable object I get no errors, but I do get warnings which I think are
> > due to non-integer entries in the weighted table. I also get from a
> > svytable() usingits first example on the help page an object that is NOT a
> > set of 2 x 2 tables in an array of the structure as expected by epi.2by2().
> > The fact that epi.2by2() will report numbers with labels for a 2 x 3 table
> > means that its error checking is weak.
> >
> > This is the output of str(dat) from one of the example on epi.2by2's help
> > page:
> >
> >> str(dat)
> > table [1:2, 1:2, 1:3] 41 13 6 53 66 37 25 83 23 37 ...
> > - attr(*, "dimnames")=List of 3
> > ..$ Exposure: chr [1:2] "+" "-"
> > ..$ Disease : chr [1:2] "+" "-"
> > ..$ Strata : chr [1:3] "20-29 yrs" "30-39 yrs" "40+ yrs"
> >
> > Notice that is is a 2 x 2 x n array. (Caveat:: from here on out I am simply
> > reading the help pages and using str() to look at the objects created to
> > get an idea regarding success or failure. I am not an experienced user of
> > either package.) I doubt that what you got from svytable is a 2 x 2
> > table. As another example you can build a 2 x 2 x n table from the built-in
> > dataset: UCBAdmissions
> >
> > DF <- as.data.frame(UCBAdmissions)
> > ## Now 'DF' is a data frame with a grid of the factors and the counts
> > ## in variable 'Freq'.
> > dat2 <- xtabs(Freq ~ Gender + Admit+Dept, DF)
> > epiR::epi.2by2(dat = dat2, method = "case.control", conf.level = 0.95,
> > units = 100, homogeneity = "breslow.day", verbose = TRUE)$OR.homog
> > #-------------
> > test.statistic df p.value
> > 1 18.82551 5 0.00207139
> >
> > Using svydesign and svytable I _think_ this is how one would go about
> > constructing a 2 x 2 table:
> >
> > tbl2<-svydesign( ~ Gender + Admit+Dept, weights=~Freq, data=DF)
> > summary(dclus1)
> > (tbl2by2 <- svytable(~ Gender + Admit+Dept, tbl2))
> > epiR::epi.2by2(dat = tbl, method = "case.control", conf.level = 0.95,
> > units = 100, homogeneity = "breslow.day", verbose = TRUE)$OR.homog
> > #-----------------------
> > test.statistic df p.value
> > 1 18.82551 5 0.00207139
> >
> > (At least I got internal consistency. I see you copied Thomas Lumley, which
> > is a good idea. I'll be happy to get corrected on any point. I'm adding the
> > maintainer of epiR to the recipients.)
> >
>
> Yes, that will give internal consistency from a data structure point
> of view. It won't give a valid test in real examples, though --
> epi.2by2 doesn't know about complex sampling, and what you're passing
> it is just an estimate of the population 2x2xK table.
>
> What would work, though it's not quite the same as the Breslow-Day
> test, is to use svyloglin() and do a Rao-Scott test comparing the
> model with all two-way interactions ~(Gender+Dept+Admit)^2 to the
> saturated model ~Gender*Dept*Admit.
>
> -thomas
>
>
> --
> Thomas Lumley
> Professor of Biostatistics
> University of Auckland
[[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.
