Verification CA Test Implementation

library(drcHelper)
#> Loading required package: drc
#> Loading required package: MASS
#> Loading required package: drcData
#> 
#> 'drc' has been loaded.
#> Please cite R and 'drc' if used for a publication,
#> for references type 'citation()' and 'citation('drc')'.
#> 
#> Attaching package: 'drc'
#> The following objects are masked from 'package:stats':
#> 
#>     gaussian, getInitial
source("../knitr-setup.R")

Comparison: Cochran-Armitage Test Implementation vs prop.trend.test

The implementation in drcHelper and prop.trend.test are mathematically equivalent but use different computational approaches. In particular, the drcHelper implementation uses actual dose values as default scoring, adds overdispersion handling by Rao-Scott correction for clustered data, the grouping by doses are done automatically, and uses the more intuitive Z-statistic. It is more comprehensive and has a cleaner interface for dose-response studies.

Here’s the detailed comparison:

Key Differences

Aspect	drcHelper Implementation	prop.trend.test
Computational Method	Direct Z-statistic calculation	Weighted linear regression
Test Statistic	Z-score (normal distribution)	Chi-squared (\(Z^2\))
Score Usage	Uses actual doses values	Uses seq_along(x) by default
Overdispersion	Rao-Scott correction available	No correction
Data Grouping	Groups by doses automatically	Assumes pre-grouped data

Mathematical Relationship

The relationship between the two approaches:

• Z-statistic: \(Z = \frac{\sum n_i(s_i - \bar{s})p_i}{\sqrt{p(1-p)\sum n_i(s_i - \bar{s})^2}}\)
• prop.trend.test Chi-squared: \(\chi^2 = Z^2\)

Where \(p(\chi^2) = 2 \times p(|Z|)\) for two-sided tests.

Verification Example ✅

# Test data
successes <- c(83, 90, 129, 70)
totals <- c(86, 93, 136, 82)

# R's function (default scores: 1,2,3,4)
result_r <- prop.trend.test(successes, totals)
result_r
#> 
#>  Chi-squared Test for Trend in Proportions
#> 
#> data:  successes out of totals ,
#>  using scores: 1 2 3 4
#> X-squared = 8.2249, df = 1, p-value = 0.004132
# drcHelper function (after correction, using doses = 1,2,3,4)
result <- cochranArmitageTrendTest(successes, totals, doses = 1:4)
result
#> 
#>  Cochran-Armitage test for trend (doses scoring)
#> 
#> data:  proportions 0.965, 0.968, 0.949, 0.854 at doses 1, 2, 3, 4
#> = -2.8679, p-value = 0.004132
#> alternative hypothesis: two.sided
# Should give: sqrt(result_r$statistic) ≈ abs(result_yours$statistic)
sqrt(result_r$statistic) - abs(result$statistic)
#>    X-squared 
#> 6.217249e-15