Michael Friendly
2010-Sep-19 17:20 UTC
[R] odds ratios for n-way tables: seeking an *apply-able method
I'm looking for a way to generalize the calculation of (log) odds
ratios for 2 x 2 x {strata} frequency
tables in vcd::oddsratio() to R x C x {strata} tables.
For an R x C table, F, odds ratios can be defined
in several ways; one simple way, particularly for tables with ordered
factors, is to calculate the
set of continuation odds ratios, based on the (R-1)x(C-1) set of 2x2
tables with elements
F[i:(i+1), j:(j+1)]
What I think I want is a method to extract all such 2x2 subtables to a
structure (list?) so that I can use
*apply or plyr methods to do the computations of (log) odds ratio (& std
error) for each, and return these
in a suitable form.
Below, I define a function OR2x2 that calculates the (log) odds ratio
for one 2x2 table, and another, OR,
that uses loops to calculate what I want, at least for a two-way table.
But the programming would be
quite messy to generalize to n-way tables where all but the first two
dimensions are considered {strata}.
Can someone help me convert this to a form using *apply methods?
# (log) odds ratio for 2x2 subtable
OR2x2 <- function(f, log=TRUE) {
if (any(f == 0)) f <- f + 0.5
f <- log(f)
lor <- f[1, 1] + f[2, 2] - f[1, 2] - f[2, 1]
if (log) lor else exp(lor)
}
# continuation (log) odds ratios for two-way table
OR <- function(f, log=TRUE) {
d <- dim(f)
R <- d[1]; C <- d[2]
res <- matrix(0, R-1, C-1)
rownames(res) <- rep("", R-1)
colnames(res) <- rep("", C-1)
for (i in 1:(R-1)) {
rownames(res)[i] <- paste(rownames(f)[i:(i+1)],
collapse=":")
for (j in 1:(C-1)) {
tab <- f[i:(i+1), j:(j+1)]
res[i,j] <- OR2x2(tab, log=log)
}
}
for (j in 1:(C-1)) {
colnames(res)[j] <- paste(colnames(f)[j:(j+1)],
collapse=":")
}
res
}
### Examples
> (HE <- apply(HairEyeColor, 1:2, sum))
Eye
Hair Brown Blue Hazel Green
Black 68 20 15 5
Brown 119 84 54 29
Red 26 17 14 14
Blond 7 94 10 16
>
> OR(HE)
Brown:Blue Blue:Hazel Hazel:Green
Black:Brown 0.87547 -0.1542 0.4769
Brown:Red -0.07658 0.2477 0.6217
Red:Blond 3.02227 -2.0466 0.4700
> # by strata
> OR(HE1 <- HairEyeColor[,,1])
Brown:Blue Blue:Hazel Hazel:Green
Black:Brown 1.00957 -0.5978 0.6931
Brown:Red 0.05827 0.3365 0.5108
Red:Blond 2.30259 -1.4351 0.4700
> OR(HE2 <- HairEyeColor[,,2])
Brown:Blue Blue:Hazel Hazel:Green
Black:Brown 0.7230 0.4287 0.1881
Brown:Red -0.1634 0.1591 0.7282
Red:Blond 3.5993 -2.5494 0.4700
>
--
Michael Friendly Email: friendly AT yorku DOT ca
Professor, Psychology Dept.
York University Voice: 416 736-5115 x66249 Fax: 416 736-5814
4700 Keele Street Web:http://www.datavis.ca
Toronto, ONT M3J 1P3 CANADA
Dennis Murphy
2010-Sep-19 22:08 UTC
[R] odds ratios for n-way tables: seeking an *apply-able method
Hi Michael:
I'm not an expert on this topic, but I tried a few things and got something
that *looks* right...
The idea is to convert the 3D table to a data frame, split on Sex, use
lapply to reconstruct the individual tables and then use lapply again to get
the results. I'm sure there's a better way to get from table to list of
tables, but simply making a list out of the original table doesn't work.
Begin the circuitous route:
d <- as.data.frame.table(HairEyeColor)
l <- lapply(split(d, d$Sex), function(df) with(df, xtabs(Freq ~ Hair +
Eye)))
lapply(l2, OR)
$Male
Brown:Blue Blue:Hazel Hazel:Green
Black:Brown 1.00957172 -0.5978370 0.6931472
Brown:Red 0.05826891 0.3364722 0.5108256
Red:Blond 2.30258509 -1.4350845 0.4700036
$Female
Brown:Blue Blue:Hazel Hazel:Green
Black:Brown 0.7230001 0.4287220 0.1880522
Brown:Red -0.1633844 0.1590647 0.7282385
Red:Blond 3.5992673 -2.5494452 0.4700036
A plyr approach taken from the data frame d:
g <- function(df) with(df, xtabs(Freq ~ Hair + Eye))
l2 <- dlply(d, 'Sex', g)
llply(l2, OR)
(Same result as above.)
HTH,
Dennis
On Sun, Sep 19, 2010 at 10:20 AM, Michael Friendly
<friendly@yorku.ca>wrote:
> I'm looking for a way to generalize the calculation of (log) odds
ratios
> for 2 x 2 x {strata} frequency
> tables in vcd::oddsratio() to R x C x {strata} tables.
>
> For an R x C table, F, odds ratios can be defined
> in several ways; one simple way, particularly for tables with ordered
> factors, is to calculate the
> set of continuation odds ratios, based on the (R-1)x(C-1) set of 2x2 tables
> with elements
>
> F[i:(i+1), j:(j+1)]
>
> What I think I want is a method to extract all such 2x2 subtables to a
> structure (list?) so that I can use
> *apply or plyr methods to do the computations of (log) odds ratio (&
std
> error) for each, and return these
> in a suitable form.
>
> Below, I define a function OR2x2 that calculates the (log) odds ratio for
> one 2x2 table, and another, OR,
> that uses loops to calculate what I want, at least for a two-way table. But
> the programming would be
> quite messy to generalize to n-way tables where all but the first two
> dimensions are considered {strata}.
> Can someone help me convert this to a form using *apply methods?
>
> # (log) odds ratio for 2x2 subtable
> OR2x2 <- function(f, log=TRUE) {
> if (any(f == 0)) f <- f + 0.5
> f <- log(f)
> lor <- f[1, 1] + f[2, 2] - f[1, 2] - f[2, 1]
> if (log) lor else exp(lor)
> }
>
> # continuation (log) odds ratios for two-way table
> OR <- function(f, log=TRUE) {
> d <- dim(f)
> R <- d[1]; C <- d[2]
> res <- matrix(0, R-1, C-1)
> rownames(res) <- rep("", R-1)
> colnames(res) <- rep("", C-1)
> for (i in 1:(R-1)) {
> rownames(res)[i] <- paste(rownames(f)[i:(i+1)],
collapse=":")
> for (j in 1:(C-1)) {
> tab <- f[i:(i+1), j:(j+1)]
> res[i,j] <- OR2x2(tab, log=log)
> }
> }
> for (j in 1:(C-1)) {
> colnames(res)[j] <- paste(colnames(f)[j:(j+1)],
collapse=":")
> }
> res
> }
>
> ### Examples
>
> > (HE <- apply(HairEyeColor, 1:2, sum))
> Eye
> Hair Brown Blue Hazel Green
> Black 68 20 15 5
> Brown 119 84 54 29
> Red 26 17 14 14
> Blond 7 94 10 16
> >
> > OR(HE)
> Brown:Blue Blue:Hazel Hazel:Green
> Black:Brown 0.87547 -0.1542 0.4769
> Brown:Red -0.07658 0.2477 0.6217
> Red:Blond 3.02227 -2.0466 0.4700
> > # by strata
> > OR(HE1 <- HairEyeColor[,,1])
> Brown:Blue Blue:Hazel Hazel:Green
> Black:Brown 1.00957 -0.5978 0.6931
> Brown:Red 0.05827 0.3365 0.5108
> Red:Blond 2.30259 -1.4351 0.4700
> > OR(HE2 <- HairEyeColor[,,2])
> Brown:Blue Blue:Hazel Hazel:Green
> Black:Brown 0.7230 0.4287 0.1881
> Brown:Red -0.1634 0.1591 0.7282
> Red:Blond 3.5993 -2.5494 0.4700
> >
>
> --
> Michael Friendly Email: friendly AT yorku DOT ca
> Professor, Psychology Dept.
> York University Voice: 416 736-5115 x66249 Fax: 416 736-5814
> 4700 Keele Street Web:http://www.datavis.ca
> Toronto, ONT M3J 1P3 CANADA
>
> ______________________________________________
> 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]]