I am trying to use the delta method to compute confidence intervals for non-linear regression models. My standard errors is an array of 5,3,3 however I think standards errors should be 5 values, not matrices. I am getting a value error when trying to add or subtract a standard_error array from my predicted values. Any help would be appreciated.
df = pd.DataFrame({
'cumsum_days': [1,2,3,4,5],
'pred': [388.259631, 368.389649, 349.754534, 332.264306, 315.836485]})
cov = np.array([[2.67918945e+04, 2.62421460e+02, 9.08452505e+00],
[2.62421460e+02, 4.31869566e+00, 1.24995272e-01],
[9.08452505e+00, 1.24995272e-01, 3.90413410e-03]])
# estimate confidence interval for predicted probabilities
gradient = np.gradient(df['pred'], df['cumsum_days'])
std_errors = np.array([np.sqrt(np.dot(np.dot(g, cov), g)) for g in gradient])
c = 1.96 # multiplier for confidence interval
upper = np.maximum(0, np.minimum(1, (df['pred'] + std_errors * c)))
lower = np.maximum(0, np.minimum(1, (df['pred'] - std_errors * c)))
ValueError here:
print(df['pred'] + std_errors)
ValueError: operands could not be broadcast together with shapes (5,) (5,3,3)
That's a shape problem, your gradient is of shape (5,) and your cov is of shape (3,3), by iterating over each element g in gradient, np.dot(g, cov) will just multiply each element of the cov matrix by a number g (1), giving a shape(3,3) and create a list of 5 elements of shape (3,3), hence the shape (5,3,3).
I don't know the details about the delta method but it seems to me that you should either have a cov matrix of shape (5,5) or have a gradient of shape (3,) by removing the first and last value and then you will be able to do std_error = np.sqrt(np.dot(np.dot(gradient, cov), gradient))