opt.maxLik: Main function of the package. Automated general to specific...

Description Usage Arguments Value Examples

View source: R/optMaxlik.R

Description

Numerical minimization of AICc by optimization of a likelihood function and automated constraining of parameters.

Usage

1
2
opt.maxLik(LL, start, initialfix = numeric(), nobs = Inf, method = "BFGS",
  penalized = numeric(), pw = 2)

Arguments

LL

a log likelihood function.

start

named vector of starting values.

initialfix

vector indicating the parameters that should be treated as constants with values supplied in the start vector.

nobs

number of observations to be used in the corrction of the AIC(c), defaults to Inf.

method

numerical algorithm. See maxLik package. Defaults to BFGS

penalized

vector indicating which parameters are penalized in your likelihood function. Will be used to compute adjustments for the AIC(c).

pw

the weight of you penalty. Supported penalties are abs(parameter value)^pw.

Value

maxLik object of final model.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
library(maxLik)
library(compiler) 
# Fill a matrix with some random data.
mydata<-matrix(rnorm(1000), ncol=10, nrow=100)

# Create you own Log likelihood function.
crossectionalARMA_44 <- function(pars, ret="LL", Y){
	 data=Y
	 Log.L <-numeric()[1:T]
	 b0    <-pars[1]

	 B.y   <-pars[2]
	 B.y2  <-pars[3]
	 B.y3  <-pars[4]
	 B.y4  <-pars[5]

	 B.e   <-pars[6]
	 B.e2  <-pars[7]
	 B.e3  <-pars[8]
	 B.e4  <-pars[9]

	 nu    <-pars[10]
	 sigma <-pars[11]

	 df = as.numeric.matrix(data)
	 T=nrow(df)
	 N=ncol(df)

	 p=4

	 e = matrix(0,T,N)
	 e[1:p,] <- 0 ## Initialize with e1 = 0

	 Log.L[1:p]<-0

	 A =N*log((gamma((nu+1)/2))/(((pi*(nu-2))^0.5)*gamma(nu/2))) - 0.5*N*log(max(0,sigma)^2)
	 A2=((nu+1)/2)
	 A3=(max(0,sigma)^2*(nu-2))  

	 for (t in (p+1):T) {
	  y = df[t,]#as.numeric(matrix(c(df[t,])))
	  ymin =  df[t-1,]#as.numeric(matrix(c(df[t-1,])))
	  ymin2 =  df[t-2,]#as.numeric(matrix(c(df[t-2,])))
	  ymin3 =  df[t-3,]#as.numeric(matrix(c(df[t-3,])))
	  ymin4 =  df[t-4,]#as.numeric(matrix(c(df[t-4,])))

	  emin=e[t-1,]
	  emin2=e[t-2,]
	  emin3=e[t-3,]
	  emin4=e[t-4,]

   MA = B.e*emin + B.e2*emin2 + B.e3*emin3 + B.e4*emin4 
	  e[t,] = y - b0 - B.y*ymin - B.y2*ymin2 - B.y3*ymin3 - B.y4*ymin4 - MA

	  }


	 Log.L <- A - A2*(rowSums(log(1+ (e)^2 / A3)))
	 Log.L[1:p]<-0

	 if(ret=="e"){return(e)}else if (ret =="LLvec") {return(Log.L)} else(sum(Log.L))

} 

# Optionally compile your function 

crossectionalARMA_44<-cmpfun(crossectionalARMA_44)

# opt.Maxlik takes only functions that have only a parameter vector as input.
# fix the data argument.

LL <- function(pars){crossectionalARMA_44(pars, ret="LL", Y=mydata)}

# create a vector of names parameter starting values.

start =c(b0=0, b1=0, b2=0, b3=0,b4=0, ma1=0, ma2=0, ma3=0, ma4=0, nu=120, sigma=sd(mydata))

# t-estimation
results <- opt.maxLik (LL=LL, start=start,initialfix=c(), nobs=dim(mydata)[1]*dim(mydata)[2])
summary(results)

estimates=coef(results)

parsToOriginal(estimates, start)

# approximate normal estiamtion by fixing the degrees of freedom
results2 <- opt.maxLik (LL=LL, start=start,initialfix=c(10), nobs=dim(mydata)[1]*dim(mydata)[2])
summary(results2)

estimates2=coef(results2)

parsToOriginal(estimates2, start)
parsToOriginal2(estimates2, start)

BPJandree/optMaxlik documentation built on Nov. 24, 2018, 6:04 p.m.