Quoting from Harrell's [_Regression Modeling Strategies_](http://www.springer.com/us/book/9783319194240), second edition, page 494:

> When correlations among predictors are mild, plots of estimated predictor transformations without adjustment for other predictors (i.e., marginal transformations) may be useful. Martingale residuals may be obtained quickly by fixing $\hat \beta$ = 0 for all predictors. Then smoothed plots of predictor against residual may be made for all predictors.

So there is nothing necessarily wrong with examining martingale residuals from a null model for linearity, provided that "correlations among predictors are mild." As noted in Table 20.3 on page 494, there are several acceptable ways to use martingale residuals in this context, depending on whether you want to estimate transformations versus checking for nonlinearity, and whether you wish to adjust for other predictors in the process.

The apparent differences among the cited references represent which way the martingale residuals were calculated and used.

If you display residuals for a model _null in the continuous predictor_ of interest, then you will display the _shape of the relation of the predictor to outcome_. So if that's a straight line you have a linear relation, as in the first example you cite. If you don't get a straight line, the shape of the curve might suggest a useful functional form to try for that predictor. That's what you do when you try to _estimate_ the functional form of the relation.

If instead you check the relation of the martingale residuals to the predictor from a model _including the estimated coefficient for that predictor_, then you hope for a _flat horizontal_ smoothed line, indicating that the single linear coefficient for the predictor is adequate. That's what you do when you _test for_ linearity in the predictor.

Nevertheless, particularly with this many data points, starting with restricted spline fits for continuous predictors provides a simple and general way to both test for and adjust for nonlinearities. If the relation of a predictor to outcome is linear, then the higher-order spline terms will have non-significant coefficients. If the relation is non-linear, you can adjust for any reasonable degree of nonlinearity by adding more knots. You ultimately should calibrate the full model (checking linearity with respect to the overall linear predictor, and correcting for optimism), as with the `calibrate()` function in Harrell's `rms` package in R.