nlcom for Risk Ratio

[Vicky Taxiarchi]

Hi, I run a univariate logistic regression model where the independent variable is a nominal with >2 categories. I use survey data and Stata MP version 16.1. My aim is to extract the risk ratio (RR) of each category to the reference and a global p-value for these.

I use the nlcom command for that, however the p-values following the nlcom are quite contradicting to the 95%CIs. Eg. for ratio 2: RR=1.04; 95%CI =(0.51, 1.57); p-value<0.001. Am I doing something wrong, or why this contradiction? And can I trust that the global p-value=0.351 is correct?

svy: logistic outcome i.cov
(running logistic on estimation sample)

note: 2.cov!= 0 predicts failure perfectly
      2.cov dropped and 13 obs not used

Survey: Logistic regression

Number of strata   =         1                  Number of obs     =        565
Number of PSUs     =       565                  Population size   = 552.706152
                                                Design df         =        564
                                                F(   2,    563)   =       0.51
                                                Prob > F          =     0.5987

-------------------------------------------------------------------------------
              |             Linearized
      outcome | Odds Ratio   Std. Err.      t    P>|t|     [95% Conf. Interval]
--------------+----------------------------------------------------------------
      cov|
       cov2 |          1  (empty)
       cov3 |   .5865927   .3142127    -1.00   0.320     .2048348    1.679846
cov4  |   1.063686   .4331284     0.15   0.880     .4780377    2.366816
              |
        _cons |   .5334486   .0505784    -6.63   0.000     .4428056    .6426463
-------------------------------------------------------------------------------
Note: _cons estimates baseline odds.

. margins cov, post

Adjusted predictions

Number of strata   =         1                  Number of obs     =        565
Number of PSUs     =       565                  Population size   = 552.706152
Model VCE    : Linearized                       Design df         =        564

Expression   : Pr(outcome), predict()

-------------------------------------------------------------------------------
              |            Delta-method
              |     Margin   Std. Err.      t    P>|t|     [95% Conf. Interval]
--------------+----------------------------------------------------------------
      cov |
       cov1  |   .3478751   .0215093    16.17   0.000     .3056269    .3901233
       cov3  |   .2383373   .0957039     2.49   0.013     .0503578    .4263168
cov4  |   .3620096   .0914604     3.96   0.000     .1823649    .5416542
-------------------------------------------------------------------------------


. nlcom (ratio1: _b[3.cov] / _b[1.cov]) (ratio2: _b[4.cov] / _b[1.cov]), post

      ratio1:  _b[3.cov] / _b[1.cov]
      ratio2:  _b[4.cov] / _b[1.cov]

------------------------------------------------------------------------------
             |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      ratio1 |   .6851231   .2783521     2.46   0.014      .139563    1.230683
      ratio2 |   1.040631   .2706705     3.84   0.000     .5101263    1.571135
------------------------------------------------------------------------------

. test _b[ratio1] = _b[ratio2]

 ( 1)  ratio1 - ratio2 = 0

           chi2(  1) =    0.87
         Prob > chi2 =    0.3510

[Andrew Musau]

the nlcom are quite contradicting to the 95%CIs. Eg. for ratio 2: RR=1.04; 95%CI =(0.51, 1.57); p-value<0.001.

What is contradictory about this?

[Daniel Schaefer]

It seems like the issue may be that in nlcom 1 is within the confidence interval for ratio 2, but these are relative risk ratios where 1 rather than 0 indicates no effect, and if 1 is within the confidence interval the ratio should not be significantly different from 1.

My understanding is that ncolm is a general command, and does not necessarily apply the correct transformation to the confidence intervals given that this is a logistic regression. This thread might provide some clues as to what is going on. Notice as well that the confidence interval is approximately the ratio +- (1.96 * the standard error). The coefficient and standard error should still be correct.

I also notice you don't use the survey prefix. I believe you should use -svy:- before nlcom here.

[Andrew Musau]

The null is that the coefficient is equal to zero, whether you define it as a risk ratio or otherwise. You are better off working with coefficients and then you take the exponents to determine the corresponding risk ratios. Or use margins that will allow you to specify the transformation.

help margins

[Maarten Buis]

You are better off working with coefficients and then you take the exponents to determine the corresponding risk ratios.

That would result in odds ratios and not risk ratios

[Andrew Musau]

Thanks. Shows that I am not familiar with these terminologies!

[Joseph Coveney]

My aim is to extract the risk ratio (RR) of each category to the reference and a global p-value for these.

If that's your aim, then why not just fit a risk-ratio regression model?

svy: glm outcome i.cov, family(binomial) link(log) eform

Unless you're having trouble with convergence that seems like a more straightforward approach.

[Tiago Pereira]

I think the suggestion offered by Joseph is way more elegant. It would be very interesting to compare it to this one, though:

https://journals.sagepub.com/doi/pdf...867X1101100208

-nlnom- should be employed in the ln scale. For example:

 nlcom (lnrr: ln(_b[1.low]/_b[0.low])), post

[Tiago Pereira]

* Example generated by -dataex-. For more info, type help dataex
clear
input byte(stage low died count)
1 0 0 50
1 0 1  5
1 1 0 10
1 1 1  2
2 0 0 57
2 0 1 17
2 1 0 13
2 1 1  9
3 0 0  6
3 0 1  9
3 1 0  2
3 1 1 12
end


expand count
logistic died i.low i.stage, nolog
margins low, post
nlcom (lnrr: ln(_b[1.low]/_b[0.low])), post eform
glm died i.low i.stage, fam(poisson) link(log) nolog vce(robust) eform

[Leonardo Guizzetti]

#9 is my preference for RR estimates since, using the log-binomial model, this rarely has convergence issues.

[Vicky Taxiarchi]

Many thanks to all for your suggestions. Tiago Pereira that was my mistake there indeed, I should have the logarithmic scale first. The results are very similar to the binomial model Joseph Coveney suggested, and quite close to the estimations of the poison model.