0
$\begingroup$

I'm trying to model the probability that somebody can do a 90 day Snapchat streak. I have some binary success data in the form of count data (the participant only reached 90 days once):

1

I've converted this count data into the probability that the person will successfully reach 90 days in future, based on the day they have reached:

2

The fitted curve was created using nls() in R, using the code below for a logistic growth curve:

# Logistic growth model
model5 <- nls(`Probability of success` ~ a/(1+exp(-(b+c*Day.number))),
          data=prob_90,
          start=list(a=1, b=1, c=1),
          upper=list(a=1), # Set upper limit on a
          algorithm='port')

However, I've since realised that because the count data is poisson or gamma distributed, a symmetrical function like the logistic curve won't work. It appears the CDF of the poisson distribution would work, however this includes a factorial term and the gamma function, which returns errors in nls().

Or am I doing this completely backwards? Should I be trying to model the count data in Figure 1, then converting that to a curve for Figure 2? Or some other method?

Thanks very much!

Improve this question
$\endgroup$
4
  • $\begingroup$ Is there a reason you don't want to use a GLM in this case? Poisson regression seems much more straightforward than nls. In fact, if you don't have any predictors and just the observed streak days, you may as well just determine by maximum likelihood the parameter(s) of the distribution (e.g. using MASS::fitdistr). If this is the case, it might be helpful to share the data as well to create a reprex. $\endgroup$
    – awhug
    Jan 10, 2022 at 4:32
  • $\begingroup$ Thank you, that was exactly the function I needed. So using fitdistr() I was able to model the distribution almost perfectly with a negative binomial model, but I've realised that I can't just convert it to a CDF to model the probabilities in the second graph. Is this where I need to switch to a GLM or something similar? Part of the problem is as you say, I don't have any predictors at the moment so it's just a pure distribution. $\endgroup$
    – Piethon
    Jan 10, 2022 at 7:09
  • $\begingroup$ It seems unlikely that doing another day is independent of the days done so far, but if it was then you could simply estimate that the probability of failure on a particular day is something like $0.133$. That does not fit the rate of immediate failure or give a substantial probability of reaching $90$ so you may prefer another model $\endgroup$
    – Henry
    Jan 10, 2022 at 9:52
  • $\begingroup$ Ah, this is true. Is there any alternative to the negative binomial distribution in cases where each trial is dependent? $\endgroup$
    – Piethon
    Jan 10, 2022 at 20:11

1 Answer 1

Reset to default
1
$\begingroup$

I don't think this is the right approach, but you have said something interesting

the probability that the person will successfully reach 90 days in future, based on the day they have reached

This is reminiscent of survival analysis. It doesn't seem like you have any other variables save time, so it might be worth just fitting a Kaplan Meir estimator in R to get the estimated survival function. I've gone ahead and copied your data and run the analysis myself

library(tidyverse)
library(survival)
days = c(0,12,3,4,5,7,8,9,10,15,17,18,21,28,29,34,90)
counts = c(21,6,6,5,12,2,1,1,2,1,1,1,1,1,1,1,1)
status = rep(1, length(days))
status[length(status)]=0

d = tibble(days, counts, status) %>% 
    uncount(counts)

km_fit = survfit(Surv(days, status)~1, data=d)

plot(km_fit, xlab = 'Days', ylab='Survival Prob')

enter image description here

This curve shows the probability of the streak will end after $t$ days. We can print out the predictions easily, which give the probability that the streak will end after the indicated time plus a confidence interval for the prediction.

> summary(km_fit, times=c(1:5, 10, 30, 90))
Call: survfit(formula = Surv(days, status) ~ 1, data = d)

 time n.risk n.event survival std.err lower 95% CI upper 95% CI
    1     43      21   0.6719  0.0587      0.56615        0.797
    2     43       0   0.6719  0.0587      0.56615        0.797
    3     43       6   0.5781  0.0617      0.46895        0.713
    4     37       5   0.5000  0.0625      0.39135        0.639
    5     32      12   0.3125  0.0579      0.21729        0.449
   10     16       6   0.2188  0.0517      0.13768        0.348
   30      2      12   0.0312  0.0217      0.00799        0.122
   90      1       1   0.0156  0.0155      0.00224        0.109

As you can see, the estimated risk very quickly goes to zero, meaning it is very unlikely many people will make it to 90 days.

Improve this answer
$\endgroup$
4
  • $\begingroup$ Thank you for that answer, this is very useful! I did consider survival analysis but haven't done it before, so this is useful to get started. The one thing I am still wondering is whether I can fit a general curve to the distribution in graph 2, perhaps using a negative binomial model or something similar (see awhug's comment thread). From what I've seen, survival curve analysis is excellent for creating an accurate non-parametric model, but it's not great for determining a general rule that fits into a nice mathematical form. Or is there a way to model survival curves? $\endgroup$
    – Piethon
    Jan 10, 2022 at 7:14
  • $\begingroup$ You have a couple of typos in your data. Perhaps days=c(0,1,2,3,4,5,7,8,9,10,15,17,18,21,28,29,34,90) and counts=c(21,6,6,4,1,2,2,1,1,2,1,1,1,1,1,1,1,1) $\endgroup$
    – Henry
    Jan 10, 2022 at 9:46
  • $\begingroup$ @Piethon Why do you have a preference for mathematics when a non-parametric estimator gives you something closer to what you actually want (i.e. the probability of the streak being broken before time $t$). $\endgroup$
    – Demetri Pananos
    Jan 10, 2022 at 14:53
  • $\begingroup$ I guess in my mind it's cleaner? I come from a physics background so it's always nice to be able to fit a general mathematical model to something. However you're right, the non-parametric model will do the job in a pinch. $\endgroup$
    – Piethon
    Jan 10, 2022 at 20:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.