R/GEV_kuma_fit_ar.R

In leonardobfn/ThesiR:

GEV_kuma_fit_ar=function(teta,BETA,alpha,LAMBDA,cov_kl=cov_kl,cov_delta=cov_delta,
                         DADOS=DADOS,link_kl=link_kl,link_delta=link_delta,
                         Lm=Lm,w=w,quantil=quantil,Esp=Esp){

  phi_ar = teta
  gr <- levels((DADOS$grupos))
  vm=logr=vector(l=length(Lm))


  for(j in 1:length(Lm)){
    #j=1
    yy = DADOS %>% dplyr::filter(grupos==gr[j])
    y_grupo = matrix(yy[["y"]])
    TT = nrow(y_grupo)/Lm[j]

    cov_delta_grupo = as.matrix((yy[c(paste0(colnames(cov_delta),"_delta"))]))
    eta_delta = cov_delta_grupo%*%LAMBDA
    delta = link_delta$linkinv(eta_delta)
    cov_kl_grupo = as.matrix((yy[c(colnames(cov_kl))]))
    id_eta_first <- seq(1,nrow(y_grupo),TT)
    id_eta_last <- seq(TT,nrow(y_grupo),TT)
    eta_1 <- cov_kl_grupo[,-1]%*%BETA[-1]
    eta_2 <- phi_ar*(link_kl$linkfun(y_grupo) - eta_1)
    eta_kl <- BETA[1] + eta_1[-c(id_eta_first)] + eta_2[-c(id_eta_last)]
    kl = link_kl$linkinv(eta_kl)


    omega.aux = rep(w[[j]],each=TT)
    aux = logr_vmteta(quantil=quantil,omega = omega.aux[-c(id_eta_first)],dados=y_grupo[-c(id_eta_last)],kl=kl,delta=delta[-c(id_eta_first)])
    vm[j] = aux[1]
    logr[j] = aux[2]
  }
  #B=Esp[,3]
  #CB = (((sin(pi*alpha*B)/sin(pi*B))^(1/(1-alpha)))*((sin(pi*(1-alpha)*B))/sin(pi*alpha*B)))
  #alpha = 1/(1+Esp[,2])
  #print(CB)
  Q1 = sum(TT*(sapply(w,log) %>% sapply(sum))-vm*Esp[,1]+logr)
  #print(Q1)
  #Q2 = sum(log(alpha/(1-alpha))-(1/(1-alpha))*Esp[,2]+log(CB) - CB*(Esp[,1]^(-alpha/(1-alpha))))
  Q2 = sum(log(alpha/(1-alpha))-(1/(1-alpha))*Esp[,2]+Esp[,3] - Esp[,4]+TT*Lm*Esp[,2])

  # deltak = ifelse(it <= c*W,1,1/(it-c*W))
  # s = (Q1+Q2)
  # L = (Lant + deltak*(s - Lant))
  #print(c(s))
  L = Q1+Q2

  return(-L)

}