When calculating the Jacobian for the family, the lpi does not include any information about the thetas, so they are missed. This causes a mismatch of dimensions between Dfx and the output of varFun.
Also need to look at the model matrix, the theta columns are set to 0 when fitting the model, but we might need them to be 1 to properly calculate the jacobian.
When looking at the linear predictor, anything related to the theta should be dropped. In particular, the columns of 0 can be dropped from the model, as the theta has no relationship with the linear predictor.
When calculating the Jacobian for the family, the lpi does not include any information about the thetas, so they are missed. This causes a mismatch of dimensions between Dfx and the output of varFun.
Also need to look at the model matrix, the theta columns are set to 0 when fitting the model, but we might need them to be 1 to properly calculate the jacobian.