Bug when IV has three categories?

[rmtrane]

I've been using your package for a project I've been working on, and I think there might be a bug. If I'm way off, I apologize for the inconvenience.

Anyway, I was running a bunch of simulations, and encountered a scenario where using bpbounds on a table of conditional probabilities claims that the IV inequalities are violated. However, I don't think that's the case.

tabp <- array(data = c(0.5279183, 0.02208171, 0.1220817, 0.3279183, 
                       0.0849975, 0.3150025, 0.4650025, 0.1349975, 
                       0.1132796, 0.3867204, 0.1867204, 0.3132796),
              dim = c(2,2,3),
              dimnames = list(x = c(0,1),
                              y = c(0,1),
                              z = c(0,1,2)))

rowSums(apply(tabp, c(1,2), max))

The last line of code above returns \max_z P(X = x, Y = 0 | Z = z) + \max_z P(X = x, Y = 1 | Z = z) for x=0,1. Both are less than 1, indicating the IV inequalities are not violated. But when I run bpbounds(as.table(tabp)), it returns Instrumental inequality: FALSE.

I'm also slightly confused by one of the examples from your vignette (which is otherwise great, by the way!):

mt3 <- c(.83, .05, .11, .01, .88, .06, .05, .01, .72, .05, .20, .03)
p3 = array(mt3,
           dim = c(2, 2, 3),
           dimnames = list(
             x = c(0, 1),
             y = c(0, 1),
             z = c(0, 1, 2)
           ))

rowSums(apply(p3, c(1,2), max))

Here, P(X = 0, Y = 0 | Z = 1) + P(X = 0, Y = 1 | Z = 2) = 1.08, i.e. the IV inequalities are violated, but this is not caught by bpbounds(as.table(p3))

I'm not sure what's going on, and again, if I'm missing something obvious, please let me know, and I apologize. But thought I would bring this to your attention in case there is indeed a bug in your package.

[remlapmot]

Thanks for this.

I will have a look and get back to you.

[sachsmc]

Hello, this may be old news but I have also recently had issues with the 3 level IV example. I think I have identified the problem:

In the code here: https://github.com/remlapmot/bpbounds/blob/master/R/bpbounds_calc_tri_z3.R#L4 the expected order of the probabilities is to have p(Y = 1, X = 0 | Z = z) coming before p(Y = 0, X = 1 | Z = z). However when using a single index on a 3-dimensional array it goes from top down then left to right. So, in the example array above, p[2] is in fact p(Y = 0, X = 1 | Z = 0).

Compare that to the code here: https://github.com/remlapmot/bpbounds/blob/master/R/bpbounds_calc_tri_z2.R#L4 which has the correct order, i.e., p[3] comes before p[2], and p[7] before p[6].

Now this wasn't apparently a problem for your Mendelian randomization example (you get the same bounds as the stata package) because the probabilities were entered in the expected order, but the dimnames of the array are wrong, i.e., mt3[2] = p(Y = 1, X = 0 | Z = 0) and mt3[3] = p(Y = 0, X = 1 | Z = 0) and likewise for Z = 1, 2 .

So, an easy fix for the example is to create the array with dimnames for x and y swapped position:

mt3 <- c(.83, .05, .11, .01, .88, .06, .05, .01, .72, .05, .20, .03)
p3 = array(mt3,
       dim = c(2, 2, 3),
       dimnames = list(
         y = c(0, 1),
         x = c(0, 1),
         z = c(0, 1, 2)
       ))

Then probabilities are labeled correctly.

Hopefully that is clear enough, and let me know if I missed something!

Best regards,
Michael