Commit 8bb9ac62 authored by Oscar Jurado's avatar Oscar Jurado
Browse files

modified gev.d.params to include the theta_zero option that comes from the gev.d.fit function

parent a1a6ad49
......@@ -247,7 +247,6 @@ gev.d.fit<-
}
class( z) <- "gev.d.fit"
invisible(z)
browser()
}
......@@ -447,7 +446,8 @@ gev.d.params <- function(fit,ydat){
npmu <- length(fit$model[[1]]) + 1
npsc <- length(fit$model[[2]]) + 1
npsh <- length(fit$model[[3]]) + 1
if(class(fit)=="gev.d.fit"){npth <- length(fit$model[[4]]) + 1}
if(class(fit)=="gev.d.fit" & !fit$theta_zero){npth <- length(fit$model[[4]]) + 1} #Including theta parameter (default)
if(class(fit)=="gev.d.fit" & fit$theta_zero){npth <- 0} #With no theta parameter, asked by user
if(class(fit)=="gev.d.fit"){npet <- length(fit$model[[5]]) + 1}
# inverse link functions
......@@ -455,7 +455,7 @@ gev.d.params <- function(fit,ydat){
mulink <- fit$link$mulink$linkinv
siglink <- fit$link$siglink$linkinv
shlink <- fit$link$shlink$linkinv
thetalink <- fit$link$thetalink$linkinv
if(!fit$theta_zero) thetalink <- fit$link$thetalink$linkinv
etalink <- fit$link$etalink$linkinv
}else{
mulink <- eval(parse(text=fit$link))[[1]]
......@@ -467,17 +467,22 @@ gev.d.params <- function(fit,ydat){
mumat <- cbind(rep(1, dim(ydat)[1]), matrix(ydat[, fit$model[[1]]],dim(ydat)[1],npmu-1))
sigmat <- cbind(rep(1, dim(ydat)[1]), matrix(ydat[, fit$model[[2]]],dim(ydat)[1],npsc-1))
shmat <- cbind(rep(1, dim(ydat)[1]), matrix(ydat[, fit$model[[3]]],dim(ydat)[1],npsh-1))
if(class(fit)=="gev.d.fit"){thmat <- cbind(rep(1, dim(ydat)[1]), matrix(ydat[, fit$model[[4]]],dim(ydat)[1],npth-1))}
if(class(fit)=="gev.d.fit" & !fit$theta_zero){thmat <- cbind(rep(1, dim(ydat)[1]), matrix(ydat[, fit$model[[4]]],dim(ydat)[1],npth-1))}
if(class(fit)=="gev.d.fit"){etmat <- cbind(rep(1, dim(ydat)[1]), matrix(ydat[, fit$model[[5]]],dim(ydat)[1],npet-1))}
# calculate parameters
mut <- mulink(mumat %*% (fit$mle[1:npmu]))
sc0 <- siglink(sigmat %*% (fit$mle[seq(npmu + 1, length = npsc)]))
xi <- shlink(shmat %*% (fit$mle[seq(npmu + npsc + 1, length = npsh)]))
if(class(fit)=="gev.d.fit"){theta <- thetalink(thmat %*% (fit$mle[seq(npmu + npsc + npsh + 1, length = npth)]))}
if(class(fit)=="gev.d.fit" & !fit$theta_zero){theta <- thetalink(thmat %*% (fit$mle[seq(npmu + npsc + npsh + 1, length = npth)]))}
if(class(fit)=="gev.d.fit"){eta <- etalink(etmat %*% (fit$mle[seq(npmu + npsc + npsh + npth + 1, length = npet)]))}
if(class(fit)=="gev.d.fit"){return(data.frame(mut=mut,sig0=sc0,xi=xi,theta=theta,eta=eta))
if(class(fit)=="gev.d.fit"){
if(!fit$theta_zero){ #When theta parameter is used (default)
return(data.frame(mut=mut,sig0=sc0,xi=xi,theta=theta,eta=eta))
}else{ #When theta parameter was not used
return(data.frame(mut=mut,sig0=sc0,xi=xi,eta=eta))
}
}else{return(data.frame(mu=mut,sig=sc0,xi=xi))}
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment