R: Ошибка с пользовательским кодом для информации о рыбаках (логистическая регрессия) - PullRequest
0 голосов
/ 25 января 2019

ниже - некоторый пользовательский код, который я создал для выполнения Логистической регрессии. Я использую тот факт, что функция оценки и информационная матрица Фишера для логистической регрессии могут быть выражены как функция строк расширенной матрицы данных.

I. h_b_x

Эта функция вычисляет: мю

h_b_x<- function(beta,x)
{
  if(length(beta) != length(x))
  {return("Incompatible length of vectors")}
  return((exp(sum(beta*x)))/(1 + exp(sum(beta*x))))
}

II. Функция счета

Эта функция вычисляет: Функция счета

score <- function(beta,X,y)
{
  # n denotes the number of observations in the data matrix X
  p <- dim(X)[2]

  # providing initial value for the score function
  sum_1 <- rep(0,p)

  for (i in 1:n) 
  {
    sum_1 <- sum_1 +  (y[i] - h_b_x(beta,X[i,])) * X[i,]
  }

  # sum_1 is a vector of length p 
  return(sum_1)
}

III. Информационная матрица рыболовов

Эта функция вычисляет: Информационная матрица Фишера

Fishers_matrix <- function(beta,X,y)
{
  # n denotes the number of observations in the data matrix X
  n <- dim(X)[1]

  # p denotes the number of parameters chosen for the logistic regression model
  p <- dim(X)[2]

  # providing initial value for the Fishers Information Matrix
  sum_1 <- matrix(rep(0,p*p), nrow = p, ncol = p)

  for (i in 1:n) 
  {
    sum_1 <- sum_1 +  (((X[i,] %*% t(X[i,])) * h_b_x(beta,X[i,]) * ( 1 - h_b_x(beta,X[i,]))))
  }

  # sum_1 is a matrix of length p * p 
  return(sum_1)
}

IV. EDA с набором данных Титаник

Этот раздел кода позволит вам скопировать EDA, который я провел:

library(titanic)
data("titanic_train")
data("titanic_test")


library(dplyr)
titanic_test$Survived <- 2
complete_data <- rbind(titanic_train, titanic_test)
complete_data$Embarked[complete_data$Embarked==""] <- "S"
complete_data$Age[is.na(complete_data$Age)] <- median(complete_data$Age,na.rm=S)
complete_data <- as.data.frame(complete_data)
titanic_data <- select(complete_data,-c(Cabin, PassengerId, Ticket, Name))

titanic_data <- titanic_data[!titanic_data$Survived == "2", ]

y_vals <- titanic_data$Survived



x <- as.data.frame(titanic_data[,-1])



X <- model.matrix(y_vals~.,data = x)

# Specifying intial value for beta
beta <- as.numeric(rep(0.01, dim(X)[2]))

V. Проблема, с которой я сталкиваюсь

Я разбил данные на 3 части -

X_1 <- X[1:297,]
X_2 <- X[298:594,]
X_3<- X[595:891,]

y_1 <- y_vals[1:297]
y_2 <- y_vals[298:594]
y_3 <- y_vals[595:891]

F_1 <- Fishers_matrix(beta,X_1,y_1)
F_2 <- Fishers_matrix(beta,X_2,y_2)
F_3 <- Fishers_matrix(beta,X_3,y_3)
F_4 <- score(beta,X_1,y_1)
F_5 <- score(beta,X_2,y_2)
F_6 <- score(beta,X_3,y_3)

Интересно, что сумма информации Фишера, хранящейся у каждой стороны, равна информации Фишера, основанной на совокупной информации. Однако то же самое не относится к функции оценки. Я нахожу невероятно трудным понять, почему это так, учитывая, что функцию оценки (как описано на 2-м изображении выше) можно вычислить как сумму строк данных (после корректировки на константы).

# Sum of each party's Fishers Info and Score function
Fisher_sum <- F_1 + F_2 + F_3
Score_sum <- F_4 + F_5 + F_6

# Fisher info and Score function of aggregate information
Fish_full <- Fishers_matrix(beta,X,y_vals)
score_full <- score(beta,X,y_vals)

sum (Fisher_sum - Fish_full)

sum(Score_sum - score_full )

Буду признателен за некоторую проницательность - я считаю, что здесь есть ошибка с моим кодом (но я рад, что исправлюсь, если моя статистика неверна.)

...