Вот моя проблема дня:
В настоящий момент я обучаю себя эконометрике и использую логистическую регрессию. У меня есть код SAS, и я хочу быть уверенным, что хорошо его понимаю, прежде чем пытаться преобразовать его в R. (у меня нет и я не знаю SAS). В этом коде я хочу смоделировать вероятность того, что один человек окажется «безработным». Под этим я подразумеваю «возраст» от 15 до 64, а «такт» = «безработный». Я хочу попытаться предсказать этот результат с помощью следующих переменных: пол, возраст и idnat (номер национальности). (При прочих равных).
Код SAS:
/* Unemployment rate : number of unemployment amongst the workforce */
proc logistic data=census;
class sex(ref="Man") age idnat(ref="spanish") / param=glm;
class tact (ref=first);
model tact = sex age idnat / link=logit;
where 15<=age<=64 and tact in ("Employee" "Jobless");
weight weight;
format age ageC. tact $activity. idnat $nat_dom. inat $nationalty. sex $M_W.;
lsmeans sex / obsmargins ilink;
lsmeans idnat / obsmargins ilink;
lsmeans age / obsmargins ilink;
run;
Это пример того, как должна выглядеть база данных:
idnat sex age tact
[1,] "english" "Woman" "42" "Employee"
[2,] "french" "Woman" "31" "Jobless"
[3,] "spanish" "Woman" "19" "Employee"
[4,] "english" "Man" "45" "Jobless"
[5,] "english" "Man" "34" "Employee"
[6,] "spanish" "Woman" "25" "Employee"
[7,] "spanish" "Man" "39" "Jobless"
[8,] "spanish" "Woman" "44" "Jobless"
[9,] "spanish" "Man" "29" "Employee"
[10,] "spanish" "Man" "62" "Retired"
[11,] "spanish" "Man" "64" "Retired"
[12,] "english" "Woman" "53" "Jobless"
[13,] "english" "Man" "43" "Jobless"
[14,] "french" "Man" "61" "Retired"
[15,] "french" "Man" "50" "Employee"
Вот какой результат я хочу получить:
Variable Modality Value ChiSq Indicator
Sex Women 56.6% 0.00001 -8.9%
Men 65.5%
Nationality
1:Spanish 62.6%
2:French 51.2% 0.00001 -11.4%
3:English 48.0% 0.00001 -14.6%
Age
<25yo 33.1% 0.00001 -44.9%
Ref:26<x<54yo 78.0%
55yo=< 48.7% 0.00001 -29.3%
(Я интерпретирую вышесказанное следующим образом: при прочих равных, вероятность трудоустройства женщин составляет -8,9% по сравнению с мужчинами, а вероятность трудоустройства у женщин младше 25 лет составляет -44,9%, чем у женщин в возрасте от 26 до 54 лет).
Итак, если я хорошо понимаю, лучшим подходом было бы использование бинарной логистической регрессии (ссылка = логит). Здесь используются ссылки «мужчина против женщины» (пол), «работник против безработного» (из переменной 'такта') ... Я предполагаю, что 'такт' автоматически преобразуется SAS в двоичную (0-1) переменную.
Вот моя первая попытка в R. Я еще не проверял (нужен собственный компьютер):
### before using multinom function
### change all predictors to factors and relevel
recens$sex <- relevel(factor(recens$sex), ref = "Man")
recens$idnat <- relevel(factor(recens$idnat), ref = "spanish")
recens$TACT <- relevel(factor(recens$TACT), ref = "employee")
### Calculations of the probabilities with function multinom,
### formatted variables, and conditions with subset
glm1 <- glm(TACT ~ sex + age + idnat, data=census,
+ weights = weight, subset=age[(15<=recens$age|recens$age<=64)] & TACT %in%
+ c("Employee","Jobless"), family=binomial())
Мои вопросы :
На данный момент кажется, что есть много функций для выполнения логистической регрессии в R, например glm, что кажется подходящим.
Однако после посещения многих форумов кажется, что многие люди рекомендуют не пытаться точно воспроизвести SAS PROC LOGISTIC, особенно функцию LSMEANS functions. Доктор Франк Харрел (автор package:rms), например.
Тем не менее, я думаю, что моя большая проблема - это LSMEANS и его варианты Obsmargins и ILINK. Даже после многократного перечитывания его описания я с трудом могу понять, как он работает.
Пока что я понимаю Obsmargin, так это то, что он учитывает структуру всего населения базы данных (т. Е. Вычисления выполняются с пропорциями от общего населения). ILINK, по-видимому, используется для получения прогнозируемого значения вероятности (уровня безработицы, уровня занятости) для каждого из предикторов (например, женщины, а затем мужчины), а не значения, найденного с помощью (экспоненциальной) модели?
Короче говоря, как это можно сделать с помощью R с rms функциями, такими как lrm?
Я действительно потерялся во всем этом. Если бы кто-то мог объяснить мне это лучше и сказать, на правильном ли я пути, это сделало бы мой день лучше.
Спасибо за вашу помощь и извините за все ошибки, мой английский немного ржавый.
Бинь
