regfac.expand.1par {RegressionFactory}

Expander Function for Single-Parameter Base Distributions

Description

This function produces the full, high-dimensional gradient and Hessian from the base-distribution derivatives for linear transformations of the arguments of a single-parameter base distribution.

Usage

regfac.expand.1par(beta, X, y, fbase1, fgh=2, ...)

Arguments

beta

Vector of coefficients in the regression model.

X

Matrix of covariates in the regression model. Note that ncol(X) must be equal to length(beta).

y

Vector of response variable in the regression model. Note that length(y) must be equal to nrow(X).

fbase1

Base distribution function fbase1(u, y, ...) for the regression model. It must return a list with elements f,g,h corresponding to the function and its first and second derivatives relative to its first argument, u.

fgh

Integer with possible values 0,1,2. If fgh=0, the function only calculates and returns the log-likelihood function. If fgh=1, it returns the log-likelihood and its gradient vector. If fgh=2, it returns the log-likelihood, the gradient vector and the Hessian matrix.

...

Other parameters to be passed to fbase1.

Value

A list with elements f,g,h corresponding to the function, gradient vector, and Hessian matrix of the function fbase(X%*%beta,y), i.e. the base function fbase(u,y) projected onto the high-dimensional space of beta through the linear transformation of its first argument (u <- X%*%beta).

Author(s)

Alireza S. Mahani, Mansour T.A. Sharabiani

References

Mahani, Alireza S. and Sharabiani, Mansour T.A. (2013) Metropolis-Hastings Sampling Using Multivariate Gaussian Tangents https://arxiv.org/pdf/1308.0657v1.pdf

See Also

regfac.expand.2par

Examples

## Not run: 
library(sns)
# simulating logistic regression data
N <- 1000 # number of observations
K <- 10 # number of variables
X <- matrix(runif(N*K, min=-0.5, max=+0.5), ncol=K)
beta <- runif(K, min=-0.5, max=+0.5)
Xbeta <- X%*%beta
y <- 1*(runif(N)<1/(1+exp(-Xbeta)))
beta.est <- rep(0,K)
# run sns in non-stochastic mode, i.e. Newton-Raphson optimization
for (i in 1:10) {
  beta.est <- sns(beta.est, regfac.expand.1par, rnd=F, X=X, y=y
    , fbase1=fbase1.binomial.logit)
}
# use glm to estimate beta and compare
beta.est.glm <- glm(y~X-1, family="binomial")$coefficients
cbind(beta.est, beta.est.glm)

## End(Not run)