# Logistic regression with imperfect measurement of the outcome

**URL:** <https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804>\
**Category:** Modeling\
**Tags:** specification\
**Created:** [September 4, 2019, 10:48pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804 "2019-09-04T22:48:21Z")\
**Posts on this page:** 12\
**Page:** 1

<div class="post-metadata">

**Author:** ![Felix\_Fischer](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/felix_fischer/32/5787_2.png) [@Felix\_Fischer](https://discourse.mc-stan.org/u/Felix_Fischer)\
**Post date:** [September 4, 2019, 10:48pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/1 "2019-09-04T22:48:21Z")

</div>

Hi,

i want to use Stan to fit a logistic regression which accounts for imperfect measurement of the outcome. The paper describing the Bayesian approach (McInturff 2004) uses discrete latent parameters in bugs. In order to perform marginalization i turned to an earlier paper describing an EM algorithm to estimate this kind of regression (Magder & Hughes, 1997). This paper gives the conditional probabilities, but my model seems to be not correct as i can’t get reasonable regression estimates even when i impose high sensitivity and specificity.

I think my error could be where I sum up the target density, but frankly i am lost. Anyone who spots my misdoing here?

best wishes, felix

```stan

    data{
      int n;
      int outcome[n]; //outcome measured with error
      int gender[n]; //a predictor
    }
      
    parameters {
      real<lower=0,upper=1> sensitivity;
      real<lower=0,upper=1> specificity;
    
      real slope;
      real intercept;
    }
    
    transformed parameters {
      real linear_predictor[n];
      real pie[n];
      real likelihood[n];
      
      for (i in 1:n){
        linear_predictor[i] = slope*gender[i]+intercept;
         pie[i] = exp(linear_predictor[i])/(1+exp(linear_predictor[i])); //probability 
      
        //marginalization --> different probabilities for true disease state given observed outcome and linear predictor
        if (outcome[i] == 0){          
          likelihood[i] = log(pie[i]*sensitivity/(pie[i]*sensitivity+pie[i]*(1-specificity)));
        } else if (outcome[i] == 1){
          likelihood[i] = log(pie[i]*(1-sensitivity)/(pie[i]*(1-sensitivity)+pie[i]*specificity));
        }
      }
    }         

model {
      sensitivity ~ beta(200,2); 
      specificity ~ beta(200,2);
      
      intercept ~ normal(0,1);
      slope ~ normal(0,1);
      
      target += log_sum_exp(likelihood);
    }

```

---

<div class="post-metadata">

**Author:** ![mike-lawrence](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/mike-lawrence/32/59_2.png) [@mike-lawrence](https://discourse.mc-stan.org/u/mike-lawrence)\
**Post date:** [September 5, 2019, 2:34pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/2 "2019-09-05T14:34:50Z")

</div>

Here’s the model description from McInturff 2004:

 ![Screenshot%20from%202019-09-05%2011-27-54](https://canada1.discourse-cdn.com/flex030/uploads/mc_stan/original/2X/2/22860ef7a87bb46f6a68a4138c02e52a9b59aa9b.png)

So, I think you have the two likelihoods flipped. Plus, possibly the algebra works out the same, but copying directly from McInturff I get:

```
        if (outcome[i] == 1){          
          likelihood[i] = log( pie[i]*sensitivity + (1-pie[i]) * (1-specificity) );
        } else if (outcome[i] == 0){
          likelihood[i] = log( pie[i]*(1-sensitivity) + (1-pie[i])*specificity );
        }

```

---

<div class="post-metadata">

**Author:** ![Felix\_Fischer](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/felix_fischer/32/5787_2.png) [@Felix\_Fischer](https://discourse.mc-stan.org/u/Felix_Fischer)\
**Post date:** [September 5, 2019, 6:10pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/3 "2019-09-05T18:10:42Z")

</div>

Thanks for looking into this. Pasting your code unfortunately still does not give reasonable parameters.

What I am calculating there (according to Magder / Hughes) is the probability of true disease state given the observed test result and the predictors (i should have not named that variable likelihood). I think this information is needed to weight the likelihood contributions.

I ultimately want to model is the relation between the predictors and the true disease state and i feel that part is still missing? Should I calculate the likelihood now for both states (true diseases yes/no) and then sum it up weighted by the probabilities?

 ![hughes](https://canada1.discourse-cdn.com/flex030/uploads/mc_stan/original/2X/4/41613206af723ad7b24116fc76fa48a416986005.png)

---

<div class="post-metadata">

**Author:** ![LucC](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/lucc/32/4475_2.png) [@LucC](https://discourse.mc-stan.org/u/LucC)\
**Post date:** [September 5, 2019, 6:47pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/4 "2019-09-05T18:47:20Z")

</div>

I can’t readily follow the math in your likelihood definition. But I do have some remarks that may help.

The expected probability p of observing outcome ==1 is: p = sens x p\_true + (1-spec) x (1-p\_true). From there it should be straightforward to work out the log-likelihood. Maybe your definition is right, but I don’t understand how you end up with a multiplication of the specificity and true prevalence in your likelihood…

General tip: define the inverse logit transformation as 1 / (1+exp(-(linear predictor))). Thats fewer terms. Not sure if it makes a difference for the autodiff cost though.

---

<div class="post-metadata">

**Author:** ![LucC](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/lucc/32/4475_2.png) [@LucC](https://discourse.mc-stan.org/u/LucC)\
**Post date:** [September 5, 2019, 6:52pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/5 "2019-09-05T18:52:24Z")

</div>

Also, if you work out the likelihood based on the expected probability that I defined above, you have already marginalised out all the possible latent states and dot’t need to do the log\_sum\_exp thing. I don’t think you needed that in the firste place…

---

<div class="post-metadata">

**Author:** ![mcol](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/mcol/32/5598_2.png) [@mcol](https://discourse.mc-stan.org/u/mcol)\
**Post date:** [September 5, 2019, 6:59pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/6 "2019-09-05T18:59:30Z")

</div>

> [@LucC](#):
>
> General tip: define the inverse logit transformation as 1 / (1+exp(-(linear predictor))). Thats fewer terms. Not sure if it makes a difference for the autodiff cost though.

Wouldn’t be better to use `inv_logit()`?

---

<div class="post-metadata">

**Author:** ![mike-lawrence](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/mike-lawrence/32/59_2.png) [@mike-lawrence](https://discourse.mc-stan.org/u/mike-lawrence)\
**Post date:** [September 5, 2019, 6:59pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/7 "2019-09-05T18:59:57Z")

</div>

Good catch; I didn’t look closely enough at the model section. Since `likelihood` is a vector filled with log-likeilhoods (might consider renaming it to `log_lik` for clarity), you don’t need to exponentiate then sum then log to get the sum log likelihood, just sum: `target += sum(likelihood)`

---

<div class="post-metadata">

**Author:** ![Felix\_Fischer](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/felix_fischer/32/5787_2.png) [@Felix\_Fischer](https://discourse.mc-stan.org/u/Felix_Fischer)\
**Post date:** [September 5, 2019, 7:12pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/8 "2019-09-05T19:12:33Z")

</div>

Thanks so much, i think i got it now!

```stan
   data{
     int n;
     int outcome[n];
     int gender[n];
   }
     
   parameters {
     real<lower=0,upper=1> sensitivity;
     real<lower=0,upper=1> specificity;
   
     real slope;
     real intercept;
   }
   
   transformed parameters {
     real linear_predictor[n];
     real pie[n];
     real loglik[n];

     
     for (i in 1:n){
       linear_predictor[i] = slope*gender[i]+intercept;
       pie[i] = 1/(1+exp(-(linear_predictor[i])));
     
       if (outcome[i] == 1){          
         loglik[i] = log( pie[i]*sensitivity + (1-pie[i]) * (1-specificity) );
       } else if (outcome[i] == 0){
         loglik[i] = log( pie[i]*(1-sensitivity) + (1-pie[i])*specificity );
       }
}
}

    
   
   model {
     sensitivity ~ beta(200,2);
     specificity ~ beta(200,2);
     
     intercept ~ normal(0,1);
     slope ~ normal(0,1);
     
    target += sum(loglik);
   }

```

---

<div class="post-metadata">

**Author:** ![mike-lawrence](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/mike-lawrence/32/59_2.png) [@mike-lawrence](https://discourse.mc-stan.org/u/mike-lawrence)\
**Post date:** [September 5, 2019, 7:17pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/9 "2019-09-05T19:17:55Z")

</div>

Be sure to double check by simulating some fake data and ensuring the simulation parameters are recovered reasonably.

---

<div class="post-metadata">

**Author:** ![Felix\_Fischer](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/felix_fischer/32/5787_2.png) [@Felix\_Fischer](https://discourse.mc-stan.org/u/Felix_Fischer)\
**Post date:** [September 5, 2019, 7:19pm UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/10 "2019-09-05T19:19:21Z")

</div>

Already on it :D

---

<div class="post-metadata">

**Author:** ![LucC](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/lucc/32/4475_2.png) [@LucC](https://discourse.mc-stan.org/u/LucC)\
**Post date:** [September 6, 2019, 7:49am UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/11 "2019-09-06T07:49:09Z")

</div>

Sure yes, if that’s an existing function in Stan, that’s probably even more stable and efficient

---

<div class="post-metadata">

**Author:** ![mcol](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/mcol/32/5598_2.png) [@mcol](https://discourse.mc-stan.org/u/mcol)\
**Post date:** [September 6, 2019, 10:17am UTC](https://discourse.mc-stan.org/t/logistic-regression-with-imperfect-measurement-of-the-outcome/10804/12 "2019-09-06T10:17:05Z")

</div>

Yes, one could write

```
pie[i] = inv_logit(linear_predictor[i]);

```

I don’t think it’s vectorized, so the loop over elements will have to stay.