Package 'ibist'

Observer Ratings for Aneurysmal Subarachnoid Hemorrhage Scales

Description

Reconstructed paired observer ratings for three clinical grading scales used after aneurysmal subarachnoid hemorrhage: WFNS, Hunt and Hess (HH), and PAASH.

Format

A data frame with 309 rows (103 paired ratings for each of 3 scales) and 4 variables:

scale

Grading scale: WFNS, HH, or PAASH.

assessment

Sequential row number within scale. These numbers are generated for the reconstructed data and are not study participant IDs.

observer1

Ordered factor giving the first observer's grade (1–5).

observer2

Ordered factor giving the second observer's grade (1–5).

Details

The paired ratings are expanded from the published cross-tabulations. They reproduce the observed cell counts for each scale and can be used to illustrate agreement summaries such as kappa. The original participant-level records are not available here; generated rows cannot be linked across scales.

Source

Degen, L. A. R., Dorhout Mees, S. M., Algra, A., and Rinkel, G. J. E. (2011). Interobserver variability of grading scales for aneurysmal subarachnoid hemorrhage. Stroke, 42(6), 1546–1549. doi:10.1161/STROKEAHA.110.601211

Examples

data(asah)
wfns <- subset(asah, scale == "WFNS")
table(wfns$observer1, wfns$observer2)
# irr::kappa2(wfns[c("observer1", "observer2")], weight = "squared")

---

Cochran-Armitage Trend Test

Description

Tests for a linear trend in a binary outcome across ordered groups using the Cochran-Armitage statistic. The test can use a normal approximation or an exact conditional distribution.

Usage

ca.trend.test(
  x,
  scores = NULL,
  alternative = c("two.sided", "greater", "less"),
  exact = FALSE
)

Arguments

x	An $R \times 2$ matrix or table of non-negative integer counts. Rows are ordered groups and columns are outcome categories. The first column is treated as the event of interest.
scores	Optional numeric scores for the ordered groups. Defaults to 1:R; scores must be finite and strictly increasing.
alternative	Alternative hypothesis: "two.sided", "greater" (increasing trend), or "less" (decreasing trend).
exact	Logical; if TRUE, compute the exact conditional p-value given both margins. Otherwise, use the asymptotic normal approximation.

Details

A positive signed statistic indicates that the event proportion tends to increase with the group scores. The exact test conditions on the row totals and the total number of events and uses the conditional distribution of $T = \sum_i x_i n_{i1}$. For a two-sided exact test, the p-value is twice the smaller one-sided tail probability, truncated at one. The exact calculation is performed by dynamic programming and can become slow when the table has many rows or large margins.

Value

An object of class "htest". The signed Cochran-Armitage statistic is returned along with its square, the score statistic $T$, the p-value, and the alternative hypothesis.

References

Cochran, W. G. (1954). Some methods for strengthening the common chi-squared tests. Biometrics, 10(4), 417–451. doi:10.2307/3001616

Armitage, P. (1955). Tests for linear trends in proportions and frequencies. Biometrics, 11(3), 375–386. doi:10.2307/3001775

Portier, C., and Hoel, D. G. (1984). Type 1 error of trend tests in proportions and the design of cancer screens. Communications in Statistics - Theory and Methods, 13(1), 1–14.

Examples

tab <- matrix(c(106, 109, 128, 130, 99, 95, 76, 74), nrow = 4)
ca.trend.test(tab)
ca.trend.test(tab, scores = c(1, 2, 4, 8), alternative = "greater")
small.tab <- matrix(c(2, 0, 1, 1, 0, 2), nrow = 3, byrow = TRUE)
ca.trend.test(small.tab, exact = TRUE)

---

CMH Linear-by-Linear Association Test

Description

Tests for a zero score correlation in an RxC contingency table using the Cochran-Mantel-Haenszel linear-by-linear statistic.

Usage

cmh.ll.test(
  x,
  measure = c("pearson", "spearman"),
  alternative = c("two.sided", "greater", "less"),
  row.scores = NULL,
  col.scores = NULL
)

Arguments

x	A two-dimensional matrix or table of non-negative integer counts. Rows and columns must be ordered from the lowest to highest category.
measure	The score correlation to test: "pearson" or "spearman". Pearson correlation uses the supplied category scores or, by default, numeric category labels and otherwise category order numbers. Spearman correlation uses marginal midrank scores.
alternative	Alternative hypothesis: "two.sided", "greater", or "less". A positive association means that larger row scores tend to occur with larger column scores.
row.scores, col.scores	Optional numeric scores for row and column categories when measure = "pearson". These must be specified together. They are not used for Spearman correlation.

Details

Let $r$ be the Pearson correlation between the assigned row and column scores (or their marginal midranks for Spearman correlation). The signed test statistic is $Z = r\sqrt{n-1}$, and the CMH statistic is $X^2 = Z^2$, which has an asymptotic chi-squared distribution with one degree of freedom under the null hypothesis. One-sided tests use the signed normal statistic. This function tests zero association only; it does not return a confidence interval. Use ordinal.assoc() for the multinomial-based Wald confidence interval for Pearson or Spearman correlation. Those intervals are not obtained by inverting this test. The Pearson and Spearman tests in ordinal.assoc() use their own multinomial null variances and are distinct from this CMH test; their p-values can differ.

Value

An object of class "htest". The statistic element contains the signed $Z$ statistic, and the chisq element contains its square, the CMH statistic. The estimate element contains the score correlation.

References

Mantel, N. (1963). Chi-square tests with one degree of freedom: extensions of the Mantel-Haenszel procedure. Journal of the American Statistical Association, 58(303), 690–700. doi:10.1080/01621459.1963.10500879

Examples

tab <- matrix(c(109, 68, 14, 72, 79, 37, 9, 42, 139), nrow = 3,
              byrow = TRUE)
cmh.ll.test(tab, measure = "pearson")
cmh.ll.test(tab, measure = "spearman", alternative = "greater")
ordinal.assoc(tab, measure = "spearman")

---

Asymptotic Inference for Cohen's Kappa

Description

Computes Cohen's kappa or weighted kappa from a square table of paired ratings, with a Wald confidence interval and an asymptotic test of $H_0: \kappa = 0$.

Usage

cohen.kappa.test(
  x,
  weights = "unweighted",
  conf.level = 0.95,
  alternative = c("two.sided", "greater", "less"),
  test.variance = c("fce", "delta")
)

Arguments

x	A square matrix or table of non-negative integer counts. Rows and columns represent the same ordered rating categories.
weights	The agreement weights: "unweighted", "linear", "quadratic", or a numeric matrix of user-specified weights. The named choices use category positions in x; the quadratic weights are Fleiss–Cohen weights. A supplied matrix must have the same dimensions as x, be symmetric, have ones on the diagonal, and contain values between zero and one.
conf.level	Confidence level for the two-sided Wald confidence interval.
alternative	Alternative hypothesis: "two.sided", "greater", or "less".
test.variance	Variance used for the test: "fce" for the Fleiss–Cohen–Everitt null variance, or "delta" for the multinomial delta-method variance evaluated at the observed proportions.

Details

The confidence interval uses the Fleiss–Cohen–Everitt (FCE) large-sample variance of the kappa estimate. This variance is algebraically equivalent to the multinomial delta-method variance used by psych::cohen.kappa(). By default, the test uses the distinct FCE variance under the null hypothesis $\kappa = 0$. Setting test.variance = "delta" instead uses the delta-method variance evaluated at the observed cell proportions. Evaluating the delta-method variance under the null gives the FCE null variance, so it is not a separate test option. The interval is an untruncated Wald interval.

This function provides asymptotic inference only. SAS PROC FREQ also provides exact tests for simple and weighted kappa; exact inference is not implemented here.

Value

An object of class "htest" containing the estimate, confidence interval, Z statistic, and p-value. It also includes the standard error used for the confidence interval, the standard error under the null hypothesis, the selected test variance method, the weight matrix, and the sample size.

References

Fleiss, J. L., Cohen, J., and Everitt, B. S. (1969). Large sample standard errors of kappa and weighted kappa. Psychological Bulletin, 72(5), 323–327. doi:10.1037/h0028106

SAS Institute Inc. (2015). Tests and measures of agreement. In Base SAS 9.4 Procedures Guide: Statistical Procedures, Third Edition. https://support.sas.com/documentation/cdl/en/procstat/67528/HTML/default/procstat_freq_details76.htm

Examples

peff <- matrix(
  c(581, 134, 12, 6, 0, 2,
    91, 73, 18, 31, 3, 1,
    4, 18, 15, 19, 3, 3,
    1, 8, 7, 24, 11, 6,
    0, 0, 2, 4, 3, 5,
    0, 0, 0, 9, 11, 13),
  nrow = 6,
  byrow = TRUE
)
cohen.kappa.test(peff)
cohen.kappa.test(peff, weights = "quadratic")
cohen.kappa.test(peff, weights = "quadratic", test.variance = "delta")

---

Demonstrate the Central Limit Theorem

Description

The demo_clt() function generates plots to illustrate the Central Limit Theorem (CLT) using a specified random number generator. The function displays standardized sampling distributions for different sample sizes and overlays the standard normal density.

Usage

demo_clt(rng, n, nrep = 10000, ..., pmean = NULL, psd = NULL)

Arguments

rng	A random number generator function taking the sample size as its first argument (e.g., runif, rnorm, rgamma).
n	A numeric vector of sample sizes (e.g., c(5, 10, 20, 40)).
nrep	The number of repetitions for generating sample means (default is 10000).
...	Additional arguments passed to the random number generator (e.g., shape and rate for rgamma).
pmean	The population mean of the distribution. If NULL, it is estimated from a large Monte Carlo sample.
psd	The population standard deviation of the distribution. If NULL, it is estimated from a large Monte Carlo sample.

Value

A ggplot2 object showing the standardized sampling distributions for different sample sizes, compared against the standard normal curve.

Examples

set.seed(123)
demo_clt(runif, n = c(5, 20), nrep = 100, min = 0, max = 1,
         pmean = 0.5, psd = sqrt(1 / 12))

demo_clt(rgamma, n = c(5, 20), shape = 2, rate = 1,
         nrep = 100, pmean = 2, psd = sqrt(2)
)

---

Depression, insomnia, and stress dataset (Kalmbach et al., 2018)

Description

This dataset is derived from a longitudinal study examining the role of insomnia and cognitive intrusions in predicting future depression. It includes participants who were not depressed at baseline and who met inclusion criteria on age, missingness, and stress exposure.

Usage

depress

Format

A data frame with 1730 observations and 9 variables:

age

Age in years.

gender

Gender of the participant (coded per original survey).

stressorA

Indicator for experiencing at least one stressor in the past year.

depress1yr

Depression status at 1-year follow-up.

depress2yr

Depression status at 2-year follow-up.

depress1or2yr

Depression at either 1 or 2 years.

IEStert

Tertile of cognitive intrusion score based on the Impact of Event Scale (IES).

IESgrp

Grouped cognitive intrusion: 0 = lowest tertile, 1 = middle tertile, 2 = highest tertile.

insom

Insomnia indicator defined as: sleep onset latency (SOL) > 30 minutes or wake after sleep onset (WASO) > 30 minutes.

Details

The dataset is constructed by applying the following criteria:

• Baseline depression score < 10.

• Age >= 18.

• Non-missing gender, cognitive intrusion (IES), and WASO.

• At least one stressor in the past year.

• Exclusion of extreme nocturnal wakefulness (> 90 minutes).

The resulting dataset (n = 1730) includes all three cognitive intrusion tertiles. In the original paper, the analytic sample was further restricted by removing observations with missing IEStert, which effectively excludes the middle tertile (coded as IESgrp = 1), yielding a smaller sample (n = 1126). This optional restriction is not applied here but can be reproduced by subsetting:

depress_paper <- subset(depress, !is.na(IEStert))

Variable definitions follow the original study documentation provided in the supplementary spreadsheet (S1 Table).

Source

Kalmbach, D. A., Pillai, V., and Drake, C. L. (2018). Nocturnal insomnia symptoms and stress-induced cognitive intrusions in risk for depression: a 2-year prospective study. PLOS ONE, 13(2), e0192088.

References

Kalmbach, D. A., Pillai, V., and Drake, C. L. (2018). Nocturnal insomnia symptoms and stress-induced cognitive intrusions in risk for depression: a 2-year prospective study. PLOS ONE, 13(2), e0192088.

Examples

data(depress)
str(depress)

# proportion developing depression within 2 years
mean(depress$depress1or2yr)

# insomnia prevalence
mean(depress$insom)

# reproduce analytic sample used in the paper
depress_paper <- subset(depress, !is.na(IEStert))

---

Dental Fluorosis Study Data

Description

These data arise from a prospective study examining the association between dental fluorosis and amoxicillin use during early childhood.

Usage

fluorosis

Format

A data frame with 24 rows and 5 variables:

om

Oral hygiene measure (factor with levels 0, 1)

fl_lev

Fluoride level (factor with levels 1, 2, 3)

amox

Amoxicillin use during early childhood (factor: 0 = no, 1 = yes)

fluorosis

Dental fluorosis status (factor: 0 = absent, 1 = present)

count

Number of subjects in the corresponding cell

Details

Dental fluorosis — a developmental defect of tooth enamel — was assessed as a binary outcome (presence or absence) when participants were approximately 9 years old. Information on fluoride intake, infections, and amoxicillin use during early childhood (birth to 32 months) was collected as part of the Iowa Fluoride Study.

The dataset is provided in aggregated form as a four-way contingency table with cell counts.

The data have a 2 x 3 x 2 x 2 factorial structure: oral hygiene (2 levels), fluoride level (3 levels), amoxicillin use (2 levels), and fluorosis status (2 levels). The counts sum to 308 subjects.

Source

Hong L, Levy SM, Warren JJ, Dawson DV, Bergus GR, Wefel JS (2005). Association of amoxicillin use during early childhood with developmental tooth enamel defects. Archives of Pediatrics & Adolescent Medicine, 159(10), 943–948.

Examples

data(fluorosis)
with(fluorosis,
     xtabs(count ~ om + fl_lev + amox + fluorosis))

---

Pearson Chi-squared Goodness-of-Fit Test for Logistic Regression

Description

Performs the Pearson chi-squared goodness-of-fit (GoF) test for a logistic regression model fitted via glm(..., family = "binomial"). The test assesses calibration by comparing observed and expected counts across covariate patterns.

Usage

logistic.gof(fit, min_n = 5, min_expected = 5, pool = TRUE)

logistic_gof(...)

Arguments

fit	A fitted glm object with family = "binomial".
min_n	Minimum total count per group for pooling. Default is 5.
min_expected	Minimum expected count per group for pooling. Default is 5.
pool	Logical; if TRUE, adjacent groups (ordered by fitted probabilities) are pooled to meet minimum thresholds. Default is TRUE.
...	Arguments passed to logistic.gof() by the deprecated logistic_gof() alias.

Details

The Pearson GoF statistic is defined as

$T_P = \sum_{\ell=1}^J \frac{(O_\ell - E_\ell)^2}{V_\ell},$

where $O_\ell$ is the observed number of events, $E_\ell = n_\ell \hat{\pi}_\ell$ is the expected number of events, and $V_\ell = n_\ell \hat{\pi}_\ell (1 - \hat{\pi}_\ell)$ is the variance for the $\ell$-th covariate pattern. Under the null hypothesis of correct model specification, $T_P$ approximately follows a chi-squared distribution with $J - (k + 1)$ degrees of freedom, where $J$ is the number of groups and $k + 1$ is the number of model parameters.

This test requires replicated covariate patterns or grouped binomial data. When groups are small, optional pooling can be applied to ensure stable chi-squared approximation.

The function supports three types of binomial responses:

• Binary outcomes (0/1)

• Proportions with weights (grouped binomial)

• Two-column matrix response via cbind(success, failure)

When no replication exists (i.e., each observation has a unique covariate pattern), the Pearson GoF test is not appropriate. In such cases, consider alternative methods such as the Hosmer–Lemeshow test.

Pooling stabilizes the test by merging adjacent groups with small counts, trading a small bias for reduced variance in the chi-squared approximation.

Value

An object of class "ibist_logistic_gof" and "htest" with components:

statistic	The Pearson chi-squared and deviance statistics.
parameter	Degrees of freedom.
p.value	The corresponding p-values.
method	Description of the test.
data.name	Description of the data.

Additional components (not printed) include:

observed	Observed counts by group.
expected	Expected counts by group.
group_n	Group sizes.
groups	Data frame of grouped summaries.

See Also

glm, chisq.test

Examples

## Example: Surfactant use and birthweight (RDS data)
data(rds, package = "ibist")

fit <- glm(death ~ surf + bwt,
  data = rds, weights = count,
  family = "binomial"
)

## Pearson GoF test with pooling (default)
(res <- logistic.gof(fit))

## Without pooling (may be unstable if small groups exist)
logistic.gof(fit, pool = FALSE)

## Inspect grouped diagnostics
res$groups

---

Check the Mantel-Haenszel Asymptotic Condition

Description

Checks whether the large-sample condition for the Mantel-Haenszel test is satisfied for a collection of stratified 2 x 2 tables.

Usage

mantelhaen.asymcheck(x)

Arguments

x	A list of 2 x 2 tables, or a 2 x 2 x K array of non-negative integer counts. Rows are the two groups and columns are the two outcome categories.

Details

Conditional on the margins of stratum $s$, the upper-left cell count $N_{s11}$ has hypergeometric variance

$V(N_{s11}) = \frac{n_{s1.}n_{s2.}n_{s.1}n_{s.2}}{n_s^2(n_s-1)}.$

The chi-squared approximation used by the Mantel-Haenszel test is considered adequate when

$\sum_s V(N_{s11}) \ge 5.$

This diagnostic applies to the asymptotic test. It is not required for an exact conditional test such as or.common.test.

Value

An object of class "ibist_mantelhaen_asymcheck" containing:

variance	the hypergeometric variance for each stratum.
variance.sum	the sum of the stratum-specific variances.
condition.met	a logical value indicating whether the sum is at least 5.
data.name	a character string describing the data.

See Also

mantelhaen.test, or.common.test

Examples

strat <- array(
  c(88, 469, 48, 653, 27, 480, 34, 1324),
  dim = c(2, 2, 2)
)
mantelhaen.asymcheck(strat)

---

Blood Pressure Cohort Data (midbp)

Description

Aggregated data of cases and person-years by diastolic blood pressure category and gender.

Usage

midbp

Format

A data frame with 16 rows and 4 variables:

dbp

Diastolic blood pressure category (1–8).

gender

Gender (Men, Women).

cases

Number of observed cases.

pyears

Person-years of follow-up.

Details

The dataset is aggregated by DBP category and gender, suitable for rate modeling (e.g., Poisson regression with offset).

The DBP categories are coded as: 1 = <70 mmHg, 2 = 70–74 mmHg, 3 = 75–79 mmHg, 4 = 80–84 mmHg, 5 = 85–89 mmHg, 6 = 90–94 mmHg, 7 = 95–99 mmHg, and 8 = >=100 mmHg.

Examples

data(midbp)
head(midbp)

---

HDL Cholesterol and Event Counts by Gender

Description

Aggregated data on event counts and person-years by HDL cholesterol category and gender. Suitable for Poisson rate modeling and incidence rate comparisons.

Usage

mihdl

Format

A data frame with 24 rows and 4 variables:

hdl

HDL category (integer, 1–12)

gender

Gender (factor: Men, Women)

cases

Number of events

pyears

Person-years at risk

Details

The dataset provides grouped count data for modeling rates using Poisson regression with an offset for log person-years.

Examples

data(mihdl)
fit <- glm(cases ~ factor(hdl) + gender,
           offset = log(pyears),
           family = poisson,
           data = mihdl)
summary(fit)

---

Non-restorative sleep and physical activity (Japan cohort study)

Description

A large observational dataset from a cohort study conducted in Japan to examine the association between non-restorative sleep (NRS) and physical activity, gender, and age. The data are used to illustrate logistic regression modeling for a binary outcome in a large-sample setting.

Usage

nrs

Format

A data frame with 91,795 subjects (90,122 complete records) on the following variables:#'

id

Subject identifier.

gender

Gender of the subject (F = Female; M = Male).

age

Age in years in 2013.

ex

Indicator of regular exercise in 2013 (integer-coded).

pa

Physical activity measure in 2013 (integer-coded).

nrs

Indicator of non-restorative sleep in 2013 (1 = presence; 0 = absence).

Details

Non-restorative sleep (NRS) is defined as a subjective feeling of lack of refreshment on awakening and reflects qualitative aspects of sleep. Hidaka et al. (2019) analyzed these data using logistic regression to assess whether the probability of NRS is associated with physical activity, gender, and age in a large cohort of adult subjects in Japan. Within this package, the dataset is provided for methodological illustration of binary regression models rather than for substantive epidemiological inference.

Missing values are represented as NA, which are present in 'ex', 'pa', and 'nrs'.

Source

Hidaka et al. (2019).

Examples

data(nrs)
summary(nrs)

---

Confidence Intervals for a 2 x 2 Odds Ratio

Description

Computes confidence intervals for the odds ratio in a 2 x 2 table.

Usage

or.2x2.ci(
  x,
  conf.level = 0.95,
  method = c("wald", "gart", "cornfield", "bp"),
  midp = FALSE,
  ...
)

Arguments

x	A 2 x 2 table of non-negative integer counts.
conf.level	Confidence level for the interval.
method	Method or methods for confidence interval. One or more of "wald", "gart", "cornfield", or "bp".
midp	Logical; should the conditional exact method be modified using mid-p probabilities? The default is FALSE.
...	Reserved for future extensions.

Details

For a table with entries $n_{11}$, $n_{12}$, $n_{21}$, and $n_{22}$, the sample odds ratio is $n_{11} n_{22} / (n_{12} n_{21})$.

The "wald" method uses the usual large-sample interval on the log odds-ratio scale. The addition of 0.5 to every cell is the Haldane–Anscombe correction, introduced independently by Haldane (1956) and Anscombe (1956). The "gart" method applies this correction using the interval construction in Gart (1966). The "cornfield" method inverts two one-sided conditional tests, while the "bp" method inverts a two-sided conditional test ordered by null probabilities. Both use the noncentral hypergeometric distribution. Setting midp = TRUE subtracts half the observed-table probability from the corresponding exact p-value. Exact intervals are conservative; mid-p intervals are not guaranteed to attain at least the nominal coverage level.

Value

An object of class "ibist_ci" containing the common estimate and a confidence-limit matrix with one row per method.

References

Anscombe, F. J. (1956). On estimating binomial response relations. Biometrika, 43(3–4), 461–464. doi:10.1093/biomet/43.3-4.461

Baptista, J., and Pike, M. C. (1977). Algorithm AS 115: Exact two-sided confidence limits for the odds ratio in a 2 x 2 table. Journal of the Royal Statistical Society. Series C, 26(2), 214–220.

Fagerland, M. W., Lydersen, S., and Laake, P. (2015). Recommended confidence intervals for two independent binomial proportions. Statistical Methods in Medical Research, 24(2), 224–254.

Gart, J. J. (1966). Alternative analyses of contingency tables. Journal of the Royal Statistical Society: Series B, 28(1), 164–179.

Haldane, J. B. S. (1956). The estimation and significance of the logarithm of a ratio of frequencies. Annals of Human Genetics, 20(4), 309–311. doi:10.1111/j.1469-1809.1955.tb01285.x

Examples

tab <- matrix(c(7, 27, 1, 33), nrow = 2, byrow = TRUE)
or.2x2.ci(tab)
or.2x2.ci(tab, method = "gart")
or.2x2.ci(tab, method = "cornfield")
or.2x2.ci(tab, method = "bp", midp = TRUE)

---

Exact Inference for Stratified 2 x 2 Tables

Description

Computes an exact conditional test or confidence interval for a common odds ratio across a collection of 2 x 2 tables.

Usage

or.common.test(
  x,
  theta = 1,
  alternative = c("two.sided", "greater", "less"),
  ...
)

or.common.ci(
  x,
  conf.level = 0.95,
  algorithm = c("network", "convolution"),
  ...
)

or.homogeneity.test(x, ...)

Arguments

x	A list of 2 x 2 tables, or a 2 x 2 x K array of non-negative integer counts. Rows are the two groups and columns are the two outcome categories.
theta	The null common odds ratio.
alternative	A character string specifying the alternative hypothesis: "two.sided", "greater", or "less".
...	Reserved for future extensions.
conf.level	Confidence level for the interval.
algorithm	Computational algorithm for the exact confidence interval. The default, "network", uses dynamic programming. The "convolution" option is retained as an independent check.

Details

Conditional on the margins of every stratum, the upper-left cell in each table has a noncentral hypergeometric distribution. The test combines the strata through the sum of those upper-left cells. The two-sided test uses the Jung et al. (2014) rule $2\min(p_\mathrm{less},p_\mathrm{greater})$, truncated at one. The confidence interval is obtained by inverting these one-sided exact tests.

or.homogeneity.test() performs Zelen's exact test of $H_0: \theta_1 = \cdots = \theta_K$. It conditions on the observed sum of the upper-left cells, which removes the unspecified common odds ratio from the null distribution. The two-sided p-value uses probability ordering.

All three procedures are exact conditional procedures; none is an asymptotic alternative. The first two procedures test or estimate a common odds ratio. The homogeneity procedure tests whether the stratum-specific odds ratios are equal, without specifying their common value.

Value

or.common.test() returns an object of class "htest". or.common.ci() returns an object of class "ibist_ci". or.homogeneity.test() returns an object of class "htest".

References

Jung, S. H., Biggerstaff, B. J., and Kim, C. (2014). A class of exact tests for a common odds ratio in stratified 2 x 2 tables. Biometrical Journal, 56(1), 97–109. doi:10.1002/bimj.201300048

Mehta, C. R., Patel, N. R., and Gray, R. (1985). On computing an exact confidence interval for the common odds ratio in several 2 x 2 contingency tables. Journal of the American Statistical Association, 80(392), 969–973. doi:10.1080/01621459.1985.10478212

Zelen, M. (1971). The analysis of several 2 x 2 contingency tables. Biometrika, 58(1), 129–137. doi:10.1093/biomet/58.1.129

Examples

tabs <- list(
  matrix(c(7, 3, 2, 8), nrow = 2, byrow = TRUE),
  matrix(c(4, 6, 1, 9), nrow = 2, byrow = TRUE)
)
or.common.test(tabs)
or.common.ci(tabs)
or.homogeneity.test(tabs)

---

Confidence Intervals for a Matched-Pair Odds Ratio

Description

Computes confidence intervals for the conditional odds ratio from a 2 x 2 table of matched pairs.

Usage

or.paired.ci(
  x,
  conf.level = 0.95,
  method = c("wald", "wilson", "clopper-pearson", "midp", "blaker")
)

Arguments

x	A 2 x 2 table of non-negative integer pair counts. Rows describe the first member of each pair and columns the second, with the positive category first in both dimensions. Thus x[1, 2] is the number positive only for the first member and x[2, 1] the number positive only for the second. In a matched case-control study, use cases in rows and controls in columns, with exposure present first.
conf.level	Confidence level for the intervals.
method	One or more of "wald" (the default), "wilson", "clopper-pearson", "midp", or "blaker".

Details

Write $b = x_{12}$ and $c = x_{21}$. The matched-pair odds ratio is estimated by $b/c$, not by the cross-product ratio of the pair-count table. Concordant pairs do not contribute information about this conditional odds ratio. Transposing the table or reversing both binary category orders reciprocates the odds ratio. Reversing only one category order changes which pairs are discordant and is not the same analysis.

The "wald" interval is $\exp\{\log(b/c) \mathbin{\pm} z_{1-\alpha/2}\sqrt{1/b+1/c}\}$. It equals the asymptotic Mantel–Haenszel interval when each pair is represented as an individual stratum. It requires both discordant counts to be positive; zero concordant counts are allowed.

Conditional on $b+c$, $B$ is binomial with success probability $p = \theta/(1+\theta)$. The remaining methods transform binomial probability limits $(L,U)$ to $(L/(1-L),U/(1-U))$. The Wilson method inverts the score test, Clopper–Pearson inverts equal-tailed exact tests, and mid-p subtracts half the observed mass from each inclusive tail. Blaker's method takes the smallest interval containing the acceptance set of its exact two-sided test. Its acceptability statistic is the smaller inclusive binomial tail, and its p-value sums probabilities of outcomes with acceptability no larger than the observed value.

Clopper–Pearson and Blaker intervals have at least nominal conditional coverage; Wald, Wilson, and mid-p intervals have no such guarantee. With only one nonzero discordant count the estimate is zero or infinity, and transformed intervals remain available. With no discordant pairs, the odds ratio is unidentified: the estimate is NA, transformed intervals are c(0, Inf), and a warning is issued.

This function estimates a conditional matched-pair odds ratio, not a marginal odds ratio or a difference of paired proportions. For the latter, use prop.paired.ci().

Value

An object of class "ibist_ci" with a matched-pair odds-ratio estimate and one confidence-limit row per requested method.

References

Fagerland, M. W., Lydersen, S., and Laake, P. (2014). Recommended tests and confidence intervals for paired binomial proportions. Statistics in Medicine, 33, 2850–2875. doi:10.1002/sim.6148

Blaker, H. (2000). Confidence curves and improved exact confidence intervals for discrete distributions. Canadian Journal of Statistics, 28, 783–798. doi:10.2307/3315916

Examples

twins <- matrix(c(51, 30, 10, 41), nrow = 2, byrow = TRUE)
or.paired.ci(twins)
or.paired.ci(twins,
             method = c("wilson", "clopper-pearson", "midp", "blaker"))
or.paired.ci(matrix(c(8, 4, 0, 6), 2, byrow = TRUE), method = "wilson")

---

Ordinal Association Measures with Asymptotic Confidence Intervals and Tests

Description

Computes ordinal association measures and their asymptotic confidence intervals from an RxC table. It also tests the null hypothesis that each association measure is zero.

Usage

ordinal.assoc(
  x,
  measure = c(
    "gamma", "tau-b", "tau-c", "somers-c|r", "somers-r|c",
    "pearson", "spearman"
  ),
  conf.level = 0.95,
  alternative = c("two.sided", "greater", "less"),
  row.scores = NULL,
  col.scores = NULL
)

Arguments

x	A two-dimensional matrix or table of non-negative integer counts. Rows and columns must be ordered from the lowest to highest category.
measure	One or more measures: "gamma", "tau-b", "tau-c", "somers-c|r", "somers-r|c", "pearson", and "spearman". By default, all seven are computed.
conf.level	Confidence level for the two-sided asymptotic confidence intervals.
alternative	Alternative hypothesis for the tests: "two.sided", "greater", or "less". Positive association means that higher row categories tend to occur with higher column categories.
row.scores, col.scores	Optional numeric scores for row and column categories when computing Pearson correlation. By default, numeric category labels are used; otherwise, category order numbers are used. Spearman correlation always uses marginal midrank scores.

Details

The procedure computes the concordant and discordant pair counts, n_C and n_D, and sets S = n_C - n_D. Its null variance is

$\mathrm{Var}_0(S) = \sum_{ij} n_{ij}(A_{ij}-D_{ij})^2 - 4S^2/n,$

where A_ij and D_ij are the concordant and discordant contributions for cell (i,j). The test statistic is Z = S/sqrt(Var_0(S)).

The concordance-based measures share a test statistic, but Pearson and Spearman correlations use their own multinomial null variances. Standard errors and confidence intervals are computed separately for each measure. Confidence intervals are two-sided regardless of the test alternative. Rows and columns with zero marginal totals are dropped before computation.

The "somers-c|r" measure is Somers' Delta(C|R), treating the row variable as explanatory; "somers-r|c" is Delta(R|C), treating the column variable as explanatory.

Pearson correlation uses the supplied category scores or the default category scores described above (numeric labels such as 1, 2, and 10 are used as-is). Spearman correlation uses marginal midranks, so observations tied within a category receive the same rank. The correlation standard errors and null variances use multinomial contingency-table formulas. The confidence intervals are Wald intervals.

Value

An object of class "ibist_ordinal_assoc" containing:

estimate

Named vector of selected association estimates.

std.error

Named vector of asymptotic standard errors.

conf.int

Matrix of lower and upper confidence limits.

statistic

Named asymptotic test statistics, one per selected measure.

p.value

Named p-values, one per selected measure.

S, var.S

The observed concordance contrast and its null variance.

concordant, discordant

Concordant and discordant pair counts.

References

Brown, M. B., and Benedetti, J. K. (1977). Sampling behavior of tests for correlation in two-way contingency tables. Journal of the American Statistical Association, 72(358), 309–315.

Goodman, L. A., and Kruskal, W. H. (1954). Measures of association for cross classifications. Journal of the American Statistical Association, 49(268), 732–764.

Kendall, M. G. (1970). Rank Correlation Methods (4th ed.). Griffin.

Somers, R. H. (1962). A new asymmetric measure of association for ordinal variables. American Sociological Review, 27(6), 799–811.

Stuart, A. (1953). The estimation and comparison of strengths of association in contingency tables. Biometrika, 40(1–2), 105–110.

Examples

tab <- matrix(c(12, 3, 1, 4, 9, 2, 1, 5, 13), nrow = 3, byrow = TRUE)
ordinal.assoc(tab)
ordinal.assoc(
  tab,
  measure = c("tau-b", "somers-c|r", "pearson", "spearman")
)

---

One-sample proportion power and sample size calculation

Description

Computes power, sample size, or detectable effect size for a one-sample binomial test. The interface mirrors stats::power.prop.test(), but for the one-sample setting.

Usage

power.p1s.test(
  n = NULL,
  p0 = NULL,
  p1 = NULL,
  power = NULL,
  sig.level = 0.05,
  alternative = c("two.sided", "less", "greater"),
  correct = FALSE,
  exact = FALSE,
  exact.method = c("quantile", "midp", "midp-rand", "cp"),
  strict = FALSE,
  tol = .Machine$double.eps^0.5,
  max_n = 1e+07,
  size.rule = c("close", "minimal"),
  alpha.min.frac = 0.9
)

Arguments

n	Sample size for the single group.
p0	Null hypothesis proportion.
p1	Alternative hypothesis proportion.
power	Desired power.
sig.level	Significance level.
alternative	Character string specifying the alternative hypothesis; one of "two.sided", "less", or "greater".
correct	Logical; if TRUE, apply continuity correction in the normal approximation. The default is FALSE, matching the uncorrected normal approximation in stats::power.prop.test().
exact	Logical; if TRUE, use an exact binomial method.
exact.method	Method used for exact binomial power calculation. "quantile" uses fixed binomial rejection regions for stable power and sample-size inversion; "midp" uses a deterministic rejection region defined by one-sided or central two-sided mid-p values; "midp-rand" uses the quantile rejection boundaries and assigns rejection probability one-half at each boundary; "cp" inverts Clopper–Pearson acceptance regions to guarantee size not exceeding sig.level.
strict	Logical; if TRUE and alternative = "two.sided", computes power accounting for both tails explicitly. Otherwise uses the standard one-sided approximation. The default is FALSE, matching stats::power.prop.test().
tol	Numerical tolerance used in root finding.
max_n	Maximum allowable sample size when solving for n.
size.rule	Character string specifying the rule used to select the sample size when exact = TRUE. Must be one of "close" (default) or "minimal". "close" selects the smallest sample size $n$ such that the achieved significance level satisfies $\alpha_{\mathrm{act}} \le \alpha$, the power meets the target, and the achieved level is not too far below the nominal level: $\alpha_{\mathrm{act}} \ge c \alpha,$ where $c =$ alpha.min.frac. This rule avoids overly conservative designs and makes more efficient use of the type I error budget. "minimal" selects the smallest sample size $n$ such that the achieved significance level satisfies $\alpha_{\mathrm{act}} \le \alpha$ and the power is at least the target value. This is the classical exact test design but may be conservative, with $\alpha_{\mathrm{act}} \ll \alpha$.
alpha.min.frac	Numeric value in (0, 1) specifying the minimum fraction of the nominal significance level required for the achieved level when size.rule = "close". The achieved significance level must satisfy $\alpha_{\mathrm{act}} \ge \text{alpha.min.frac} \times \alpha.$ The default value 0.9 requires the achieved level to be at least 90% of the nominal level. Smaller values allow more conservative designs, while larger values enforce closer agreement with the nominal level. See also achieved.sig.level in the returned value. For exact.method = "midp-rand", this component is the randomized null rejection probability, including one-half of the mass at each rejection boundary.

Details

Power can be computed using a normal approximation or exact binomial methods. When exact = TRUE, the exact test is determined by exact.method.

For exact.method = "quantile", the one-sided critical boundary is the most extreme binomial quantile whose null tail probability does not exceed sig.level. For a two-sided test, up to half of sig.level is allocated to each tail.

For exact.method = "midp", the lower and upper mid-p tails at an observed count $x$ are, respectively,

$P_0(X < x) + \tfrac{1}{2}P_0(X = x)$

and

$P_0(X > x) + \tfrac{1}{2}P_0(X = x).$

A one-sided test rejects when the applicable mid-p tail is no greater than sig.level. The central two-sided test implemented here rejects when twice the smaller mid-p tail is no greater than sig.level. Other two-sided mid-p tests may use different orderings of the sample space. Power and achieved.sig.level are the full probabilities of this deterministic rejection region under $p_1$ and $p_0$, respectively. Unlike the other deterministic exact methods, a mid-p test need not control its actual size at the nominal level for every $p_0$.

For exact.method = "midp-rand", the quantile method's critical boundaries are retained, but an observation on either boundary is assigned rejection probability one-half. Thus, each reported probability consists of the probability strictly beyond the boundaries plus one-half of the boundary mass. This randomized construction explains its connection to mid-p tail probabilities.

For exact.method = "cp", rejection is determined by the one-sided exact binomial p-value or, for a two-sided test, twice the smaller exact tail probability. This is equivalent to inversion of the corresponding Clopper–Pearson acceptance region and controls size at sig.level.

Exactly one of n, p0, p1, power, or sig.level must be NULL; the missing quantity is solved numerically.

Value

An object of class "power.htest" containing the computed quantity and test specifications. For exact calculations, achieved.sig.level is the null rejection probability. With exact.method = "midp-rand", it includes the one-half randomized contribution at each boundary.

Examples

## Normal approximation (default)
power.p1s.test(n = 50, p0 = 0.1, p1 = 0.25)

## Exact binomial power (quantile-based, default exact method)
power.p1s.test(n = 50, p0 = 0.1, p1 = 0.25, exact = TRUE)

## Exact power for the deterministic mid-p test
power.p1s.test(
  n = 50, p0 = 0.1, p1 = 0.25,
  exact = TRUE, exact.method = "midp"
)

## Randomized mid-p power at the quantile rejection boundaries
power.p1s.test(
  n = 50, p0 = 0.1, p1 = 0.25,
  exact = TRUE, exact.method = "midp-rand"
)

## Exact Clopper--Pearson power (guaranteed size control)
power.p1s.test(
  n = 50, p0 = 0.1, p1 = 0.25,
  exact = TRUE, exact.method = "cp"
)

---

Power Calculation for Two-Sample Proportion Test with Unequal Group Sizes

Description

Computes power, sample size, detectable proportions, or significance level for a two-sample test comparing proportions, allowing unequal group sizes through a user-specified allocation ratio.

Usage

power.p2s.test(
  n = NULL,
  p1 = NULL,
  p2 = NULL,
  sig.level = 0.05,
  power = NULL,
  group.rate = 1,
  alternative = c("two.sided", "one.sided"),
  correct = FALSE,
  strict = FALSE,
  tol = .Machine$double.eps^0.25
)

Arguments

n	Sample size in the first group. The second group size is n * group.rate.
p1	Proportion in the first group.
p2	Proportion in the second group.
sig.level	Significance level (type I error rate), must be in [0, 1].
power	Target power of the test.
group.rate	Ratio of sample sizes between the second and first groups, defined as $n_2 / n_1$. Defaults to 1 (equal group sizes).
alternative	Character string specifying the alternative hypothesis, either "two.sided" or "one.sided".
correct	Logical; if TRUE, apply continuity correction in the normal approximation. The default is FALSE, matching the uncorrected normal approximation in stats::power.prop.test().
strict	Logical; if TRUE and alternative = "two.sided", computes power accounting for both tails explicitly. Otherwise uses the standard one-sided approximation. The default is FALSE, matching stats::power.prop.test().
tol	Tolerance for numerical root-finding when solving for an unknown parameter.

Details

Exactly one of n, p1, p2, power, or sig.level must be NULL, and that parameter is solved from the others.

The calculation is based on the normal approximation to the binomial distribution. Let $n_1 = n$ and $n_2 = n \times group.rate$. The test statistic uses a pooled variance under the null and an unpooled variance under the alternative.

When strict = TRUE for a two-sided test, power is computed as the sum of tail probabilities. Otherwise, a one-sided approximation is used.

The code is adapted from stats::power.power.test() by allowing different group sizes and continuity correction.

Value

An object of class "power.htest" with components:

n	Sample size in the first group.
p1	Proportion in the first group.
p2	Proportion in the second group.
sig.level	Significance level.
power	Power of the test.
alternative	Type of alternative hypothesis.
note	Clarifies that n refers to the first group size and states whether continuity correction was applied.
method	Description of the method.

Note

This function generalizes stats::power.prop.test() by allowing unequal group sizes through group.rate. When group.rate = 1, the two functions are equivalent.

Examples

# Power with unequal group sizes (n2 = 2 * n1)
power.p2s.test(n = 100, p1 = 0.3, p2 = 0.5, group.rate = 2)

# Required sample size in group 1
power.p2s.test(p1 = 0.3, p2 = 0.5, power = 0.8, group.rate = 1.5)

---

Print an ibist Confidence Interval Object

Description

Print an ibist Confidence Interval Object

Usage

## S3 method for class 'ibist_ci'
print(x, digits = getOption("digits"), ...)

Arguments

x	An object of class "ibist_ci".
digits	Number of significant digits to print.
...	Reserved for future extensions.

---

Confidence Intervals for Paired Proportion Difference

Description

Computes confidence intervals for the difference of two paired proportions.

Usage

prop.paired.ci(
  b,
  c,
  n,
  conf.level = 0.95,
  method = c("score", "wald", "waldcc", "agresti-min", "wang"),
  ...
)

Arguments

b	Number of discordant pairs favorable to group 1.
c	Number of discordant pairs favorable to group 2.
n	Total number of matched pairs.
conf.level	Confidence level for the interval.
method	Method or methods for confidence interval. One or more of "score", "wald", "waldcc", "agresti-min", or "wang".
...	Additional arguments for the selected method. For "score", these include tol; for "wang", they include precision, grid.one, and grid.two.

Details

For a matched-pair table with discordant counts $b$ and $c$, the estimated paired difference is $(b - c) / n$. The "wald" and "waldcc" methods use the paired Wald standard error, with "waldcc" adding the continuity-correction term $1 / n$. The "agresti-min" method applies the Wald calculation to $b + 1/2$, $c + 1/2$, and $n + 2$. The "score" method inverts Tango's score statistic by numerically maximizing the multinomial likelihood subject to the candidate difference. The "wang" method computes Wang's exact interval based on an inductive ordering of the paired sample space. This method can be time-consuming when the sample size is moderate because it searches over nuisance parameters and refines the confidence limits numerically.

Value

An object of class "ibist_ci" containing the common estimate and a confidence-limit matrix with one row per method.

Method-specific arguments

For method = "score", tol is the numerical tolerance used when solving for the confidence limits; its default is 1e-8. For method = "wang", precision controls the numerical precision of the confidence limits (default 1e-5), while grid.one and grid.two control the sizes of the two nuisance-parameter search grids (defaults 30 and 20, respectively). Increasing these values can improve numerical accuracy at the cost of longer computation.

References

Tango, T. (1998). Equivalence test and confidence interval for the difference in proportions for the paired-sample design. Statistics in Medicine, 17, 891–908.

Agresti, A., and Min, Y. (2005). Simple improved confidence intervals for comparing matched proportions. Statistics in Medicine, 24, 729–740.

Wang, W. (2012). An inductive order construction for the difference of two dependent proportions. Statistics & Probability Letters, 82, 1623–1628.

Examples

prop.paired.ci(b = 8, c = 25, n = 180)
prop.paired.ci(b = 8, c = 25, n = 180, method = "waldcc")
prop.paired.ci(b = 8, c = 25, n = 180, method = "agresti-min")
prop.paired.ci(b = 3, c = 0, n = 4, method = "wang",
                precision = 0.0001)

---

Confidence Intervals for a One-Sample Poisson Rate

Description

Computes confidence intervals for a Poisson rate parameter $\lambda = x / T$ using several methods.

Usage

rate.1s.ci(
  x,
  exposure = 1,
  conf.level = 0.95,
  method = c("exact", "score", "wh", "wald", "log"),
  correct = TRUE,
  ...
)

Arguments

x	Non-negative integer. Observed number of events.
exposure	Positive numeric. Exposure time (or total time at risk).
conf.level	Confidence level. Default is 0.95.
method	Method or methods for confidence interval. One or more of: "exact": Exact (Garwood) interval "score": Score interval (inversion of score test) "wh": Wilson-Hilferty interval "wald": Wald interval (normal approximation) "log": Log-transformed interval
correct	Logical. Apply continuity correction for "score", "wh", and "wald" methods. Default is TRUE.
...	Reserved for future extensions.

Details

This function provides several confidence intervals for the Poisson rate $\lambda = x / T$.

Exact (Garwood): Based on inversion of the Poisson test using chi-square quantiles:

$\left[ \frac{1}{2T} \chi^2_{2x, \alpha/2}, \frac{1}{2T} \chi^2_{2(x+1), 1-\alpha/2} \right].$

Score: Obtained by inverting the score test. Computed from solving quadratic equations for boundaries.

Wilson-Hilferty: The Wilson–Hilferty (WH) interval is based on a cube-root transformation of a chi-square approximation, which reduces skewness and yields an approximately normal pivot, leading to the cubic form of the limits. The lower limit is

$\frac{x}{T} \left(1 - \frac{1}{9x} + \frac{z_{\alpha/2}}{3\sqrt{x}}\right)^3$

and the upper limit is

$\frac{x + 1}{T} \left(1 - \frac{1}{9(x + 1)} + \frac{z_{1-\alpha/2}}{3\sqrt{x + 1}}\right)^3.$

Wald: For $x > 0$,

$\hat{\lambda} \pm z_{\alpha/2} \sqrt{\hat{\lambda}/T}.$

For $x = 0$, lower limit 0, upper limit $-\log(\alpha/2)/T$.

Log: For $x > 0$, applies normal approximation to $\log(\lambda)$ and transforms back. For $x = 0$, lower limit 0, upper limit $-\log(\alpha/2)/T$.

When correct = TRUE, continuity correction is applied on the count scale for methods that support it.

Value

An object of class "ibist_ci" containing the common estimate and a confidence-limit matrix with one row per method.

Examples

rate.1s.ci(5, 10)
rate.1s.ci(5, 10, method = "score")
rate.1s.ci(5, 10, method = c("exact", "score", "log"))
rate.1s.ci(0, 10, method = "exact")

---

Confidence Intervals for a Two-Sample Poisson Rate Ratio

Description

Computes confidence intervals for the rate ratio $\rho = \lambda_1 / \lambda_2$ from two independent Poisson counts with known exposures.

Usage

rate.2s.ci(
  x,
  exposure = c(1, 1),
  conf.level = 0.95,
  method = c("log", "score", "exact"),
  ...
)

Arguments

x	Length-2 vector of non-negative integer event counts.
exposure	Length-2 positive numeric vector of exposures.
conf.level	Confidence level. Default is 0.95.
method	Method or methods for confidence interval. One or more of: "log": Log-Wald interval "score": Transformed Wilson score interval "exact": Transformed exact Clopper-Pearson interval
...	Reserved for future extensions.

Details

For two independent Poisson variables $X_i \sim \mathrm{Pois}(\lambda_i T_i)$, the point estimate of the rate ratio is

$\hat\rho = \frac{X_1/T_1}{X_2/T_2}.$

The "log" method uses the large-sample interval

$\log(\hat\rho) \pm z_{1-\alpha/2}\sqrt{1/X_1 + 1/X_2}$

and exponentiates the endpoints. This method requires both counts to be positive.

The "score" and "exact" methods use the conditional binomial representation. Given $n = X_1 + X_2$, $X_1 \sim \mathrm{Binom}(n, \pi)$, where

$\pi = \frac{T_1\rho}{T_1\rho + T_2}.$

Confidence limits $L$ and $U$ for $\pi$ are transformed to confidence limits for $\rho$ as

$\left( \frac{L T_2}{(1 - L)T_1}, \frac{U T_2}{(1 - U)T_1} \right).$

The "score" method uses Wilson score limits for $\pi$; the "exact" method uses the exact Clopper-Pearson limits returned by binom.test.

Value

An object of class "ibist_ci" containing the rate estimates and a confidence-limit matrix with one row per method.

Examples

rate.2s.ci(c(151, 55), exposure = c(57518.1, 74573.5))
rate.2s.ci(c(151, 55), exposure = c(57518.1, 74573.5), method = "score")
rate.2s.ci(c(151, 55), exposure = c(57518.1, 74573.5),
           method = c("log", "score"))
rate.2s.ci(c(9, 12), exposure = c(1817.6, 7496.3), method = "exact")

---

Test of Poisson Rates Using Normal Approximation

Description

Performs a large-sample (normal) test for one or two Poisson rates with known exposures. The test is carried out on the observed count scale and then translated to rates for reporting.

Usage

rate.test(
  x,
  exposure = 1,
  r = 1,
  alternative = c("two.sided", "less", "greater"),
  conf.level = 0.95,
  correct = TRUE
)

Arguments

x	a vector of event counts. A single value specifies a one-sample test; a vector of length two specifies a two-sample comparison.
exposure	a vector of exposures corresponding to x (e.g., person-time). Must have the same length as x.
r	a positive number specifying the null rate per unit exposure, $\lambda_0$, for a one-sample test. Default to 1.0. Ignored for two-sample tests.
alternative	a character string specifying the alternative hypothesis, one of "two.sided", "less", or "greater".
conf.level	confidence level for the confidence interval. For one-sample tests, the interval is the score confidence interval for the rate. For two-sample tests, the interval is the log-Wald confidence interval for the rate ratio.
correct	logical; if TRUE, applies a continuity correction on the count scale. For one-sample tests, a $\pm 0.5$ correction is applied to $x$ in the direction determined by alternative. For "two.sided", the correction uses the sign of $x-\mu_0$.

Details

One-sample test: Assume $X \sim \mathrm{Pois}(\mu)$ with $\mu = \lambda T$. The null hypothesis is $H_0:\lambda=\lambda_0$, where $\lambda_0=r$. Let $\mu_0=\lambda_0 T$. The normal approximation gives

$Z = \frac{(x - c) - \mu_0}{\sqrt{\mu_0}} \approx N(0,1),$

where $c$ is 0 if correct = FALSE and is a continuity correction on the count scale if correct = TRUE. Two-sided p-values use $2\{1-\Phi(|Z|)\}$.

Two-sample test: Assume independent $X_1 \sim \mathrm{Pois}(\lambda_1 T_1)$ and $X_2 \sim \mathrm{Pois}(\lambda_2 T_2)$ with known exposures $T_1$ and $T_2$. The null hypothesis is $H_0:\lambda_1=\lambda_2$.

For hypothesis testing, the function uses the exact conditional representation under $H_0$:

$X_1 \mid (X_1+X_2=n) \sim \mathrm{Binom}\!\left(n,\; \frac{T_1}{T_1+T_2}\right),$

and applies the same large-sample normal approximation as prop.test (with optional Yates continuity correction) to obtain the test statistic and p-value.

For estimation and confidence intervals, the inference target is the rate ratio, $\rho = \lambda_1 / \lambda_2$. The point estimate is $\hat\rho = (X_1/T_1) / (X_2/T_2)$, and the confidence interval is constructed on the log scale using estimated standard error

$\sqrt{1/X_1 + 1/X_2}.$

The continuity correction affects the hypothesis test but not the confidence interval for $\rho$.

Confidence intervals: For one-sample tests, the confidence interval is the score interval for $\lambda$, as computed by rate.1s.ci with method = "score". For two-sample tests, the confidence interval is for the rate ratio $\lambda_1/\lambda_2$. The log-Wald interval requires both event counts to be positive.

Value

An object of class "htest" containing:

statistic	the standardized normal statistic z for one-sample tests or the chi-square statistic from prop.test for two-sample tests.
parameter	degrees of freedom (df = 1).
p.value	the p-value.
conf.int	a confidence interval for the rate in one-sample tests or for the rate ratio, $\lambda_1/\lambda_2$, in two-sample tests.
estimate	estimated rate (one-sample) or estimated rates and rate ratio (two-sample).
null.value	the null rate (one-sample) or the null rate ratio $\lambda_1/\lambda_2=1$ (two-sample).
alternative	the alternative hypothesis.
method	a character string describing the test.
data.name	a character string describing the data.

See Also

rate.2s.ci, prop.test, poisson.test

Examples

## One-sample test: compare observed rate to a reference unit rate
rate.test(x = 411, exposure = 25800, r = 0.0119)
rate.test(x = 411, exposure = 25800, r = 0.0119, correct = FALSE)

## Two-sample test: compare two Poisson rates
rate.test(x = c(12, 5), exposure = c(100, 80))
rate.test(x = c(12, 5), exposure = c(100, 80), correct = FALSE)

## One-sided alternative
rate.test(x = 411, exposure = 25800, r = 0.0119,
          alternative = "greater")

---

Neonatal Respiratory Distress Syndrome Data

Description

Aggregated counts of neonatal outcomes by birth weight, surfactant use, and survival status.

Usage

rds

Format

A data frame with 16 rows and 4 variables:

bwt

Birth weight category in grams (factor with 4 levels).

surf

Surfactant use (Yes/No).

death

Outcome (1 = death; 0 = alive).

count

Number of infants in each group.

Details

The data are presented in aggregated form. Each row corresponds to a combination of birth weight category, surfactant use, and outcome, with count indicating the number of observations.

Analyses should account for the grouped structure, e.g., via weighted models or binomial responses.

Source

Classical neonatal RDS dataset (exact source to be specified).

---