Régression probit contrainte dans R
Je cherche à exécuter un modèle probit dans R définissant certains coefficients égaux les uns aux autres.
Prenons l'exemple simple où quatre équipes s'affrontent une fois à domicile et une fois sur la route:
Home <- c('NY','NY','NY','LA','LA','LA','BOS','BOS','BOS','CHI','CHI','CHI')
Away <- c('LA','CHI','BOS','NY','CHI','BOS','LA','CHI','NY','LA','NY','BOS')
HomeWin <- c(1,1,0,1,0,1,0,1,0,0,0,1)
results <- data.frame(Home,Away,HomeWin)
Supposons que je veuille exécuter un modèle probit dans lequel j'inclus des variables fictives pour l'équipe à domicile et l'équipe à l'extérieur.
model <- glm(HomeWin ~ as.factor(Home) + as.factor(Away), family = binomial(link="probit"), data = results)
Le résultat du modèle fournit des estimations de coefficient pour trois des équipes à domicile (par rapport à une équipe à domicile exclue) et trois des équipes à l'extérieur (par rapport à une équipe à l'extérieur exclue). Supposons que je veuille définir le modèle de telle sorte que l'estimation du coefficient du domicile pour NY soit égale à l'estimation du coefficient d'absence pour NY (et la même chose pour les autres villes). Comment pourrais-je faire ça? Mes données complètes contiennent 30 de ces groupes et avec beaucoup plus de variables.
Réponses
Si je comprends bien la question, ce que vous recherchez en fait, c'est d'avoir homeet awayd'avoir des effets opposés. Par exemple. beta_{home=NY} = - beta_{away=NY}. Ce n'est cependant pas tout à fait clair. Mais un moyen simple d'y parvenir serait de concevoir manuellement vos variables fictives, de telle sorte que vous ayez une variable fictive pour NY_home_or_awayavec home=1et away=-1. Dans ce cas, il beta_NY_home_or_awayserait basé à la fois à domicile et à l'extérieur, mais aurait un signe négatif.
library(dplyr)
competitors <- unique(unlist(results[, c('Home', 'Away')]))
new_cols <- lapply(competitors, function(x){
home <- results[['Home']] == x
away <- results[['Away']] == x
case_when(home ~ 1,
away ~ -1,
TRUE ~ 0)
})
names(new_cols) <- competitors
results_wide <- bind_cols(results, new_cols)
fit <- glm(HomeWin ~ NY + LA + CHI + BOS, data = results_wide, family = binomial('probit'))
summary(fit)
Call:
glm(formula = HomeWin ~ NY + LA + CHI + BOS, family = binomial("probit"),
data = results_wide)
Deviance Residuals:
Min 1Q Median 3Q Max
-1.64597 -0.73997 0.01633 1.19731 1.19731
Coefficients: (1 not defined because of singularities)
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.927e-02 3.823e-01 -0.077 0.939
NY 6.786e-01 6.676e-01 1.017 0.309
LA 6.786e-01 6.676e-01 1.017 0.309
CHI -2.898e-16 6.527e-01 0.000 1.000
BOS NA NA NA NA
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 16.636 on 11 degrees of freedom
Residual deviance: 14.537 on 8 degrees of freedom
AIC: 22.537
Number of Fisher Scoring iterations: 5
Notez que maintenant, le signe dépend du signe de savoir si l'équipe est Awayet Homecomme Away=-1. De plus, tout test statistique devrait probablement être effectué avec une certaine prudence après avoir effectué une telle transformation, car leur interprétation et leur validité dépendront d'autres variables. Notez également qu'une équipe recevra des NAestimations, car les variables factices sont linéairement dépendantes.
Vous pouvez créer des variables factices pour chaque nom d'équipe répertorié comme Domicile ou Absent et utiliser ces variables dans la régression.
(L'exemple ci-dessous peut fonctionner numériquement de manière étrange étant donné les exemples de données que vous avez fournis, mais il devrait fonctionner avec les données réelles.)
library(dplyr)
library(fastDummies)
teams <- results$Home %>% unique()
# function to add a dummy for a given team is either Home or Away
add_HoA <- function(df, team) {
HoA_str <- paste0('HoA_',team)
HoA <- ensym(HoA_str)
df <- df %>% mutate(!!HoA := (Home ==team | Away==team) %>% as.integer())
return (df)
}
for (team in teams) {
results <- add_HoA(results, team)
}
# using HoA_ variables for all teams
model2 <- glm(HomeWin ~ ., family = binomial(link="probit"),
data = results %>% dplyr::select(HomeWin, starts_with('HoA_')))
summary(model2)
results <- fastDummies::dummy_cols(results, select_columns = c('Home','Away'))
# using HoA_ variables for NY
model3 <- glm(HomeWin ~ ., family = binomial(link="probit"),
data = results %>%
dplyr::select(HomeWin, HoA_NY, starts_with('Home_'), starts_with('Away_')) %>%
dplyr::select(-Home_NY, -Away_NY))
summary(model3)
# using HoA_ variables for BOS
model4 <- glm(HomeWin ~ ., family = binomial(link="probit"),
data = results %>%
dplyr::select(HomeWin, HoA_BOS, starts_with('Home_'), starts_with('Away_')) %>%
dplyr::select(-Home_BOS, -Away_BOS))
summary(model4)