convert_para: Convert parameters in the abk parametrization to pms...

In sqyu/ZiDAG: Directed Graphical Models and Causal Discovery for Zero-Inflated Data

Description

Gradient of log_dhurdle1d_abk() (with scalar v and y) with respect to a, b, k.

Usage

convert_para(paras, pms_to_abk)

Arguments

paras	A vector of 3 numbers. If pms_to_abk == TRUE, paras should contain the log odds and the conditional mean and variance for the Gaussian part. If pms_to_abk == FALSE, it should contain the a, b and k parameters.
pms_to_abk	A logical.

Details

The log density of a sample y from the 1-d Hurdle model with abk parametrization with respect to the sum of the Lebesgue measure and a point mass at 0 is a*v+b*y-y^2*k/2-log(1+sqrt(2pi/k)*exp(a+b^2/(2k))) The pms parametrization composes of the log odds (p) and the conditional mean (m) and variance (s) of the Gaussian part, i.e. the density is 1-p if v == 0, or p*dnorm(x, m, sqrt(s)) if v == 1. The two can be related as s = 1 / k, m = b / k, p = a + b^2/2/k - log(k/2/pi)/2. Also see Page 7 of McDavid (2019).

Value

A vector of 3 numbers.

Examples

set.seed(1)
convert_para(c(0, 0, 1), TRUE) # (-log(2*pi)/2, 0, 1)
convert_para(c(-log(2*pi)/2, 0, 1), FALSE) # (0, 0, 1)
paras <- c(rnorm(1), rnorm(1), abs(rnorm(1)))
convert_para(convert_para(paras, TRUE), FALSE) - paras
convert_para(convert_para(paras, FALSE), TRUE) - paras

n <- 100
v <- stats::rbinom(100, 1, 0.8)
y <- stats::rnorm(100) * v
a <- paras[1]; b <- paras[2]; k <- paras[3]
tmp <- convert_para(paras, FALSE)
log_odds <- tmp[1]; mu <- tmp[2]; sigmasq <- tmp[3]
prob <- exp(log_odds) / (exp(log_odds)+1)
likelihoods <- numeric(100)
likelihoods[v == 0] <- 1 - prob
likelihoods[v == 1] <- prob * dnorm(y[v == 1], mu, sqrt(sigmasq))
mean(log(likelihoods)) - log_dhurdle1d_abk(v, y, a, b, k)

sqyu/ZiDAG documentation built on Jan. 19, 2021, 4:11 p.m.