Skip to content

Instantly share code, notes, and snippets.

@josephsdavid
Last active June 20, 2019 03:42
Show Gist options
  • Save josephsdavid/0f29eca6e59066ae2e5a6c20bc42228e to your computer and use it in GitHub Desktop.
Save josephsdavid/0f29eca6e59066ae2e5a6c20bc42228e to your computer and use it in GitHub Desktop.
Calculating AR(p) psi weights with R

Calculating psi weights of an AR(p) time series in R

Test Case

phis <- c(1.7,-.72)

Control

library(tswge)
psi.weights.wge(phi = phis, lag.max = 4)
# [1] 1.7000 2.1700 2.4650 2.6281

Goal function:

psi weight recursive definition

phis <- c(1.7,-.72)
l = 1
psi = as.numeric(1)
psi
# [1] 1
n <- phis[1:l] * psi[1:l]
# [1] 1.7
psi <- rlang::prepend(psi,n)
# [1] 1.7 1.0
l = 2
n  <- phis[1:l] * psi[1:l]
n <- sum(n)
# [1] 2.17
psi <- rlang::prepend(psi,n)
# [1] 2.17 1.70 1.00
l = 3
n  <- c(phis,0) * psi[1:l]
n <- sum(n)
# [1] 2.465
psi <- rlang::prepend(psi,n)
# [1] 2.465 2.170 1.700 1.000

First attempt

Failing to grow vector

psii <- function(phi,l){
	if (l>length(phi)){
		phi <- append(phi,rep(0,l-length(phis)))
	}

	if (l==1)	return(phi[1])
	else 		return(sum(phi[1:l]*(rlang::prepend(1,psii(phi,l-1)))))
}

Test:

psip <- function(phi,l){
	as.numeric(lapply(1:l, psii, phi = phi))
}
psip(phis,4)
# [1] 1.7000 2.1700 2.9690 4.3273

Why I think it broke: We are not making a vector of psi's, we are only successfully getting the first two values

second attempt: failed vectorization

For some reason, this code fails to output a vector. I am not sure why because i think we are growing a vector:

phitest <- function(phi,l){
	if (l>length(phi)){
		 return( append(phi,rep(0,l-length(phis))))
	}
	else{
		return(phi)
	}
}
psic <- function(phi,l, psi = as.numeric(1)){
	phi <- phitest(phi,l)
	if (l==1)	return(phi[1:l] * psi[1:l])
	else 		return(psic(phi, l - 1, rlang::prepend(psi, phi[1:l] * psi[1:l]))
	)
}

psic(phis,1)
# [1] 1.7
psic(phis,2)
# [1] 2.89

How to proceed????

Working for loop

We dont want to be doing this for loop because for loop bad, we want to vectorize or recursive

multfun <- function(phi,l,psi){
	sum(phi[1:l] * psi[1:l])
}

phitest <- function(phi,l){
	if (l>length(phi)){
		 return( append(phi,rep(0,l-length(phis))))
	}
	else{
		return(phi)
	}
}
psiloop <- function(phi,l){
	phi <- phitest(phi,l)
	psi <- as.numeric(1)
	for (i in 2:l){
		psi[i] <-  multfun(phi,i-1,rev(psi))
	}
psi
}
psiloop(phis,5)
# [1] 1.0000 1.7000 2.1700 2.4650 2.6281

Successful AR(p,q) loop:

recursive definition

multfun <- function(phi,l,psi){
	sum(phi[1:l] * psi[1:l])
}
phitest <- function(phi,l){
	if (l>length(phi)){
		 return( append(phi,rep(0,l-length(phis))))
	}
	else{
		return(phi)
	}
}
psiloop <- function(phi,l, theta = as.numeric(0)){
	phi <- phitest(phi,l)
	psi <- as.numeric(1) 
	theta <- c(0,theta)
	theta <- phitest(theta,l)
	for (i in 2:l){
		psi[i] <-  multfun(phi,i-1,rev(psi)) - theta[i]
	}
psi
}
@carson-drake
Copy link

def

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment