Generuj losowo skorelowane dane między zmienną binarną a ciągłą

23

Chcę wygenerować dwie zmienne. Jedna to zmienna wyniku binarnego (powiedzmy sukces / porażka), a druga to wiek w latach. Chcę, aby wiek był pozytywnie skorelowany z sukcesem. Na przykład powinno być więcej sukcesów w wyższych segmentach wiekowych niż w niższych. Idealnie powinienem być w stanie kontrolować stopień korelacji. W jaki sposób mogę to zrobić?

Dzięki

użytkownik333
źródło

Odpowiedzi:

20

Podejście @ ocram na pewno zadziała. Pod względem właściwości zależności jest to jednak nieco restrykcyjne.

Inną metodą jest użycie kopuły do ​​uzyskania wspólnego rozkładu. Możesz określić rozkład krańcowy dla sukcesu i wieku (jeśli masz istniejące dane, jest to szczególnie proste) oraz rodzinę kopuł. Różnicowanie parametrów kopuły da różne stopnie zależności, a różne rodziny kopuł dadzą różne zależności (np. Silna zależność górnej części ogona).

Najnowszy przegląd robienia tego w R za pomocą pakietu copula jest dostępny tutaj . Zobacz także dyskusję w tym dokumencie dla dodatkowych pakietów.

Jednak niekoniecznie potrzebujesz całego pakietu; oto prosty przykład z wykorzystaniem kopuły Gaussa, marginalnego prawdopodobieństwa sukcesu 0,6 i wieku rozproszenia gamma. Zmieniaj r, aby kontrolować zależność.

r = 0.8 # correlation coefficient
sigma = matrix(c(1,r,r,1), ncol=2)
s = chol(sigma)
n = 10000
z = s%*%matrix(rnorm(n*2), nrow=2)
u = pnorm(z)

age = qgamma(u[1,], 15, 0.5)
age_bracket = cut(age, breaks = seq(0,max(age), by=5))
success = u[2,]>0.4

round(prop.table(table(age_bracket, success)),2)

plot(density(age[!success]), main="Age by Success", xlab="age")
lines(density(age[success]), lty=2)
legend('topright', c("Failure", "Success"), lty=c(1,2))

Wydajność:

Stół:

           success
age_bracket FALSE TRUE
    (0,5]    0.00 0.00
    (5,10]   0.00 0.00
    (10,15]  0.03 0.00
    (15,20]  0.07 0.03
    (20,25]  0.10 0.09
    (25,30]  0.07 0.13
    (30,35]  0.04 0.14
    (35,40]  0.02 0.11
    (40,45]  0.01 0.07
    (45,50]  0.00 0.04
    (50,55]  0.00 0.02
    (55,60]  0.00 0.01
    (60,65]  0.00 0.00
    (65,70]  0.00 0.00
    (70,75]  0.00 0.00
    (75,80]  0.00 0.00

wprowadź opis zdjęcia tutaj

JMS
źródło
Świetna odpowiedź! Kopuły są pięknym, choć niedocenianym narzędziem. Model probitowy (z marginesem Gaussa na zmiennej ciągłej) jest szczególnym przypadkiem modelu kopuły Gaussa. Ale to jest znacznie bardziej ogólne rozwiązanie.
jpillow
1
@JMS: +1 Tak, Copulas są bardzo atrakcyjne. Powinienem spróbować przestudiować je bardziej szczegółowo!
ocram
@jpillow Rzeczywiście; Modele kopuły gaussowskiej obejmują wszelkiego rodzaju wielowymiarowe modele probitowe. Poprzez mieszanie skali obejmują one także kopie t / logistyczne i modele logit / robit. Tres cool :)
JMS
@ocram Do! Istnieje wiele otwartych pytań w mieszanych kontekstach danych (podczas używania ich jako modeli, a nie tylko czerpania z nich), które ludzie tacy jak ja chcieliby zobaczyć rozwiązani ...
JMS
@JMS Doskonała odpowiedź!
user333,
28

Możesz symulować model regresji logistycznej .

Dokładniej, możesz najpierw wygenerować wartości zmiennej wieku (na przykład stosując rozkład równomierny), a następnie obliczyć prawdopodobieństwo sukcesu za pomocą

π(x)=exp(β0+β1x)1+exp(β0+β1x)

β0β1β1

π

Przykład ilustrujący w R:

n <- 10
beta0 <- -1.6
beta1 <- 0.03
x <- runif(n=n, min=18, max=60)
pi_x <- exp(beta0 + beta1 * x) / (1 + exp(beta0 + beta1 * x))
y <- rbinom(n=length(x), size=1, prob=pi_x)
data <- data.frame(x, pi_x, y)
names(data) <- c("age", "pi", "y")
print(data)

         age        pi y
 1  44.99389 0.4377784 1
 2  38.06071 0.3874180 0
 3  48.84682 0.4664019 1
 4  24.60762 0.2969694 0
 5  39.21008 0.3956323 1
 6  24.89943 0.2988003 0
 7  51.21295 0.4841025 1
 8  43.63633 0.4277811 0
 9  33.05582 0.3524413 0
 10 30.20088 0.3331497 1
ocram
źródło
3
Nicea odpowiedź, choć z estetycznego punktu widzenia ( nie praktyczny jeden) model regresji probit może być jeszcze lepszy. Model probit jest równoważny z dwuwymiarowym RV Gaussa i progowaniem jednego z nich (do zera lub 1). Naprawdę polega to tylko na zastąpieniu funkcji skumulowanej normalnej Gaussa („probit”) logitem stosowanym w regresji logistycznej. Praktycznie powinno to dać taką samą wydajność (i obliczeniowo jest wolniejsze, ponieważ normcdf jest kosztowny do oceny (1 + e ^ x) ^ - 1), ale miło jest pomyśleć o Gaussa z jedną ze zmiennych ocenzurowanych („zaokrągloną”).
jpillow
@jpillow: Dziękujemy za komentarz. Zastanowię się jak najszybciej!
ocram
1
Zaletą modelu probula / kopuła Gaussa jest to, że parametry przyjmują formę macierzy kowariancji między dwiema wielkościami (z których jedna jest następnie binarna na 0 i 1). Jest to więc miłe z punktu widzenia interpretacji (ale nie tak miłe z punktu widzenia wygody obliczeniowej).
jpillow