Add Anderson rubin inference

[kv9898]

This pull request adds support for the Anderson-Rubin (AR) weak-instrument robust test for instrumental variable (IV) models in the fixest package. It introduces two new functions, ar_test() and ar_confint(), along with their print methods, and includes comprehensive tests to ensure correct behavior across a variety of scenarios.

This is requested as part of #442.

New IV Robustness Features:

• Added the Anderson-Rubin (AR) test for IV models via the new ar_test() function, which provides weak-instrument robust inference for coefficients of endogenous regressors. The related ar_confint() function computes AR confidence intervals by inverting the test.
• Exported the new ar_test and ar_confint functions in the NAMESPACE to make them available for users.
• Registered S3 print methods for the new result classes fixest_ar and fixest_ar_confint in the NAMESPACE.

Testing and Validation:

• Added extensive tests for the new AR test and confidence interval functions in tests/fixest_tests.R, covering correct output structure, distribution selection, error handling for invalid input, and integration with model features like clustering, weights, and multiple instruments.

[kylebutts]

Happy to spend a few hours reviewing this code, but could you first create a small working example showing this matches the results from ivDiag? https://yiqingxu.org/packages/ivDiag/articles/iv_tutorial.html

Also related to #593

[kv9898]

Happy to spend a few hours reviewing this code, but could you first create a small working example showing this matches the results from ivDiag? https://yiqingxu.org/packages/ivDiag/articles/iv_tutorial.html

Also related to #593

Hi @kylebutts . During my internal testing, I benchmarked it against both ivmodel and ivDiag. However, I believe that ivDiag lacks precision: for when a closed form solution is possible (iid se with no fixed effects), ivmodel (which has an exact solution algorithm) and my PR's exact and numeric solution yield same results, but ivDiag yields a different one. I will present a reproducible example shortly

My PR is superior to

• ivDiag, for mine is more precise (to the 8 decimal places by default), and allows for closed-form solutions when possible. Mine is also not conceivably slower than ivDiag (I don't feel the difference)
• ivmodel, as it allows for different vcov choices, fixed effects. It is also more user-friendly by also taking an existing feols object.

[kv9898]

Here's a simple repex:

library(fixest)
library(ivmodel)
library(ivDiag)

# Data generation ----
set.seed(123)

gamma <- 0.5
beta <- 2
obs <- 100

e_1 <- rnorm(obs, mean = 0, sd = 1)
z <- rnorm(obs, mean = 0, sd = 1)
x <- gamma * z + e_1
w <- rnorm(obs, mean = 0, sd = 1)
e_2 <- rnorm(obs, mean = 0, sd = 1) + 0.8 * e_1
y <- beta * x + w + e_2

df <- data.frame(y, x, z, w)

# ivmodel (B0=0) ----
ivm <- ivmodel(y, x, z, w) # This uses an exact closed-form solution
summary(ivm) # AR: F=10.79273, df1=1, df2=97, p=0.0014187
# CI: [1.2562, 2.5130]

ivdiag <- AR_test(df, "y", "x", "z", "w")
summary(ivdiag) # AR: F=16.2581, CI: [1.4889, 2.4556] NOT PRECISE!

reg_iv <- feols(y ~ w | x ~ z, data = df) # our PR
(test0 <- ar_test(reg_iv)) # AR: F=10.79, df1=1, df2=97, p=0.001419
# CI: [1.256, 2.513]
(test0n <- ar_test(reg_iv, ci = "numeric"))  # try our grid search algorithm
# CI: [1.256, 2.513], the numeric results are designed to be accurate to the 8th decimal place

# ivmodel (B0=1) ----
ivm_1 <- ivmodel(y, x, z, w, beta0 = 1)
summary(ivm_1) # AR: F=5.647947, df1=1, df2=97, p=0.019439
# CI: [1.2562, 2.5130]

# fixest (IV) ----
(test1 <- ar_test(reg_iv, beta0 = 1)) # must use my custom function
(test2 <- ar_test(reg_iv, beta0 = 1, ci = "numeric"))
# AR: F=5.647947, df1=1, df2=97, p=0.019439
# CI: [1.2562, 2.5130], the numeric results are designed to be accurate to the 8th decimal place

[kv9898]

@kylebutts Just realized I have probably seen you in Positron discussions? It's nice to see you in another repo again. I'm also interested in adding valid t-ratio inference from by Lee et al. (2022) and Montiel Olea &
Pueger's (2013) effective first stage F-test. But they will be separate PRs and I want to get this done first.

[lrberge]

Hi @kv9898 and thanks for this PR!
I'm not sure creating new functions is that necessary (although I understand the appeal!).
I'd like also to think to the integration with, eg, etable. I'll have a look in January.
Thanks again and in the meantime merry Christmas!

[kylebutts]

I'm working through the code and have some changes I'll commit to this branch when I finish.

I agree with @lrberge that we should think more about how users would use this. On the one-hand, we form a p-value or a confidence interval for the endogenous variable. We could imagine plugging in AR-based p-value or CI in that row of the regression table. Would it make sense to have a different kind of CI for the endogenous variable only?

Alternatively, it could be treated similar to the other fitstats and go below the coefficients

What do we think?

[grantmcdermott]

Alternatively, it could be treated similar to the other fitstats and go below the coefficients

Without having looked at the code, this would have been my first inclination.

[lrberge]

Without having looked at the code, this would have been my first inclination.

+1!

Full integration with streamlined user experience needs to be thought carefully.
This would change the coefficient table, it can be a chance to make a template for future similar additions.

[kv9898]

@kylebutts Happy New Year! Any updates?

[lrberge]

Hi @kv9898, to do it correctly it actually implies big and well thought internal changes. I'm sorry but this won't be direct. I beg for your patience.

[kv9898]

Hi @kv9898, to do it correctly it actually implies big and well thought internal changes. I'm sorry but this won't be direct. I beg for your patience.

Ohh I thought we were just planning for some nice wrappers around the function. No rush, take your time! Thanks for the quick reply, though.