`fit1()` currently returns just estimates and SEs. Gary requested having it return the full variance-covariance matrix for the coefficient estimates.