Q_function: Q(theta, theta*) objective function in the EM algorithm
In hlennon/copulaIVTS: Dependence modelling of a digitised Gaussian ARMA(p,q) copula model for Integer-Valued Time Series

Description Usage Arguments Details Value Note Author(s) References Examples

Evaluates the Q function in EM algorithm for the Gaussian copula model with NB marginals involving the terms E[YiYj | X=x]

Computes the Conditional likelihood given the observed data Computes e2/t by simulating yy (EXACTLY GHK) and averaging Arguements to be used: para (true parameter) & para.star (current parameter)

1	Q_function(para, para.star, m, low, up)

`para`	true parameters, theta, of the function
`para.star`	current parameters, theta star, of the function. Note theta* are the current values of the parameters in the iterative algorithm and therefore are fixed values once specified
`low`	lower bound of values, dependent only on data and marginal distribution the Negative Binomial
`up`	upper bound of values, dependent only on data and marginal distribution the Negative Binomial
`m`	The number of Monte Carlo samples, m, to approximate the conditional expectation

Objective Q funciton of the MCEM Algorithm

p - order of the ar parameters in ARMA(p,q), FIXED TO p=1 q - order of the ma parameters in ARMA(p,q), FIXED TO q=1

A numeric value of the conditional expectation of the complete likelihood function given the observed data

Requires: LD, ensure_stationarity_causality, GHK

Hannah Lennon

Dempster et al, (1977)

library(polynom)

para      <- c(0.9, 0.5)
para.star <- c(0.8, 0.5)

low  <- c(rep(-1, 10))
up   <- c(rep( 1, 10))

# Q_function(para, para.star, m=1, low, up)




## The function is currently defined as

Q_function <-function(para, para.star, m, low, up){
            p          <- 1
            q          <- 1
            N          <- length(low)
            para.star[1:p] <- ensure_causality_invertibility(para.star[1:p])         # Ensures parameters and current pa
            a.star     <- para.star[1:p]                      # ar current parameters
            b.star     <- para.star[(p+1):(p+q)]              # ma current parameters
            r.star     <- as.numeric(ARMAacf(ar=a.star, ma=b.star, lag.max=N, pacf=FALSE))   #r.star <- arma_cov(a.star, b.star, N); r.star <- r.star/r.star[1]
            P          <- toeplitz(r.star[1:N])                # Covariance matrix for Mult Truncated Normal
            para[1:p]  <- ensure_causality_invertibility(para[1:p])
            a          <- para[1:p]
            b          <- para[(p+1):(p+q)]
            r          <- as.numeric(ARMAacf(ar=a, ma=b, lag.max=N, pacf=FALSE))                     # Theoretical autocovariances for parameters
            output     <- LD(r,N-1)
            PHI        <- t(output[, 1:N])
            tau        <- output[,(N+1)]
            sim.data   <- GHK(m, P, low, up)
            xxx       <<- sim.data
            YY         <- sim.data *t(sim.data)/m
            tau.mat    <- diag(1/tau)                           # Matix with tau^-1 in diagonal elements
            A          <- PHI *tau.mat*t(PHI)                # Matrix form of equations
            like       <- sum(A*YY) + sum(log(tau))
            return(like)
}