1
$\begingroup$

General question

When I perform a logistic regression using lrm and specify weights for the observations, I get the following warning message:

Warning message: In lrm(Tag ~ DLL, weights = W, data = tagdata, x = TRUE, y = TRUE) : currently weights are ignored in model validation and bootstrapping lrm fits

My interpretation is that everything that the rms package will tell me regarding goodness-of-fit, notably using the residuals.lrm tool, is wrong. Is this correct?

Specific example

To be more specific, I have working example. All the code and output can be found in this GitHub repository. I have two CSV tables of data, toystudy.csv and realstudy.csv. There are three columns in each:

  1. The binomial response $y$ (0 or 1) [called Tag in code]
  2. The predictor $x$ [called DLL in code]
  3. The weight for the observation [called W in code]

The former is simulated data, where all the weights are unity and where a logistic regression $log(\pi) = \theta_0 + \theta_1 x$ should fit the data perfectly. The latter is real data from my analysis, where the validity of this simple model is in question. The real data has weighted observations. (Some of the weights are negative, but there is a well-defined reason for this). The analysis code in contained completely in regressionTest.R; the meat of the code is

library(rms)
fit <- lrm(Tag ~ DLL, weights = W, data = tagdata, x=TRUE, y=TRUE)
residuals(fit,"gof")

Here are the results for the two tables of data.

Case 1: Toy data

The goodness-of-fit claimed by lrm (which is something called the le Cessie-van Houwelingen-Copas-Hosmer test, I understand?) is very good:

enter image description here

This is confirmed by grouping the data into 20 quantiles of the predictor and overlaying the predicted success rate over the average actual success rate:

enter image description here

Case 2: Real data

In this case, the goodness-of-fit reported by lrm is horrendous:

enter image description here

However, I don't think it should be that bad. Again grouping the data into quantiles, and taking into account the weights when computing the average values in each bin:

enter image description here

Comparing the prediction to the observed values and their standard errors, I don't think this is that bad (the error bars here depend on how the standard error on a weighted mean is computed, so they might not be 100% right, but should at least be close). On the other hand, if I produce the same plot while ignoring the weights:

enter image description here

I can definitely imagine this fit being as poor as the goodness-of-fit test says.

Conclusion

So, is residuals.rm simply ignoring the weights when it calculates its goodness-of-fit statistic? And if so, is there any R package that will do this correctly?

$\endgroup$
  • $\begingroup$ I should add: The reason I am using lrm and not glm is that the latter rejects negative weights, which are crucial for my analysis. Negative weights may be suspicious, but they are well defined here. The observations belong to two categories, S(ignal) and B(ackground). The weights are calculated using the sPlot procedure [arxiv.org/pdf/physics/0602023.pdf] and are a way of unfolding the two components so that I can analyze the S component only. $\endgroup$ – jwimberley Jul 9 '15 at 13:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.