Error Message for `contrast = "RR", stratified = TRUE, weighting = "MN"`

[lovestat]

The below code will produce an error message "Error in while (wtdiff > MNtol) { : missing value where TRUE/FALSE needed".

I am not sure if it could be a potential bug. Could you please check it a bit?

library(ratesci)
  ## Example data
  n.1.1 <-  199
  n.1.2 <-  200
  n.2.1 <-  197
  n.2.2 <-  202

  r.1.1 <-  1
  r.1.2 <-  2
  r.2.1 <-  0
  r.2.2 <-  0

  fit <- scoreci(x1 = c(r.1.1, r.1.2),
                 n1 = c(n.1.1, n.1.2),
                 x2 = c(r.2.1, r.2.2),
                 n2 = c(n.2.1, n.2.2),
                 contrast = "RR",
                 stratified = TRUE,
                 weighting = "MN",
                 skew = FALSE)

[lovestat]

Thank you, Pete. It works great now.

Based on this modification, I found two different behaviors of scoreci () with contrast = "RR" and contrast = "OR":

contrast = "OR" notifies that

[1] "Note: at least one stratum contributed no information and was removed" [1] "Note: only one stratum!"

## Example data
n.1.1 <-  199
n.1.2 <-  200
n.2.1 <-  197
n.2.2 <-  202

r.1.1 <-  3
r.1.2 <-  0
r.2.1 <-  0
r.2.2 <-  0

fit <- scoreci(x1 = c(r.1.1, r.1.2),
               n1 = c(n.1.1, n.1.2),
               x2 = c(r.2.1, r.2.2),
               n2 = c(n.2.1, n.2.2),
               contrast = "OR",
               stratified = TRUE,
               weighting = "MN",
               skew = FALSE)

contrast = "RR" is silent.

## Example data
n.1.1 <-  199
n.1.2 <-  200
n.2.1 <-  197
n.2.2 <-  202

r.1.1 <-  3
r.1.2 <-  0
r.2.1 <-  0
r.2.2 <-  0

fit <- scoreci(x1 = c(r.1.1, r.1.2),
               n1 = c(n.1.1, n.1.2),
               x2 = c(r.2.1, r.2.2),
               n2 = c(n.2.1, n.2.2),
               contrast = "RR",
               stratified = TRUE,
               weighting = "MN",
               skew = FALSE)

Is this expected? I suppose that the empty strata will not contribute information in both cases. Thank you!

[PeteLaud]

Yes, this is expected.

For a stratified analysis, there are certain weighting schemes (set as default for each contrast) that produce analysis which is consistent with the Cochran-Mantel-Haenszel test for superiority. (That is, MH weighting for RD and RR, and INV weighting for OR. See Laud 2021 for more details, DOI: 10.1002/pst.2144.) For these specific analyses it is possible to retain all strata as input data to the function, which is desirable in order to satisfy the intention-to-treat principle.

“Double-zero” strata contribute no information to the test or confidence interval for these CMH-consistent analyses of OR and RR. Specifically, for the OR score, the variance of S(𝜃) is infinite, resulting in zero stratum weight with INV weighting, and for the MH-weighted RR score, the stratum’s contributions to the numerator and denominator of the score are both zero. (Interestingly, even for the RD contrast, the superiority test p-value is unchanged when the double-zero strata are removed, even though such strata do contain relevant information which affects the estimated difference and its confidence interval.)

Hence the same results are obtained whether double-zero strata are retained or omitted, but it is preferable to retain them, to avoid explicitly removing patients from an ITT analysis. Furthermore, such strata would not need to be removed for the CMH test. (See example code below.)

For the alternative weighting schemes, such as MH or MN weighting for OR, the infinite stratum variance produces a zero combined score for all values of 𝜃, so it becomes necessary to remove the double-zero strata.

Having said that, double-zero strata do pose a difficulty for the evaluation of heterogeneity (they are removed for the random-effects TDAS method), and I recognise that excluding these strata for some weighting schemes but not others leads to inconsistency in the degrees of freedom used in the heterogeneity tests. I plan to update the heterogeneity test code so that double-zero strata are not considered (for RR/OR), and to include details of the degrees of freedom in the Qtest output object for transparency.

The dropzeros argument is available to force the exclusion of double-zero strata, but this may be redundant once I’ve made the above changes to the heterogeneity test output.

library(ratesci)
# Example dataset with two "double-zero" strata
x1 <- c(0, 0, 2, 0, 3, 5)
n1 <- c(10, 10, 10, 10, 10, 12)
x2 <- c(0, 0, 0, 1, 1, 1)
n2 <- c(10, 10, 10, 10, 10, 10)
CMHarray <- aperm(array(c(x1, n1 - x1, x2, n2 - x2), dim = c(6, 2, 2)), c(2, 3, 1))

# Same p-value for RD, RR and OR: p = 0.04948997
scoreci(x1, n1, x2, n2, contrast = "RD", skew = FALSE, stratified = TRUE) # Default weighting = 'MH'
scoreci(x1, n1, x2, n2, contrast = "RR", skew = FALSE, stratified = TRUE) # Default weighting = 'MH'
scoreci(x1, n1, x2, n2, contrast = "OR", skew = FALSE, stratified = TRUE) # Default weighting = 'INV'
# $estimates (for OR)
#        Lower      MLE   Upper 
# [1,] 1.002706 3.738034 14.4833
# with skewness correction:
# [1,] 1.005936 3.633184 18.03531 

# Matches CMH test, with or without the double-zero strata: p = 0.04949
mantelhaen.test(CMHarray, correct = FALSE)
mantelhaen.test(CMHarray[, , 3:6], correct = FALSE)
# (Note the contradictory CI for OR produced by this function - 95% CI includes OR = 1 when p-value < 0.05)
# common odds ratio: 
# 3.762646 
# 95 percent confidence interval:
# 0.9533434 14.8503718

[PeteLaud]

I have now updated the code so that double-zero strata need not be excluded for weighting = "MN" - I found this also produces a test consistent with the CMH test (for all binomial contrasts). I have also made the changes mentioned above for the heterogeneity test output.