Jakie jest prawdopodobieństwo tego procesu?

Pacjent zostaje przyjęty do szpitala. Ich długość pobytu zależy od 2 rzeczy: ciężkości urazu i wysokości ubezpieczenia, jaką jest skłonny zapłacić, aby zatrzymać je w szpitalu. Niektórzy pacjenci odejdą przedwcześnie, jeśli ich ubezpieczenie zdecyduje się przestać płacić za pobyt.

Załóż, że:

1) Długość pobytu jest rozłożona poissonem (załóżmy, że może to być realistyczne założenie), z parametrem .$\lambda$

2) Różne plany ubezpieczenia obejmują pobyty na 7, 14 i 21 dni. Wielu pacjentów odejdzie po 7,14 lub 21 dniach pobytu (ponieważ kończy się ich ubezpieczenie i muszą wyjechać).

Gdybym miał uzyskać dane z tego procesu, mogłoby to wyglądać następująco:

Jak widać, w punkcie 7, 14 i 21 dnia występują skoki. Są to pacjenci, którzy wyjeżdżają po wygaśnięciu ubezpieczenia.

Oczywiście dane można modelować jako mieszaninę. Trudno mi zapisać prawdopodobieństwo tej dystrybucji. To jest jak poison z napompowaniem zerowym, ale inflacja wynosi 7, 14 i 21.

Jakie jest prawdopodobieństwo tych danych? Jaki jest proces myślowy leżący u podstaw prawdopodobieństwa?

W tym przypadku uważam, że istnieje ścieżka do rozwiązania, jeśli założymy nasz kapelusz analizy przeżycia. Zauważ, że nawet jeśli w tym modelu nie ma ocenzurowanych podmiotów (w tradycyjnym znaczeniu), nadal możemy korzystać z analizy przeżycia i mówić o zagrożeniach ze strony podmiotów.

Musimy modelować trzy rzeczy w tej kolejności: i) skumulowane zagrożenie, ii) zagrożenie, iii) prawdopodobieństwo dziennika.

i) Zrobimy część i) etapami. Jakie jest skumulowane zagrożenie zmiennej losowej Poissona? W przypadku dystrybucji dyskretnej istnieją dwa sposoby jej zdefiniowania¹, ale użyjemy definicji . Tak więc skumulowane zagrożenie dla wynosi$H(t)$$H(t) = -\log{S(t)}$$T \sim Poi(\lambda)$

$$H_T(t) = -\log{(1 - Q(t, \lambda))} = -\log{P(t, \lambda)}$$

gdzie jest odpowiednio górną, dolną regularyzowaną funkcją gamma.$Q, P$

Teraz chcemy dodać „zagrożenia” związane z brakującym ubezpieczeniem. Zaletą kumulatywnych zagrożeń jest to, że są one addytywne, dlatego musimy po prostu dodać „ryzyko” w czasach 7, 14, 21:

$$H_{T'}(t) = -\log{P(t, \lambda)} + a\cdot\mathbb{1}_{(t>7)} + b\cdot\mathbb{1}_{(t>14)} + c\cdot\mathbb{1}_{(t>21)}$$

Heurystycznie pacjent jest narażony na podstawowe ryzyko „Poissona”, a następnie ryzyko punktowe w punktach 7, 14 i 21. (Ponieważ jest to ryzyko skumulowane , kumulujemy to ryzyko punktowe, stąd ). nie wiemy, co to są i , ale później połączymy je z prawdopodobieństwem braku ubezpieczenia.$>$$a, b$$c$

W rzeczywistości, ponieważ wiemy, że 21 jest górną granicą, a następnie wszyscy pacjenci są usuwani, możemy ustawić na nieskończoność.$c$

$$H_{T'}(t) = -\log{P(t, \lambda)} + a\cdot\mathbb{1}_{(t>7)} + b\cdot\mathbb{1}_{(t>14)} + \infty \cdot\mathbb{1}_{(t>21)}$$

ii) Następnie wykorzystujemy skumulowane zagrożenie, aby uzyskać zagrożenie, . Wzór na to jest następujący:$h(t)$

$$h(t) = 1 - \exp{(H(t) - H(t+1))}$$

Podłączenie naszego skumulowanego zagrożenia i uproszczenie:

$$h_{T'}(t) = 1 - \frac{P(t+1, \lambda)}{P(t, \lambda)} \exp(-a\cdot\mathbb{1}_{(t=7)} - b\cdot\mathbb{1}_{(t=14)} - \infty \cdot\mathbb{1}_{(t=21)})$$

iii) Wreszcie, napisanie prawdopodobieństwa dziennika dla modeli przeżycia (bez cenzury) jest bardzo łatwe, gdy mamy zagrożenie i ryzyko skumulowane:

$$ll(\lambda, a, b \;|\; t) = \sum_{i=1}^N \left(\log h(t_i) - H(t_i)\right)$$

I oto jest!

Istnieją zależności, które łączą nasze punktowe współczynniki ryzyka i prawdopodobieństwa długości ubezpieczenia: .$a = -\log(1 - p_a), b = -\log(1 - p_a - p_b) - \log(1 - p_a), p_c = 1 - (p_a + p_b)$

---

Dowód jest w budyniu. Zróbmy kilka symulacji i wnioskowania przy użyciu semantyki niestandardowego modelu linii życia .

from lifelines.fitters import ParametericUnivariateFitter
from autograd_gamma import gammaincln, gammainc
from autograd import numpy as np

MAX = 1e10

class InsuranceDischargeModel(ParametericUnivariateFitter):
    """
    parameters are related by
    a = -log(1 - p_a)
    b = -log(1 - p_a - p_b) - log(1 - p_a)
    p_c = 1 - (p_a + p_b)
    """
    _fitted_parameter_names = ["lbd", "a", "b"]
    _bounds = [(0, None), (0, None), (0, None)]

    def _hazard(self, params, t):
        # from (1.64c) in http://geb.uni-giessen.de/geb/volltexte/2014/10793/pdf/RinneHorst_hazardrate_2014.pdf
        return 1 - np.exp(self._cumulative_hazard(params, t) - self._cumulative_hazard(params, t+1))

    def _cumulative_hazard(self, params, t):
        lbd, a, b = params
        return -gammaincln(t, lbd) + a * (t > 7) + b * (t > 14) + MAX * (t > 21)


def gen_data():
    p_a, p_b = 0.4, 0.2
    p = [p_a, p_b, 1 - p_a - p_b]
    lambda_ = 18
    death_without_insurance = np.random.poisson(lambda_)
    insurance_covers_until = np.random.choice([7, 14, 21], p=p)
    if death_without_insurance < insurance_covers_until:
        return death_without_insurance
    else:
        return insurance_covers_until


durations = np.array([gen_data() for _ in range(40000)])
model = InsuranceDischargeModel()
model.fit(durations)
model.print_summary(5)
"""
<lifelines.InsuranceDischargeModel: fitted with 40000 observations, 0 censored>
number of subjects = 40000
  number of events = 40000
    log-likelihood = -78845.10392
        hypothesis = lbd != 1, a != 1, b != 1

---
        coef  se(coef)  lower 0.95  upper 0.95      p  -log2(p)
lbd 18.05026   0.03353    17.98455    18.11598 <5e-06       inf
a    0.50993   0.00409     0.50191     0.51794 <5e-06       inf
b    0.40777   0.00557     0.39686     0.41868 <5e-06       inf
"""

---

¹ patrz sekcja 1.2 tutaj