Задача нахождения оценок максимального правдоподобия параметров обобщенного нормального распределения в R - PullRequest
Новая
       30

Задача нахождения оценок максимального правдоподобия параметров обобщенного нормального распределения в R

0 голосов
/ 14 апреля 2019

Сначала я оценил параметры возведенного в степень обобщенного нормального распределения, используя набор данных (данные 1: сгенерированный из нормального распределения). Затем я использовал реальный набор данных (data2) и попытался найти оценки максимального правдоподобия (MLE), используя пакеты AdequacyModel . Это дает сообщение об ошибке (упоминается прямо под кодом).

  1. Почему один и тот же код оценивает первый набор данных (данные1), но не оценивает второй набор данных (данные2)?

  2. Является ли мой способ включения базового (нормального) распределения в обобщенный класс распределений F (x) = [1- (1-G (x)) ^ бета] ^ гамма введено Кордейру и др. (2013) в среде R верно? Правильно ли я определяю cdf и pdf обобщенного нормального распределения в R? 

# For first dataset
data1 <- rnorm(100)
pdf_exps <- function(par,x){
beta = par[1]
gamma= par[2]
mean = par[3]
sd = par[4]
     ( beta*gamma* ((1-(1-(pnorm(x,mean,sd)))^beta)^(gamma-1)) * ((1-   
     (pnorm(x,mean,sd)))^(beta-1)) ) * (dnorm(x,mean,sd)) 
}

cdf_exps <- function(par,x){
beta = par[1]
gamma= par[2]
mean = par[3]
sd = par[4]
( 1-(1-(pnorm(x,mean,sd)))^beta)^gamma 
}

set.seed(1)
result = goodness.fit(pdf = pdf_exps, cdf = cdf_exps,
starts = c(1,1,1,1),data = data1  , method = "BFGS",
domain = c(-Inf,Inf), lim_inf = c(0,0,0,0),
lim_sup = c(2,2,2,2), S = 250, prop=0.1, N=50)
result$mle
[1]   0.5976990 0.1541794 1.0073006 0.4689187

 # For second data set
 data2 <-c( 20.56, 20.67, 21.86, 21.88, 18.96, 21.04, 21.69, 20.62, 22.64, 19.44, 25.75, 21.20,
  22.03, 25.44, 22.63, 21.86, 22.27, 21.27, 23.47, 23.19, 23.17, 24.54, 22.96, 19.76,
  23.36, 22.67, 24.24, 24.21, 20.46, 20.81, 20.17, 23.06, 24.40, 23.97, 22.62, 19.16,
  21.15, 21.40, 21.03, 21.77, 21.38, 21.47, 24.45, 22.63, 22.80, 23.58, 20.06, 23.01,
  24.64, 18.26, 24.47, 23.99, 26.24, 20.04, 25.72, 25.64, 19.87, 23.35, 22.42, 20.42,
  22.13, 25.17, 23.72, 21.28, 20.87, 19.00, 22.04, 20.12, 21.35, 28.57, 26.95, 28.13,
  26.85, 25.27, 31.93, 16.75, 19.54, 20.42, 22.76, 20.12, 22.35, 19.16, 20.77, 19.37,
  22.37, 17.54, 19.06, 20.30, 20.15, 25.36, 22.12, 21.25, 20.53, 17.06, 18.29, 18.37,
  18.93, 17.79, 17.05, 20.31, 22.46, 23.88, 23.68, 23.15, 22.32, 24.02, 23.29, 25.11,
 22.81, 26.25, 21.38, 22.52, 26.73, 23.57, 25.84, 24.06, 23.85, 25.09, 23.84, 25.31,
 19.69, 26.07, 25.50, 23.69, 26.79, 25.61, 25.06, 24.93, 22.96, 20.69, 23.97, 24.64,
 25.93, 23.69, 25.38, 22.68, 23.36, 22.44, 22.57, 19.81, 21.19, 20.39, 21.12, 21.89,
 29.97, 27.39, 23.11, 21.75, 20.89, 22.83, 22.02, 20.07, 20.15, 21.24, 19.63, 23.58,
 21.65, 25.17, 23.25, 32.52, 22.59, 30.18, 34.42, 21.86, 23.99, 24.81, 21.68, 21.04,
 23.12, 20.76, 23.13, 22.35, 22.28, 23.55, 19.85, 26.51, 24.78, 33.73, 30.18, 23.31,
 24.51, 25.37, 23.67, 24.28, 25.82, 21.93, 23.38, 23.07, 25.21, 23.25, 22.93, 26.86,
 21.26, 25.43, 24.54, 27.79, 23.58, 27.56, 23.76, 22.01, 22.34, 21.07)
pdf_exps <- function(par,x){
beta = par[1]
gamma= par[2]
mean = par[3]
sd = par[4]
           ( beta*gamma* ((1-(1-(pnorm(x,mean,sd)))^beta)^(gamma-1)) * 
   ((1-(pnorm(x,mean,sd)))^(beta-1)) ) * (dnorm(x,mean,sd)) 
 }

cdf_exps <- function(par,x){
beta = par[1]
gamma= par[2]
mean = par[3]
sd = par[4]
( 1-(1-(pnorm(x,mean,sd)))^beta)^gamma 
}

set.seed(1)
result = goodness.fit(pdf = pdf_exps, cdf = cdf_exps,
starts = c(1,1,1,1),data = data2  , method = "BFGS",
domain = c(-Inf,Inf), lim_inf = c(0,0,0,0),
lim_sup = c(2,2,2,2), S = 250, prop=0.1, N=50)
result$mle

 Error in optim(par = starts, fn = likelihood, x = data, method = "BFGS",  
:non-finite finite-difference value [1]
...