Menggunakan metode momen umum (GMM) untuk menghitung parameter regresi logistik

13

Saya ingin menghitung koefisien untuk regresi yang sangat mirip dengan regresi logistik (Sebenarnya regresi logistik dengan koefisien lain: ketika bisa diberikan). Saya berpikir untuk menggunakan GMM untuk menghitung koefisien, tapi saya tidak yakin kondisi apa yang harus saya gunakan.

A1+e(b0+b1x1+b2x2+),
A

Dapatkan seseorang menolong saya dengan itu?

Terima kasih!

pengguna5497
sumber
Ketika Anda mengatakan " dapat diberikan", apakah maksud Anda ditentukan oleh pengguna atau diperkirakan oleh model? SEBUAH
Makro
bagaimanapun juga. Saya dapat meletakkannya sebagai input (egA = 0,25) atau menjadi salah satu koefisien yang dapat ditemukan
user5497
Apakah ini berbeda dari satu subjek ke subjek lainnya (yaitu data) atau apakah konstanta tetap di semua pengamatan?
Makro
diperbaiki pada semua pengamatan (seperti b0, b1, ...)
user5497
2
Mengapa tidak menggunakan kemungkinan maksimum sebagai ganti GMM?
Makro

Jawaban:

6

Dengan asumsi , model ini memiliki variabel respon Bernoulli Y i denganA1Yi

Pr(Yi=1)=A1+eXib,

di mana (dan mungkin A , tergantung pada apakah itu diperlakukan sebagai konstanta atau parameter) adalah koefisien yang dipasang dan X i adalah data untuk pengamatan i . Saya berasumsi istilah intersepsi ditangani dengan menambahkan variabel dengan nilai konstan 1 ke matriks data.bAXii

Kondisi saat ini adalah:

E[(YiA1+eXib)Xi]=0.

N

m=1Ni=1N[(YiA1+eXib)Xi]=0

mmb

A

dat <- as.matrix(cbind(data.frame(IsVersicolor = as.numeric(iris$Species == "versicolor"), Intercept=1), iris[,1:4]))
head(dat)
#      IsVersicolor Intercept Sepal.Length Sepal.Width Petal.Length Petal.Width
# [1,]            0         1          5.1         3.5          1.4         0.2
# [2,]            0         1          4.9         3.0          1.4         0.2
# [3,]            0         1          4.7         3.2          1.3         0.2
# [4,]            0         1          4.6         3.1          1.5         0.2
# [5,]            0         1          5.0         3.6          1.4         0.2
# [6,]            0         1          5.4         3.9          1.7         0.4

Berikut adalah koefisien yang dipasang menggunakan regresi logistik:

summary(glm(IsVersicolor~., data=as.data.frame(dat[,-2]), family="binomial"))
# Coefficients:
#              Estimate Std. Error z value Pr(>|z|)    
# (Intercept)    7.3785     2.4993   2.952 0.003155 ** 
# Sepal.Length  -0.2454     0.6496  -0.378 0.705634    
# Sepal.Width   -2.7966     0.7835  -3.569 0.000358 ***
# Petal.Length   1.3136     0.6838   1.921 0.054713 .  
# Petal.Width   -2.7783     1.1731  -2.368 0.017868 *  

(YiA1+eXib)Xii

moments <- function(b, X) {
  A <- 1
  as.vector(X[,1] - A / (1 + exp(-(X[,-1] %*% cbind(b))))) * X[,-1]
}

b

init.coef <- lm(IsVersicolor~., data=as.data.frame(dat[,-2]))$coefficients
library(gmm)
fitted <- gmm(moments, x = dat, t0 = init.coef, type = "iterative", crit = 1e-19,
              wmatrix = "optimal", method = "Nelder-Mead",
              control = list(reltol = 1e-19, maxit = 20000))
fitted
#  (Intercept)  Sepal.Length   Sepal.Width  Petal.Length   Petal.Width  
#      7.37849      -0.24536      -2.79657       1.31364      -2.77834  
# 
# Convergence code =  0 

Kode konvergensi 0 menunjukkan prosedur konvergen, dan parameternya identik dengan yang dikembalikan oleh regresi logistik.

momentEstim.baseGmm.iterativegmm:::.obj1mmoptimgmm

gmm.objective <- function(theta, x, momentFun) {
  avg.moment <- colMeans(momentFun(theta, x))
  sum(avg.moment^2)
}
optim(init.coef, gmm.objective, x=dat, momentFun=moments,
      control = list(reltol = 1e-19, maxit = 20000))$par
#  (Intercept) Sepal.Length  Sepal.Width Petal.Length  Petal.Width 
#    7.3784866   -0.2453567   -2.7965681    1.3136433   -2.7783439 
Josliber
sumber