Obliczanie wielkości próby dla jednoczynnikowej regresji logistycznej

Jak obliczyć wielkość próby potrzebną do badania, w którym grupa badanych będzie miała jedną ciągłą zmienną zmierzoną w czasie operacji, a następnie dwa lata później zostaną one zaklasyfikowane jako wynik funkcjonalny lub wynik pogorszony.

Chcielibyśmy sprawdzić, czy pomiar ten mógł przewidzieć zły wynik. W pewnym momencie możemy chcieć uzyskać punkt odcięcia w zmiennej ciągłej, powyżej którego próbujemy interweniować, aby zmniejszyć prawdopodobieństwo pogorszenia wyniku.

Jakieś pomysły? Dowolna implementacja R.

Obliczenia wielkości próby regresji logistycznej są złożone. Nie będę tutaj próbował tego streszczać. Racjonalnie dostępne rozwiązania tego problemu można znaleźć w:

Hsieh FY. Przykładowe tabele rozmiarów dla regresji logistycznej. Statystyka w medycynie. 1989 lipiec; 8 (7): 795-802.

Hsieh FY i in. Prosta metoda obliczania wielkości próby dla regresji liniowej i logistycznej. Statystyka w medycynie. 30 lipca 1998; 17 (14): 1623–34.

Dostępną dyskusję na temat zagadnień z przykładowymi obliczeniami można znaleźć w ostatnim rozdziale (Rozdział 8.5 str. 339-347) Stosowanej regresji logistycznej Hosmer & Lemeshow .

Zwykle prowadzenie symulacji jest dla mnie łatwiejsze i szybsze. Artykuły długo czytają, rozumieją i ostatecznie dochodzą do wniosku, że nie mają zastosowania w szczególnym przypadku, który jest zainteresowany.

Dlatego wybrałbym tylko kilka tematów, symulowałem zmienną towarzyszącą, którą jesteś zainteresowany (rozłożony tak, jak ci się wydaje,), symulował dobre / złe wyniki w oparciu o postawioną funkcjonalną formę (efekty progowe współzmiennej? Nieliniowość?) z minimalnym (klinicznie) znaczącym rozmiarem efektu, który chcesz wykryć, przeprowadź wynik przez analizę i sprawdź, czy efekt zostanie znaleziony na poziomie alfa. Uruchom to ponownie 10 000 razy i sprawdź, czy uzyskałeś efekt w 80% symulacji (lub jakiejkolwiek innej mocy, której potrzebujesz). Dostosuj liczbę tematów, powtarzaj, aż uzyskasz moc, z której jesteś zadowolony.

Ma to tę zaletę, że jest bardzo ogólny, więc nie jesteś ograniczony do określonej formy funkcjonalnej lub określonej liczby lub rozkładu zmiennych towarzyszących. Możesz dołączyć przerywniki, patrz komentarz chl powyżej, losowo lub pod wpływem zmiennej towarzyszącej lub wyniku. Zasadniczo kodujesz analizę, którą zamierzasz wykonać na końcowej próbce, co czasem pomaga skupić moje myślenie na projekcie badania. I łatwo to zrobić w R (wektoryzacja!).

Kontynuując post Stephana Kolassy (nie mogę tego dodać jako komentarza), mam alternatywny kod do symulacji. Wykorzystuje tę samą podstawową strukturę, ale eksploduje nieco bardziej, więc być może jest trochę łatwiejszy do odczytania. Opiera się również na kodzie Kleinmana i Hortona do symulacji regresji logistycznej.

nn jest liczbą w próbce. Zmienna towarzysząca powinna być stale rozkładana normalnie i znormalizowana do wartości 0 i sd 1. Do wygenerowania tego używamy rnorm (nn). Wybieramy iloraz szans i przechowujemy go w nieparzystych. Ratio. Wybieramy również numer do przechwytywania. Wybór tej liczby decyduje o tym, jaka część próby doświadcza „zdarzenia” (np. 0,1, 0,4, 0,5). Musisz grać z tym numerem, aż uzyskasz odpowiednią proporcję. Poniższy kod podaje proporcję 0,1 przy wielkości próbki 950 i OR 1,5:

nn <- 950
runs <- 10000
intercept <- log(9)
odds.ratio <- 1.5
beta <- log(odds.ratio)
proportion  <-  replicate(
              n = runs,
              expr = {
                  xtest <- rnorm(nn)
                  linpred <- intercept + (xtest * beta)
                  prob <- exp(linpred)/(1 + exp(linpred))
                  runis <- runif(length(xtest),0,1)
                  ytest <- ifelse(runis < prob,1,0)
                  prop <- length(which(ytest <= 0.5))/length(ytest)
                  }
            )
summary(proportion)

podsumowanie (proporcja) potwierdza, że ​​proporcja wynosi ~ 0,1

Następnie przy użyciu tych samych zmiennych moc oblicza się na 10000 przebiegów:

result <-  replicate(
              n = runs,
              expr = {
                  xtest <- rnorm(nn)
                  linpred <- intercept + (xtest * beta)
                  prob <- exp(linpred)/(1 + exp(linpred))
                  runis <- runif(length(xtest),0,1)
                  ytest <- ifelse(runis < prob,1,0)
                  summary(model <- glm(ytest ~ xtest,  family = "binomial"))$coefficients[2,4] < .05
                  }
            )
print(sum(result)/runs)

Myślę, że ten kod jest poprawny - porównałem go z przykładami podanymi w Hsieh, 1998 (tabela 2) i wydaje się zgadzać z trzema podanymi tam przykładami. Przetestowałem go również na przykładzie na str. 342 - 343 Hosmer i Lemeshow, gdzie znaleziono moc 0,75 (w porównaniu do 0,8 w Hosmer i Lemeshow). Być może w niektórych okolicznościach takie podejście nie docenia władzy. Jednak kiedy uruchomiłem ten sam przykład w tym kalkulatorze internetowym , okazało się, że zgadza się on ze mną, a nie z wynikami w Hosmer i Lemeshow.

Jeśli ktoś mógłby nam powiedzieć, dlaczego tak jest, byłbym zainteresowany.

$\theta = 1$$_0: \theta = 0$ . wydaje się, że nie określasz żadnego kryterium wyboru wielkości próby.

w rzeczywistości brzmi to tak, jakby twoje badania były prowadzone sekwencyjnie. w takim przypadku opłaca się uczynić z tego wyraźną część eksperymentu. sekwencyjne próbkowanie może często być bardziej wydajne niż eksperyment o ustalonej wielkości próby [średnio potrzeba mniej obserwacji].

farrel: dodaję to w odpowiedzi na twój komentarz.

aby uzyskać wielkość próby, zwykle określa się jakieś kryterium precyzji dla oszacowania [takiego jak długość CI] LUB mocy przy określonej alternatywie testu, który należy przeprowadzić na danych. zdaje się, że wymieniłeś oba te kryteria. w zasadzie nie ma w tym nic złego: wystarczy wykonać dwa obliczenia wielkości próbki - jedno dla uzyskania pożądanej dokładności oszacowania - i drugie dla uzyskania pożądanej mocy przy podanej alternatywie. wymagany jest większy z dwóch rozmiarów próbek. [btw - inaczej niż mówiąc o 80% mocy - wydaje się, że nie wspomniałeś o tym, jaki test planujesz wykonać - ani o alternatywie, przy której chcesz 80% mocy.]

jeśli chodzi o stosowanie analizy sekwencyjnej: jeśli badani są zapisani do badania jednocześnie, wówczas ustalony rozmiar próby ma sens. ale jeśli przedmiotów jest niewielu i to daleko od siebie, uzyskanie wymaganej liczby może zająć rok lub dwa [lub więcej]. w ten sposób proces może trwać trzy lub cztery lata [lub więcej]. w takim przypadku schemat sekwencyjny oferuje możliwość zatrzymania się wcześniej - jeśli efekt, którego szukasz, stanie się statystycznie znaczący na wcześniejszym etapie badania.