Package 'infometrics'

Default Support Spaces from OLS

Description

Constructs default signal (Z) and error (V) support spaces for GME estimation using OLS estimates as reference points. The signal support for each coefficient is centered on the OLS estimate with half-range equal to signal_scale times the absolute OLS estimate. The error support is symmetric around zero with half-range equal to error_scale times the OLS residual standard deviation.

Usage

default_supports(
  y,
  X,
  M_signal = 5L,
  M_error = 5L,
  signal_scale = 2,
  error_scale = 3
)

Arguments

y	Numeric vector of length T. The dependent variable.
X	Numeric matrix of dimension T x K. The regressor matrix.
M_signal	Integer. Number of support points per coefficient (default 5).
M_error	Integer. Number of support points for error terms (default 5).
signal_scale	Positive numeric. Half-range multiplier for the signal support (default 2: support spans +/- 2 * |OLS estimate|).
error_scale	Positive numeric. Half-range multiplier for the error support (default 3: support spans +/- 3 * OLS residual SD, the "three-sigma rule" of Pukelsheim, 1994).

Value

A list with two elements:

Z

K x M_signal matrix. Row k is the support for beta_k.

V

Numeric vector of length M_error. The error support, symmetric around zero.

References

Golan, A. (2008). Information and Entropy Econometrics. Foundations and Trends in Econometrics, 2(1-2), 1-145. Section 6.1.

Pukelsheim, F. (1994). The three sigma rule. The American Statistician, 48, 88-91.

Examples

set.seed(42)
n <- 50; K <- 3
X <- cbind(1, matrix(rnorm(n * (K-1)), n, K-1))
y <- X %*% c(1, 2, -1) + rnorm(n)
sp <- default_supports(y, X)
sp$Z   # 3 x 5 signal support matrix
sp$V   # 5-point error support

---

Information-Theoretic (Fano) Error Bounds

Description

Generic for Fano's-inequality error bounds on a fitted probability matrix.

Usage

fano_bounds(object, ...)

Arguments

object	A fitted model object.
...	Passed to methods.

Value

A method-specific object of error bounds.

See Also

fano_bounds.multinomial_gce

---

Fano Error Bounds for a Pure Inverse-Problem (Cross-Entropy) Fit

Description

The recovered p is a single distribution over the $K$ states, so Golan's (2008, sec 7.5) Fano bound applies to it directly and is the natural companion to the normalized entropy S: reading S as an irreducible prediction error, a modal classifier that guesses $\arg\max_k p_k$ errs with probability $pe = 1 - \max_k p_k$, and Fano's inequality lower-bounds this by

$pe \ge S(p) - \log 2/\log K, \qquad S(p) = H(p)/\log K.$

The S used here is the uniform-reference $H(p)/\log K$ the Fano bound requires – distinct from inverse_ce's prior-relative reported $S = H(p)/H(p_0)$. This is an information-theoretic error bound on prediction accuracy, not a sampling standard error (for the identification SEs of the multipliers see vcov.inverse_ce / summary).

Usage

## S3 method for class 'inverse_ce'
fano_bounds(object, ...)

Arguments

object	An inverse_ce object.
...	Unused.

Value

A one-row data frame with columns p_max, pe (modal error), H (entropy, nats), S (uniform-normalized entropy), and pe_lower (Fano weak lower bound). An "overall" attribute holds mean_pe, mean_pe_lower, and S_system (here identical to the single row).

References

Golan, A. (2008). Information and Entropy Econometrics. Foundations and Trends in Econometrics, 2(1-2), 1-145 (sec 7.5); Fano, R. (1961). Transmission of Information.

See Also

fano_bounds, inverse_ce

---

Fano Error Bounds for a Noisy Inverse-Problem (GME/GCE) Fit

Description

The recovered p_hat is a single distribution over the $K$ states, so Golan's (2008, sec 7.5) Fano bound applies to it directly (a parallel of fano_bounds.inverse_ce): a modal classifier that guesses $\arg\max_k p_k$ errs with probability $pe = 1 - \max_k p_k$, and $pe \ge S(p) - \log 2/\log K$ with $S(p) = H(p)/\log K$ (uniform-reference, distinct from the prior-relative S/S_w the object reports). An information-theoretic error bound, not a sampling SE (for those see vcov.inverse_noise / summary).

Usage

## S3 method for class 'inverse_noise'
fano_bounds(object, ...)

Arguments

object	An inverse_noise object.
...	Unused.

Value

A one-row data frame with columns p_max, pe, H, S, and pe_lower; an "overall" attribute mirrors the row.

References

Golan, A. (2008). Information and Entropy Econometrics. Foundations and Trends in Econometrics, 2(1-2), 1-145 (sec 7.5); Fano, R. (1961). Transmission of Information.

See Also

fano_bounds, inverse_noise

---

Fano Error Bounds for a Mixed GCE Fit

Description

For each observation (a row of p_hat, a distribution over the $J$ categories) the modal classifier predicts $\arg\max_j p_{ij}$ with error $pe_i = 1 - \max_j p_{ij}$, lower-bounded (Golan 2008, sec 7.5) by $S(p_i) - \log 2/\log J$. This is an information-theoretic error bound, not a sampling standard error (for those see margins with se = TRUE).

Usage

## S3 method for class 'mixed_gce'
fano_bounds(object, ...)

Arguments

object	A mixed_gce object.
...	Unused.

Value

A data frame with one row per observation (p_max, pe, H, S, pe_lower) and an "overall" attribute.

See Also

fano_bounds, margins

---

Fano Error Bounds for a Multinomial GCE Fit

Description

For each observation (a row of p_hat, a distribution over the $J$ alternatives) the modal classifier predicts $\arg\max_j p_{ij}$ with error probability $pe_i = 1 - \max_j p_{ij}$. Fano's inequality (Golan 2008, sec 7.5) lower-bounds this error by the normalized entropy:

$pe_i \ge S(p_i) - \log 2/\log J,\qquad S(p_i)=H(p_i)/\log J.$

This is an information-theoretic error bound on prediction accuracy, not a sampling standard error (for the latter see margins with se = TRUE).

Usage

## S3 method for class 'multinomial_gce'
fano_bounds(object, ...)

Arguments

object	A multinomial_gce object.
...	Unused.

Value

A data frame with one row per observation and columns p_max, pe (modal error), H (entropy, nats), S (normalized entropy), and pe_lower (Fano weak lower bound). An "overall" attribute holds mean_pe, mean_pe_lower, and S_system.

References

Golan, A. (2008). Information and Entropy Econometrics. Foundations and Trends in Econometrics, 2(1-2), 1-145 (sec 3.6, 7.5); Fano, R. (1961). Transmission of Information.

See Also

margins

---

Pure Inverse Problem via Cross-Entropy (Formula Interface)

Description

Estimates the K-dimensional distribution p closest (in Kullback-Leibler divergence) to a prior p0 subject to the T pure moment constraints $y = X p$, by solving the unconstrained dual (concentrated) model of Golan (2008, Section 4.2; see also Golan, Judge & Miller, 1996, Ch. 3). This is the lm-style formula interface to the maximum-entropy / cross-entropy problem, where a uniform prior p0 gives Maximum Entropy; a non-uniform p0 gives Cross-Entropy (hence the name).

Usage

inverse_ce(formula, data, p0 = NULL, subset, na.action, control = list(), ...)

## S3 method for class 'inverse_ce'
coef(object, ...)

## S3 method for class 'inverse_ce'
fitted(object, ...)

## S3 method for class 'inverse_ce'
residuals(object, ...)

## S3 method for class 'inverse_ce'
vcov(object, ...)

## S3 method for class 'inverse_ce'
print(x, digits = max(3L, getOption("digits") - 3L), ...)

## S3 method for class 'inverse_ce'
summary(object, ...)

## S3 method for class 'summary.inverse_ce'
print(x, digits = max(3L, getOption("digits") - 3L), ...)

Arguments

formula	A model formula. The response is the moment vector y; each RHS term is a state. Usually written without an intercept, e.g. y ~ x1 + x2 + x3 - 1.
data	A data frame (or environment) with one row per moment.
p0	Numeric vector of length K (number of model-matrix columns): strictly-positive prior probabilities for the Cross-Entropy formulation. Defaults to the uniform distribution, which makes CE reduce to ME. (Named p0 to match the GME/GCE signal-prior convention; for inverse_ce it is a K-vector prior over the states, not a K-by-M support-point prior.)
subset	Optional vector specifying a subset of moments (rows) to use, as in lm.
na.action	A function indicating how to handle NAs in the model frame, as in lm.
control	A named list of control parameters merged over the defaults and passed to optim (BFGS): maxit (default 500) and reltol (default 1e-10).
...	Currently unused.
object, x	An inverse_ce object.
digits	Number of significant digits to print.

Details

The dual concentrated objective (Golan 2008, Eq. 4.4-4.5), minimised here to match the package's dual sign convention, is

$\ell(\lambda) = -\sum_t \lambda_t y_t + \log \Omega(\lambda),$

with partition function $\Omega(\lambda) = \sum_k p_{0,k} \exp(\sum_t \lambda_t x_{tk})$. Minimising this convex objective is equivalent to maximising the entropy of p subject to the moment constraints. With uniform p0 this is Maximum Entropy; with a non-uniform prior it is Cross-Entropy. The estimated probabilities are recovered as

$\hat p_k = p_{0,k} \exp(\textstyle\sum_t \hat\lambda_t x_{tk}) / \Omega(\hat\lambda)$

(Golan 2008, Eq. 4.3), and the analytic dual Hessian is $\nabla^2 \ell = X \,\mathrm{diag}(p)\, X' - (Xp)(Xp)' = \mathrm{Cov}_p(\text{moments})$ (Eq. 4.7), which is positive semidefinite.

Value

An object of class c("inverse_ce", "infometrics"): a list with components

p_hat

Named K-vector of estimated probabilities.

lambda_hat

Named T-vector of Lagrange multipliers (dual variables).

fitted.values

Fitted moments $X p$ (length T).

residuals

Moment residuals $y - X p$ (length T).

hessian

Analytic T-by-T dual Hessian = information matrix $I(\lambda)=\mathrm{Cov}_p(\text{moments})$ (Eq. 4.7), positive semidefinite.

vcov_lambda

$\mathrm{Var}(\lambda)=I^{-1}(\lambda)$ (T-by-T), or NULL when $I$ is singular.

se_lambda

Information-matrix SEs of $\lambda$ (length T), or NA when $I$ is singular.

se_p

Delta-method SEs of the probabilities p (length K), or NA when $I$ is singular.

H_signal

Shannon entropy $H(\hat p)$.

S

Normalized entropy $H(\hat p)/H(p_0)$ in [0, 1] (pseudo-R2 = 1 - S).

objective

Dual objective at the optimum.

p0

The prior distribution used.

convergence

Raw optim convergence code.

method

Solver used ("dual").

terms, model, call

Standard model components.

Methods (by generic)

• coef(inverse_ce): Estimated probabilities p.

• fitted(inverse_ce): Fitted moments $X p$.

• residuals(inverse_ce): Moment-fitting residuals $y - X p$.

• vcov(inverse_ce): Information-matrix covariance of the multipliers $\mathrm{Var}(\lambda)=I^{-1}(\lambda)$ (T-by-T), or NULL when $I$ is singular.

• print(inverse_ce): Compact printed display.

• summary(inverse_ce): lm()-style summary (information-matrix SE tables + Fano).

Functions

• print(summary.inverse_ce): Print method for the summary object.

Information-matrix standard errors

The Hessian above is the Fisher information matrix of the multipliers, $I(\lambda) = \mathrm{Cov}_p(\text{moments})$, so its inverse gives their covariance $\mathrm{Var}(\lambda) = I^{-1}(\lambda)$ (Golan 2008, p. 59; Cover & Thomas, 2006, Ch. 17). inverse_ce() reports the resulting se_lambda, the delta-method probability SEs se_p ($\mathrm{Var}(p) = J\,I^{-1}J'$ with $J_{kt} = \partial p_k/\partial\lambda_t = p_k(x_{tk}-E_p[x_t])$), and the full vcov_lambda (via vcov). These are information-matrix / curvature quantities: they measure how well the moments identify $\lambda$ – following the multiplier's reading as the marginal information content of a moment (a moment $t$ with larger $\mathrm{Var}_p(x_t)$ carries more information, hence a smaller $SE(\lambda_t)$; $\lambda_t\approx 0$ means moment $t$ adds little) – and are not frequentist sampling standard errors, because a pure inverse problem is deterministic. $I(\lambda)$ is singular when a moment row is constant/collinear or $T \ge K$ (its rank is at most $K-1$); there the affected SEs are returned as NA rather than a pseudo-inverse's misleadingly finite value. For sampling-based inference use the stochastic-moment inverse_noise.

Formula orientation

Unlike a regression, the response is the vector of moments and each right-hand-side term is a state. The model matrix built from formula and data is therefore the T-by-K matrix $X$ with rows = moment constraints (T) and columns = states / support points (K); the data frame has one row per moment. The estimated coefficients (coef) are the K probabilities p, one per column. Because columns are states, an intercept would add a spurious extra state; drop it with - 1. A warning is issued if an intercept is present.

References

Golan, A. (2008). Information and Entropy Econometrics - A Review and Synthesis. Foundations and Trends in Econometrics, 2(1-2), 1-145. Section 4.2. Cover, T. M. & Thomas, J. A. (2006). Elements of Information Theory, 2nd ed., Ch. 17. Golan, A., Judge, G. and Miller, D. (1996). Maximum Entropy Econometrics: Robust Estimation with Limited Data. Wiley. Chapter 3.

See Also

inverse_noise for the stochastic-moment (GME/GCE) estimator that reports sampling standard errors; fano_bounds for the Fano error bound on the recovered p.

Examples

# Recover a 5-state distribution from 2 exact moment constraints.
# The data frame has one row per moment (T = 2) and one column per
# state (K = 5). Because p_true is the maximum-entropy distribution for
# its own moments, it is recovered exactly.
X <- rbind(c(1, 2, 3, 4, 5),
           c(1, 4, 9, 16, 25))           # 2 moments x 5 states
p_true <- exp(as.vector(crossprod(X, c(0.2, -0.05))))
p_true <- p_true / sum(p_true)
dat <- data.frame(y  = as.vector(X %*% p_true),
                  s1 = X[, 1], s2 = X[, 2], s3 = X[, 3],
                  s4 = X[, 4], s5 = X[, 5])
fit <- inverse_ce(y ~ s1 + s2 + s3 + s4 + s5 - 1, data = dat)
fit
coef(fit)          # ~ p_true
summary(fit)       # info-matrix SEs of lambda / p + a Fano line
vcov(fit)          # Var(lambda) = I^{-1}
fano_bounds(fit)   # Fano error bound for the recovered p

---

Noisy Inverse Problem via Generalized Cross-Entropy (Formula Interface)

Description

Estimates a K-dimensional signal distribution p together with per-moment noise distributions w for the stochastic-moment problem $y = X p + e$, $e_t = \sum_j v_j w_{tj}$, by solving the unconstrained dual (concentrated) model of Golan (2008, Section 6.1). This is the noisy-moment counterpart of inverse_ce.

Usage

inverse_noise(
  formula,
  data,
  v = NULL,
  nu = 0.5,
  p0 = NULL,
  w0 = NULL,
  se_method = c("sandwich", "delta", "bootstrap", "none"),
  B = 1000L,
  subset,
  na.action,
  control = list(),
  ...
)

## S3 method for class 'inverse_noise'
coef(object, ...)

## S3 method for class 'inverse_noise'
fitted(object, ...)

## S3 method for class 'inverse_noise'
residuals(object, ...)

## S3 method for class 'inverse_noise'
vcov(object, type = c("sandwich", "delta", "bootstrap"), B = 1000L, ...)

## S3 method for class 'inverse_noise'
print(x, digits = max(3L, getOption("digits") - 3L), ...)

## S3 method for class 'inverse_noise'
summary(object, se_method = NULL, B = 1000L, ...)

## S3 method for class 'summary.inverse_noise'
print(x, digits = max(3L, getOption("digits") - 3L), ...)

Arguments

formula	A model formula. The response is the moment vector; each RHS term is a state. Usually written without an intercept, e.g. y ~ x1 + x2 + x3 - 1.
data	A data frame (or environment) with one row per moment.
v	Noise support. One of: NULL (default 3-point symmetric support of half-width $\max(3, \sqrt{2\log T})\,s$ with $s = sd(y)$); a single whole number M >= 2 (build an M-point symmetric support); or an explicit numeric vector of support points (should be symmetric around 0 for zero-mean noise). The default widens with T rather than using a fixed three-sigma rule, because the largest of T noise draws grows like $\sigma\sqrt{2\log T}$; see the Support feasibility section.
nu	Entropy weight in $(0,1)$ on the noise term; default 0.5 (equal-weight GME/GCE).
p0	Prior signal probabilities, length K (= model-matrix columns); strictly positive, default uniform (GME). A non-uniform prior gives GCE.
w0	Prior noise probabilities, a T x M matrix (one row per moment); strictly positive rows summing to 1, default uniform. Its column count must match length(v).
se_method	For summary: SE method to report ("sandwich", "delta", or "bootstrap"); defaults to the stored method.
B	Integer number of bootstrap resamples when se_method = "bootstrap" (or vcov/summary with type/se_method = "bootstrap"); default 1000.
subset	Optional vector specifying a subset of moments (rows) to use, as in lm.
na.action	A function indicating how to handle NAs in the model frame, as in lm.
control	A named list merged over the defaults and passed to optim (BFGS): maxit (default 500) and reltol (default 1e-10).
...	Currently unused.
object, x	An inverse_noise object.
type	Covariance method for vcov: "sandwich" (default), "delta", or "bootstrap".
digits	Number of significant digits to print.

Details

Unlike inverse_ce, the moments are treated as noisy: the constraint is $y = X p + e$ rather than $y = X p$. The estimator therefore shrinks p toward its prior p0; how much depends on the width of the noise support v. As v collapses toward 0 the solution approaches the pure-CE solution of inverse_ce.

The dual objective is minimised here (convex) to match the package convention used by inverse_ce:

$\ell(\lambda) = -\sum_t \lambda_t y_t + (1-\nu)\log\Omega(\lambda) + \nu \sum_t \log\Psi_t(\lambda),$

with signal partition function $\Omega = \sum_k p_{0,k} \exp(\sum_t \lambda_t x_{tk}/(1-\nu))$ and per-moment noise partition functions $\Psi_t = \sum_j w_{0,tj}\exp(\lambda_t v_j/\nu)$. The weight $\nu \in (0,1)$ splits entropy between signal and noise; $\nu = 0.5$ gives the standard equal-weight GME/GCE. Minimising this convex objective is equivalent to maximising the joint signal-plus-noise entropy subject to the noisy moment constraints. The estimates are recovered as $\hat p_k \propto p_{0,k}\exp(\sum_t \hat\lambda_t x_{tk}/(1-\nu))$ and $\hat w_{tj} \propto w_{0,tj}\exp(\hat\lambda_t v_j/\nu)$, with $\hat e_t = \sum_j v_j \hat w_{tj}$. Uniform priors (p0 = NULL, w0 = NULL) give the GME version; user-supplied priors give GCE.

Value

An object of class c("inverse_noise", "infometrics"): a list with components

p_hat

Named K-vector of signal probabilities.

foc_residual

$\max|y - Xp - e|$ at the optimum; near 0 for a healthy fit (see Support feasibility).

lambda_hat

Named T-vector of Lagrange multipliers.

w_hat

T x M matrix of noise probabilities.

fitted.values

Signal-fitted moments $X p$ (length T).

residuals

Moment residuals $y - X p$, equal to the estimated noise $e_t = \sum_j v_j w_{tj}$ at the optimum.

hessian

Analytic T-by-T dual Hessian (positive definite).

vcov_lambda, se_lambda, se_p

Covariance of $\lambda$ and standard errors of $\lambda$/$p$ by se_method (NULL when se_method = "none").

se_method

The SE method stored.

H_signal, H_error, H_error0

Shannon entropies of the signal $H(\hat p)$, the noise $H(\hat w)$, and the noise prior $H(w_0)$.

S, S_p

Normalized signal entropy $H(\hat p)/H(p_0)$ in [0, 1] (pseudo-R2 = 1 - S).

S_w

Normalized noise entropy $H(\hat w)/H(w_0)$.

objective

Dual objective at the optimum.

nu, prior, noise_prior, support

The entropy weight, signal prior p0, noise prior w0, and the noise support v.

X, y, control

Stored for SE recomputation and the bootstrap.

convergence

Raw optim convergence code.

method, terms, model, call

Solver and standard model components.

Methods (by generic)

• coef(inverse_noise): Estimated signal probabilities p.

• fitted(inverse_noise): Fitted (signal) moments $X p$.

• residuals(inverse_noise): Residuals $y - X p$ (the estimated noise).

• vcov(inverse_noise): Covariance of the Lagrange multipliers by method. type = "sandwich" (default) is the robust $H^{-1}\mathrm{diag}(\hat e^2)H^{-1}$; "delta" is the naive $H^{-1}$ (overstates ~10x); "bootstrap" refits B residual resamples. Note: the default changed from the raw $H^{-1}$ to the sandwich, since $H^{-1}$ is not a valid sampling covariance here.

• print(inverse_noise): Compact printed display.

• summary(inverse_noise): lm()-style summary; prints a p coefficient table with standard errors and t-stats, the signal/noise normalized entropies, and a Fano line. se_method selects the SE method (default: the one stored on the object).

Functions

• print(summary.inverse_noise): Print method for the summary object.

Standard errors

Structurally this is a GME regression $y = Xp + e$ with p on the simplex, so the $T$ moment rows are the effective sample size and sampling SEs are meaningful. Because $\hat\lambda$ solves $Xp(\lambda)+e(\lambda)-y=0$, the implicit-function theorem gives $\partial\hat\lambda/\partial y = H^{-1}$, so $\mathrm{Var}(\hat\lambda) = H^{-1}\Sigma_e H^{-1}$ and $\mathrm{Var}(\hat p) = J\,\mathrm{Var}(\hat\lambda)\,J'$ with $J_{kt} = \partial p_k/\partial\lambda_t = (1-\nu)^{-1} p_k(x_{tk}-(Xp)_t)$. Three methods are provided via se_method / vcov(type=) / summary(se_method=):

"sandwich" (default)

robust HC0 $H^{-1}\,\mathrm{diag}(\hat e^2)\,H^{-1}$ (meat = squared fitted noise). Accurate for p (validated against a Monte-Carlo SD and the bootstrap); per-element lambda SEs are noisier (each rests on one $\hat e_t^2$).

"delta"

the naive $H^{-1}$. Kept for comparison but it overstates the sampling SE by roughly ten-fold (it is not a sampling covariance) – do not use it for inference.

"bootstrap"

residual bootstrap: keep the rows fixed, resample the fitted noise $\hat e$, refit. Gives stable, aligned SEs for both lambda and p; the gold-standard cross-check.

All three mildly understate (~15-20%) because the finite-support entropy penalty shrinks the fitted noise below the true noise dispersion (the bootstrap shows the same, so it is a property of the estimator and of an over-wide support, not a bug).

Support feasibility

The dual is bounded below only if y can be represented as $Xp + e$ with $p$ on the simplex and every $e_t$ inside the noise support v. If it cannot, the dual is unbounded: the multipliers diverge, the softmax saturates, and p_hat collapses to a vertex of the simplex (a 0/1 corner) even though optim reports convergence. inverse_noise() therefore checks the first-order condition $y - Xp - e = 0$ at the optimum, reports it as foc_residual, and warns when it is not satisfied. The remedy is a wider noise support v.

Formula orientation

Unlike a regression, the response is the vector of moments and each right-hand-side term is a state. The model matrix built from formula and data is therefore the T-by-K matrix $X$ with rows = moment constraints (T) and columns = states / support points (K); the data frame has one row per moment. The estimated coefficients (coef) are the K probabilities p, one per column. Because columns are states, an intercept would add a spurious extra state; drop it with - 1. A warning is issued if an intercept is present.

References

Golan, A. (2008). Information and Entropy Econometrics - A Review and Synthesis. Foundations and Trends in Econometrics, 2(1-2), 1-145. Section 6.1.

See Also

inverse_ce for the exact-moment (pure) version; linreg for the GME/GCE regression estimator; fano_bounds for the Fano error bound on p_hat.

Examples

set.seed(123)
X <- matrix(runif(30), nrow = 10, ncol = 3)   # T = 10 moments, K = 3 states
p_true <- c(0.5, 0.2, 0.3)
y <- as.vector(X %*% p_true) + rnorm(nrow(X), 0, 0.3)
dat <- data.frame(y = y, x1 = X[, 1], x2 = X[, 2], x3 = X[, 3])
fit <- inverse_noise(y ~ x1 + x2 + x3 - 1, data = dat)
fit
summary(fit)                 # p table (SE + t), signal/noise entropy, Fano
vcov(fit)                    # sandwich Cov(lambda)
fano_bounds(fit)             # Fano error bound for p_hat

---

Generalized Cross-Entropy Linear Regression (Formula Interface)

Description

Fits the linear model $y = X\beta + e$ by Generalized Cross-Entropy (Golan 2008, Section 6.1, Eq. 6.5; see also Golan, Judge & Miller, 1996, Ch. 6). Each coefficient is reparameterized on a bounded support, $\beta_k = \sum_m z_{km} p_{km}$, and each error on a noise support, $e_t = \sum_j v_j w_{tj}$; the signal and noise entropies (relative to priors p0, w0) are maximized subject to the data constraints, with weight $\nu \in (0,1)$ on the noise. Estimation uses the concentrated dual in the T Lagrange multipliers.

Usage

linreg(
  formula,
  data,
  Z = NULL,
  nu = 0.5,
  p0 = NULL,
  v = NULL,
  w0 = NULL,
  subset,
  na.action,
  control = list(),
  ...
)

## S3 method for class 'linreg'
coef(object, ...)

## S3 method for class 'linreg'
fitted(object, ...)

## S3 method for class 'linreg'
residuals(object, ...)

## S3 method for class 'linreg'
vcov(object, ...)

## S3 method for class 'linreg'
predict(object, newdata, ...)

## S3 method for class 'linreg'
print(x, digits = max(3L, getOption("digits") - 3L), ...)

## S3 method for class 'linreg'
summary(object, ...)

## S3 method for class 'summary.linreg'
print(
  x,
  digits = max(3L, getOption("digits") - 3L),
  signif.stars = getOption("show.signif.stars"),
  ...
)

Arguments

formula, data, subset, na.action	Standard model-frame arguments, as in lm. The intercept is kept by default.
Z	Coefficient support. Either a length-M numeric vector (a common support recycled across all coefficients) or a K-by-M matrix whose rows correspond, in order, to the columns of the model matrix (so row 1 is the intercept's support when an intercept is present). Default: a 5-point symmetric support spanning +/- a multiple of sd(y).
nu	Entropy weight in $(0,1)$ on the noise term, hence the weight on the signal is $1-\nu$; default 0.5 (equal weight on both the signal and the noise).
p0	Prior signal probabilities, a K-by-M matrix matching Z; strictly positive rows summing to 1. Default uniform.
v, w0	Noise support and prior, as in inverse_noise: v is NULL (3-point default of half-width $\max(3, \sqrt{2\log T})\,sd(y)$, which widens with the sample size – see Support feasibility), a whole-number support count, or an explicit vector; w0 is a strictly-positive T-by-J matrix with rows summing to 1 (default uniform).
control	A named list merged over the defaults and passed to optim (BFGS): maxit (default 500) and reltol (default 1e-12).
...	Currently unused.
object, x	A linreg object.
newdata	Optional data frame of new observations.
digits	Significant digits to print.
signif.stars	Logical; show significance stars (default from options).

Details

This fits the Generalized Cross-Entropy (GCE) linear regression estimator with an lm-style interface, adding a predict method, an lm-style coefficient table, and r.squared. The dual is minimised here (convex) to match the package convention used by inverse_ce and inverse_noise; this is equivalent to maximising the joint signal-plus-noise entropy. Uniform priors (p0 = NULL, w0 = NULL) give GME; user priors give GCE.

summary() reports a per-coefficient Entropy-Ratio (ER) test of $H_0\!:\beta_k = 0$ (Golan 2008, Sec. 6.4/6.6) in place of a Wald z test. To do so it refits the model $K+1$ times (one per coefficient plus one joint test of all coefficients = 0), each warm-started from the fitted multipliers; the resolved Z, p0, v, w0, X, y are stored on the object to make these refits faithful.

Value

An object of class c("linreg", "infometrics"), which is a list containing the following components:

coefficients

Numeric vector of length K: the estimated coefficients $\hat\beta = Z\hat p$, named by the model-matrix columns. Extracted by coef.

p, p_hat

K-by-M matrix of estimated signal probabilities (two names for the same object).

lambda, lambda_hat

Numeric vector of length T: the estimated Lagrange multipliers (two names for the same object).

w, w_hat

T-by-J matrix of estimated noise probabilities (two names for the same object).

e

Numeric vector of length T: the estimated noise, $e_t = \sum_j v_j w_{tj}$.

fitted.values

Numeric vector of length T: $X\hat\beta$. Extracted by fitted.

residuals

Numeric vector of length T: $y - X\hat\beta$, which equals e at a converged optimum. Extracted by residuals.

vcov

K-by-K asymptotic covariance matrix of $\hat\beta$ (Golan 2008, p. 96). Extracted by vcov.

hessian

T-by-T dual Hessian at the optimum; positive definite.

r.squared

Single number: the ordinary coefficient of determination.

H_signal

Numeric vector of length K: the per-coefficient signal entropies $H(\hat p_k)$.

S

Single number: the normalized signal entropy $H(\hat p)/H(p_0)$. With a uniform prior this reduces to $H/(K\log M)$ and lies in [0, 1]; with an informative prior it is measured relative to that prior and may exceed 1.

objective, value

Single number: the minimised dual objective (two names for the same value).

foc_residual

Single number: $\max|y - X\hat\beta - \hat e|$ at the optimum. Near 0 for a healthy fit; see Support feasibility.

converged

Logical: TRUE when optim converged and the first-order condition is satisfied.

convergence

Integer: the raw optim convergence code.

nu

Single number: the entropy weight used.

method

Character string naming the solver, "dual".

Z, p0, v, w0, X, y, control

The resolved inputs, retained so that summary() can refit the model for the entropy-ratio test.

terms, model, xlevels, call

The standard model components, as returned by lm.

Methods (by generic)

• coef(linreg): Estimated coefficients $\beta$.

• fitted(linreg): Fitted values $X\beta$.

• residuals(linreg): Residuals $y - X\beta$ (the estimated noise e).

• vcov(linreg): Asymptotic covariance of $\beta$, Golan (2008) p. 96 (coincides with OLS as the supports widen).

• predict(linreg): Predictions for newdata (or fitted values).

• print(linreg): Compact display.

• summary(linreg): lm()-style summary with an entropy-ratio coefficient table. For each coefficient $k$ an Entropy-Ratio (ER) test of $H_0\!:\beta_k = 0$ is computed by refitting with row $k$ of the signal support collapsed to zero and forming $\mathrm{ER}_k = 2[H^*_{\mathrm{unrestricted}} - H^*_{\mathrm{restricted},k}]$ with $H^* = \sum H(p) + \sum H(w)$; under $H_0$, $\mathrm{ER}_k \sim \chi^2_1$ (Golan 2008, Sec. 6.4/6.6). An overall ER test of $H_0\!:$ all coefficients (intercept and slopes) $= 0$ is also reported. Note this performs $K+1$ warm-started refits.

Functions

• print(summary.linreg): Print method for the summary object.

Support feasibility

The dual is bounded below only if y can be represented as $X\beta + e$ with every $\beta_k$ inside its signal support Z and every $e_t$ inside the noise support v. If it cannot, the dual is unbounded: the multipliers diverge, the signal softmax saturates and $\beta$ pins to a vertex of Z even though optim reports convergence. This is why the default v widens with $T$ – the largest of $T$ errors grows like $\sigma\sqrt{2\log T}$, so a fixed three-sigma support becomes infeasible in large samples. linreg() checks the first-order condition $y - X\beta - e = 0$ at the optimum, reports it as foc_residual, and warns (setting converged = FALSE) when it is not satisfied. The remedy is a wider v and/or Z.

References

Golan, A. (2008). Information and Entropy Econometrics - A Review and Synthesis. Foundations and Trends in Econometrics, 2(1-2), 1-145. Section 6.1. Golan, A., Judge, G. and Miller, D. (1996). Maximum Entropy Econometrics: Robust Estimation with Limited Data. Wiley. Chapter 6.

See Also

inverse_ce, inverse_noise for related information-theoretic estimators.

Examples

set.seed(1)
d <- data.frame(x1 = rnorm(100), x2 = rnorm(100))
d$y <- 2 - 1.5 * d$x1 + 0.8 * d$x2 + rnorm(100, sd = 0.7)
fit <- linreg(y ~ x1 + x2, data = d, Z = seq(-20, 20, length.out = 5))
summary(fit)

---

Stochastic-Moments GME Estimator for IV Regression

Description

Fits the linear model $y = X\beta + e$ by the (relaxed) stochastic-moments Generalized Maximum Entropy estimator of Golan (2008, pp. 89-91), identified by instrument moments $IV'(y - X\beta - e) = 0$. Each coefficient is reparameterized on a bounded signal support, $\beta_k = \sum_m z_{km} p_{km}$, and each error on a noise support, $e_i = \sum_j v_j w_{ij}$; the signal and noise entropies are maximized subject to the instrument moments, with weight $\nu \in (0,1)$. It is the instrumental-variables sibling of linreg.

Usage

linreg_iv(
  y,
  X,
  IV,
  Z,
  p0 = NULL,
  v = NULL,
  w0 = NULL,
  nu = 0.5,
  se_method = c("sandwich", "delta", "bootstrap", "none"),
  control = list(),
  boot = 200L
)

Arguments

y	Numeric response vector of length N.
X	N-by-K design matrix (include an intercept column if wanted).
IV	N-by-P instrument matrix, P >= K. For exogenous regressors, use the regressor as its own instrument.
Z	Coefficient (signal) support: a K-by-M matrix whose row k is the support for $\beta_k$ (as in linreg).
p0	Optional K-by-M signal prior (rows strictly positive, summing to 1); default uniform.
v	Optional error support: NULL (default 3-point symmetric grid on $\pm 3\,\mathrm{sd}(y)$, the three-sigma rule), a single whole number giving the number of support points, or an explicit numeric vector. Unlike linreg, this default does not need to widen with the sample size; see Choosing the supports.
w0	Optional N-by-J error prior (rows strictly positive, summing to 1); default uniform.
nu	Entropy weight in $(0,1)$ on the noise term, hence the weight on the signal is $1-\nu$; default 0.5 (equal weight on both the signal and the noise).
se_method	Standard-error method: "sandwich" (default, robust), "delta" (classical), "bootstrap" (pairs resampling), or "none" (skip; se_beta/vcov left NULL).
control	Named list merged over defaults and passed to optim (BFGS): maxit (default 1000) and reltol (default 1e-12). fnscale is forced to -1 (the dual is maximized).
boot	Number of bootstrap resamples when se_method = "bootstrap" (default 200).

Details

The concentrated dual is maximized over the Lagrange multipliers $\lambda$ (one per instrument moment): with $p_{km} \propto p0_{km}\exp(z_{km} (X'IV\lambda)_k/\nu)$ and $w_{ij} \propto w0_{ij}\exp((IV\lambda)_i v_j/(1-\nu))$, the gradient is the instrument moment $IV'(y - X\beta - e)$. Supports just-identified (ncol(IV) == ncol(X)) and over-identified (ncol(IV) > ncol(X)) systems. The instruments are standardized internally for numerical conditioning ($\beta$ is invariant to instrument scaling); lambda is reported on the original scale. coef() returns $\beta$.

As the signal support widens the GME estimate approaches the exact (2SLS-type) IV solution; narrower supports shrink $\beta$ toward the support centers.

Standard errors. Because $\hat\lambda$ solves the instrument moments, $\hat\beta = Zp(\hat\lambda)$ is a Z-estimator with sandwich covariance $\mathrm{Var}(\hat\beta) = J_\beta A^{-1}\Omega A^{-1} J_\beta'$, where $A$ is the dual Hessian and $J_\beta = \partial\beta/\partial\lambda$. The "meat" $\Omega$ is set by se_method: "sandwich" (default, robust HC0 $IV'\mathrm{diag}(r^2)IV$), "delta" (classical homoskedastic $\hat\sigma^2 IV'IV$), or "bootstrap" (pairs resampling). The robust sandwich matches a Monte-Carlo sampling SD and reduces to the 2SLS robust SE as the support widens; the classical delta assumes homoskedasticity. These are asymptotic SEs: the sampling distribution of $\hat\beta$ is support-bounded and can be skewed, so symmetric Wald intervals are approximate with weak instruments or small $n$.

Value

An object of class c("linreg_iv", "infometrics"), which is a list containing the following components:

coefficients, b_hat

Numeric vector of length K: the estimated coefficients $\hat\beta = Z\hat p$ (two names for the same object). Extracted by coef.

lambda_hat

Numeric vector of length P: the estimated Lagrange multipliers, one per instrument moment, reported on the original instrument scale.

p_hat

K-by-M matrix of estimated signal probabilities.

w_hat

N-by-J matrix of estimated noise probabilities.

e_hat, e

Numeric vector of length N: the estimated noise, $e_i = \sum_j v_j w_{ij}$ (two names for the same object).

fitted.values

Numeric vector of length N: $X\hat\beta$. Extracted by fitted.

residuals

Numeric vector of length N: $y - X\hat\beta$. Extracted by residuals.

se_beta

Numeric vector of length K: standard errors of $\hat\beta$, or NULL when se_method = "none".

vcov

K-by-K covariance matrix of $\hat\beta$, or NULL when se_method = "none". Extracted by vcov.

se_method, boot

The standard-error method used and, for the bootstrap, the number of resamples.

H_signal

Numeric vector of length K: the per-coefficient signal entropies.

S

Single number: the normalized signal entropy.

objective, value

Single number: the maximized dual objective (two names for the same value).

moment_resid

Single number: the largest absolute instrument moment at the optimum. Not normalized by $N$ – see Note on moment_resid.

converged

Logical: TRUE when optim reported convergence.

convergence

Integer: the raw optim convergence code.

method

Character string naming the solver, "dual".

Z, p0, v, w0, nu, X, y, IV

The resolved inputs.

N, K, P, M

Integers: the numbers of observations, coefficients, instruments and signal support points.

call

The matched call.

Choosing the supports

The signal support Z is the consequential input here. It must contain the true coefficients: if it does not, $\hat\beta$ is pinned to the boundary of Z (for example, a support of $\pm 1$ returns $\hat\beta_k = 1$ for a coefficient whose true value is 1.5). Widening Z moves the estimate toward the exact (2SLS-type) IV solution; narrowing it shrinks $\beta$ toward the support centers.

The error support v is far less critical, and – unlike linreg and inverse_noise – its default does not need to grow with the sample size. Those estimators impose one dual condition per observation, $y_t = x_t'\beta + e_t$, so every realized error must fit inside v; because the largest of $T$ errors grows like $\sigma\sqrt{2\log T}$, a fixed three-sigma support eventually becomes infeasible and their duals become unbounded. linreg_iv() instead imposes only the $P$ aggregate instrument moments $IV'(y - X\beta - e) = 0$, so the noise never has to absorb individual residuals and no per-observation feasibility condition arises. In practice $\hat\beta$ is nearly invariant to the width of v, and large samples pose no difficulty: the multipliers stay small and the estimate tracks the exact IV solution.

Note on moment_resid

moment_resid is $\max_p |IV'(y - X\beta - e)|_p$ evaluated on the standardized instruments. It is a sum over the $N$ observations and is not divided by $N$, so its magnitude grows with the sample size even when the fit is excellent; judge it relative to $N$ (or compare fits of the same size) rather than against a fixed threshold.

References

Golan, A. (2008). Information and Entropy Econometrics - A Review and Synthesis. Foundations and Trends in Econometrics, 2(1-2), 1-145. Pages 89-91. Golan, A., Judge, G. and Miller, D. (1996). Maximum Entropy Econometrics: Robust Estimation with Limited Data. Wiley.

See Also

linreg for the (non-IV) GME/GCE regression.

Examples

set.seed(1)
n  <- 200L
z  <- rnorm(n)                         # instrument
u  <- rnorm(n)                         # structural error
xe <- 0.7 * z + u + rnorm(n)           # endogenous regressor (corr. with u)
y  <- 1 + 1.5 * xe + u                 # true slope 1.5; OLS biased
X  <- cbind(1, xe)                     # intercept + endogenous regressor
IV <- cbind(1, z)                      # intercept is its own instrument
Z  <- matrix(c(-10, 0, 10), nrow = 2, ncol = 3, byrow = TRUE)  # signal support
fit <- linreg_iv(y, X, IV, Z)
coef(fit)                              # slope ~ 1.5 (vs OLS biased upward)

---

Construct a Symmetric Support Space

Description

Creates an M-point symmetric support vector centered at zero (or a specified center) over a given half-range. Support spaces are required inputs for the GME and GCE estimators, bounding the parameter vector beta and the error vector epsilon.

Usage

make_support(half_range, M = 5L, center = 0)

Arguments

half_range	Positive numeric. The half-width of the support interval. The support spans from center - half_range to center + half_range.
M	Positive integer. Number of support points. Must be >= 2. The default of 5 (a 5-point support) is commonly used in the GME literature.
center	Numeric. Center of the support interval (default 0).

Details

In the GME framework (Golan, 2008, Section 6.1), each parameter beta_k is reparameterized as the expected value of a random variable defined on a bounded support z_k = (z_{k1}, ..., z_{kM}). The support must contain the true parameter value. A common default is a symmetric 5-point grid with half-range set to 3 standard deviations of the OLS residuals (for the error support V) or 2 times the absolute OLS estimate (for the signal support Z).

The complexity of the GME dual model is invariant to M (the number of support points), since the real parameters are the Lagrange multipliers whose dimension equals T (the number of observations), not K*M.

Value

A numeric vector of length M giving the support points, evenly spaced from center - half_range to center + half_range.

References

Golan, A. (2008). Information and Entropy Econometrics. Foundations and Trends in Econometrics, 2(1-2), 1-145. Section 6.1.

Examples

# 5-point support for error terms: [-3, -1.5, 0, 1.5, 3]
make_support(half_range = 3, M = 5)

# 3-point support centered at a non-zero value: [0, 5, 10]
make_support(half_range = 5, M = 3, center = 5)

# 7-point support for a coefficient bounded in [-2, 2]
make_support(half_range = 2, M = 7)

---

Marginal Effects

Description

Generic for average (or per-observation) marginal effects of a fitted model.

Usage

margins(object, ...)

Arguments

object	A fitted model object.
...	Passed to methods.

Value

A method-specific object of marginal effects.

See Also

margins.multinomial_gce

---

Marginal Effects for a Mixed GCE Fit

Description

Partial derivatives of the signal probabilities with respect to each design slot, holding $\lambda$ and $\rho$ fixed:

$\partial p_{ij}/\partial x_{ijk} = (\lambda_{kj}/\nu)\,\mathrm{Var}_\theta(s)_{ij},\qquad \mathrm{Var}_\theta(s)_{ij} = \sum_m s_m^2\theta_{ijm} - p_{ij}^2.$

Usage

## S3 method for class 'mixed_gce'
margins(
  object,
  average = TRUE,
  se = FALSE,
  se_method = c("sandwich", "bootstrap"),
  B = 500L,
  ...
)

Arguments

object	A mixed_gce object.
average	Logical (default TRUE). If TRUE, return the K x J matrix of marginal effects averaged over observations; if FALSE, the full N x J x K array of per-observation derivatives.
se	Logical (default FALSE). If TRUE (requires average = TRUE), return a margins_gce object with estimates, standard errors, z-values, and p-values.
se_method	"sandwich" (default) or "bootstrap"; see Details.
B	Integer number of bootstrap resamples when se_method = "bootstrap" (default 500).
...	Unused.

Details

With se = TRUE (and average = TRUE) standard errors for the average marginal effects are attached. These are a sampling quantity, unrelated to the Fano fano_bounds (which bound classification error). se_method = "sandwich" (default) uses the robust covariance $(-H)^{-1}\hat V(-H)^{-1}$ of the dual multipliers, with $\hat V=\sum_i g_i g_i'$ the per-observation score outer products, sandwiched with a numerical Jacobian of the average effects. This matches a row bootstrap; the naive Hessian-inverse would overstate by ~7-9x for this model (the support reparameterization plus per-observation $\rho$ block make solve(-H) a poor sampling-covariance estimate), so it is not used. se_method = "bootstrap" resamples observations, refits, and takes the across-resample SD.

Value

With se = FALSE: a K x J matrix (average effects) or an N x J x K array. With se = TRUE: a margins_gce object.

See Also

fano_bounds for information-theoretic error bounds.

---

Marginal Effects for a Multinomial GCE Fit

Description

Partial derivatives of the signal probabilities with respect to each column of X:

$\partial p_{ij}/\partial x_{ik} = (p_{ij}/\nu)\,(\lambda_{kj} - \textstyle\sum_l p_{il}\lambda_{kl}).$

The reference alternative's column is included like any other; its "intercept" effect (derivative w.r.t. a constant) is rarely interpreted.

Usage

## S3 method for class 'multinomial_gce'
margins(
  object,
  average = TRUE,
  se = FALSE,
  se_method = c("sandwich", "delta", "bootstrap"),
  B = 500L,
  ...
)

Arguments

object	A multinomial_gce object.
average	Logical (default TRUE). If TRUE, return the K x J matrix of marginal effects averaged over observations; if FALSE, the full N x J x K array of per-observation derivatives.
se	Logical (default FALSE). If TRUE (requires average = TRUE), return a margins_gce object carrying the estimates, standard errors, z-values, and p-values.
se_method	"sandwich" (default), "delta", or "bootstrap"; see Details.
B	Integer number of bootstrap resamples when se_method = "bootstrap" (default 500).
...	Unused.

Details

With se = TRUE (and average = TRUE) standard errors for the average marginal effects are attached. These are a sampling quantity, unrelated to the Fano fano_bounds (which bound classification error, not sampling variance). All methods sandwich a numerical Jacobian $J$ of the average effects in the free multipliers with a covariance of $\hat\lambda$:

• se_method = "sandwich" (default) uses the robust covariance $(-H)^{-1}\hat V(-H)^{-1}$ with $\hat V=\sum_i g_i g_i'$ the per-observation dual scores. This matches a row bootstrap and is the accurate choice under the GME's finite-support regularization.

• se_method = "delta" uses the naive Hessian inverse $V=$vcov(object) (the classical Golan 2008 sec. 7.3 observed-information SE); being tied to $(-H)^{-1}$ it is conservative (~1.35x the bootstrap here).

• se_method = "bootstrap" resamples observations, refits, and takes the across-resample SD (accurate, slower).

Value

With se = FALSE: a K x J matrix (average effects) or an N x J x K array. With se = TRUE: a margins_gce object.

See Also

fano_bounds for information-theoretic error bounds.

---

Cross-Entropy Estimation of a Markov Transition/Balancing Matrix (Panel)

Description

Estimates a first-order, stationary row-stochastic matrix $P$ ($P_{kj} = k \to j$) linking simplex-valued state vectors over time, from a balanced long panel, by Cross-Entropy (Golan 2008, Section 7.7.1, eqs. 7.25-7.27). States may be one-hot membership indicators or compositional shares (a panel/time-indexed extension of matrix balancing). Data are supplied in long form (one row per unit-period) and reshaped internally; transitions $t-1 \to t$ are formed within each unit.

Usage

markov_ce(
  data,
  id,
  time,
  states,
  covariates = NULL,
  p0 = NULL,
  control = list()
)

## S3 method for class 'markov_ce'
coef(object, ...)

## S3 method for class 'markov_ce'
fitted(object, ...)

## S3 method for class 'markov_ce'
fano_bounds(object, ...)

## S3 method for class 'markov_ce'
margins(
  object,
  average = TRUE,
  se = FALSE,
  se_method = c("sandwich", "bootstrap"),
  B = 500L,
  ...
)

## S3 method for class 'markov_ce'
residuals(object, ...)

## S3 method for class 'markov_ce'
print(x, digits = max(3L, getOption("digits") - 3L), ...)

## S3 method for class 'markov_ce'
summary(object, ...)

Arguments

data	A balanced long-format data.frame, one row per unit-period.
id, time	Names of the unit-identifier and time columns.
states	Character vector of state columns; each row over them is non-negative and sums to 1. K = length(states).
covariates	Optional character vector of covariate columns (S = length(covariates)); default uses the lagged-state indicator.
p0	Optional K-by-K prior transition matrix, strictly positive with rows summing to 1. Defaults to uniform.
control	A named list merged over the defaults and passed to optim (BFGS): maxit (default 1000) and reltol (default 1e-12). fnscale is forced to -1 (the dual is maximized).
object, x	A markov_ce object.
...	Unused.
average	Logical; if TRUE (default) return the origin-weighted average marginal effect (S x K); if FALSE the full K x K x S array of per-origin effects.
se	Logical (default FALSE); if TRUE (with average = TRUE) attach standard errors.
se_method	"sandwich" (default, robust and analytic) or "bootstrap" (resample units, re-solve). See margins.markov_gce.
B	Integer bootstrap resamples when se_method = "bootstrap" (default 500).
digits	Significant digits to print.

Details

Optional covariates z enter lead/lagged like the states: per transition the lagged states $X$ pair with lagged covariates $Z_2$ and the lead states $Y$ with lead covariates $Z_1$, forming the cross moments $A = \sum X' Z_2$ (K x S) and $M = \sum Z_1' Y$ (S x K), $S =$number of covariates. The concentrated dual

$\ell(\lambda) = \sum_{s,j} M_{sj}\lambda_{sj} - \sum_k \log\Omega_k(\lambda), \quad \hat P_{kj} = p0_{kj}\exp(\textstyle\sum_s A_{ks}\lambda_{sj})/\Omega_k$

is maximized over the S x K multipliers. With no covariates the conditioning is the lagged-state indicator ($S = K$, transition counts).

When the conditioning is full rank (A = X'X}, the common indicator case), the moment condition \eqn{A'P = M exactly determines $P$ as the empirical/regression transition matrix, so the prior is irrelevant for identified origin states. Only under-identified rows ($\mathrm{rank}(A) < K$, e.g. a state never observed as an origin) fall back to the prior (maximum entropy when p0 is uniform); see the rank warning and the normalized entropy S in summary().

Standard errors of $\lambda$ (se_lambda) use the analytic Hessian of the CE log-partition (Golan 2008, eq. 4.7), $H_{(s,j),(s',j')} = \sum_k A_{ks} A_{ks'} P_{kj}(\delta_{jj'} - P_{kj'})$ – the covariance of the moment functions under $\hat P$. Note that $H^{-1}$ alone overstates the SE (~$\sqrt{n}$), because the moments $A, M$ are cross-moment sums ($H = n\,\hat V$); the reported SE is the robust sandwich $H^{-1}\hat V H^{-1}$ (Golan 2008, Sec. 3.3), which matches a unit bootstrap. Because $\lambda$ is identified only up to a per-conditioning additive shift, se_lambda is relative to a reference destination (lambda_ref = 1, that column NA); a (near-)boundary fit ($P_{kj}\in\{0,1\}$) yields NA (an infinite log-odds SE).

Value

An object of class c("markov_ce", "infometrics"), which is a list containing the following components:

p, p_hat

K-by-K row-stochastic transition matrix, $\hat P_{kj} = k \to j$ (two names for the same object). Extracted by coef.

lambda, lambda_hat

S-by-K matrix of estimated Lagrange multipliers, reported on the original covariate scale (two names for the same object).

se

K-by-K matrix of large-sample multinomial standard errors of $\hat P$.

se_lambda

S-by-K matrix of analytic-Hessian sandwich standard errors of $\lambda$. The reference-destination column is NA, as are entries at a (near-)boundary fit.

lambda_ref

Integer: the index of the reference destination state used to normalize se_lambda.

vcov_lambda

Covariance matrix of the free (non-reference) elements of $\lambda$, or NULL if the sandwich could not be formed.

n_from

Numeric vector of length K: the total weight of each origin state across transitions.

states

Character vector of length K: the state labels.

se_data

List of the cross moments, lead/lagged design blocks and scaling used by margins.markov_ce, or NULL when no covariates were supplied.

n_transitions

Integer: the number of unit-period transitions, $N(T-1)$.

state_type

Character string, "indicator" or "shares", describing the supplied state columns.

fitted.values

Matrix of fitted lead states, $X\hat P$, one row per transition. Extracted by fitted.

residuals

Matrix of the same dimension: $Y - X\hat P$. Extracted by residuals.

entropy, H_signal

Single number: the mean row entropy of $\hat P$ (two names for the same value).

S

Single number: the prior-relative normalized entropy, the mean row entropy of $\hat P$ divided by that of p0.

moment_residual

Single number: $\max|M - A'\hat P|$ at the optimum.

objective, value

Single number: the maximized dual objective (two names for the same value).

prior

The resolved K-by-K prior p0.

converged

Logical: TRUE when optim reported convergence.

convergence

Integer: the raw optim convergence code.

method

Character string naming the solver, "dual".

call

The matched call.

References

Golan, A. (2008). Information and Entropy Econometrics - A Review and Synthesis. Foundations and Trends in Econometrics, 2(1-2), 1-145. Section 7.7.1.

See Also

matrix_ce, matrix_gce

Examples

set.seed(1)
P_true <- matrix(c(.6,.3,.1, .2,.5,.3, .1,.3,.6), nrow = 3, byrow = TRUE)
N <- 150L; Tn <- 5L; K <- 3L
s  <- sample(K, N, replace = TRUE)
id <- rep(seq_len(N), each = Tn); tm <- rep(seq_len(Tn), N)
state <- integer(N * Tn)
for (t in seq_len(Tn)) {
  state[tm == t] <- s
  s <- vapply(s, function(k) sample(K, 1L, prob = P_true[k, ]), integer(1))
}
oh <- diag(K)[state, ]
panel <- data.frame(i = id, t = tm, y1 = oh[, 1], y2 = oh[, 2], y3 = oh[, 3])
fit <- markov_ce(panel, id = "i", time = "t", states = c("y1", "y2", "y3"))
coef(fit)

---

Generalized Cross-Entropy Estimation of a Markov Matrix with Noisy Moments

Description

Estimates a row-stochastic $K \times K$ matrix $P$ linking simplex- valued state vectors over time, conditioning on covariates, by Generalized Cross-Entropy from a balanced long panel (Golan 2008, Section 7.7.1, eqs. 7.25a and 7.28). Unlike markov_ce (which imposes the moments exactly), each conditional moment carries an additive error $\varepsilon_{itj} = \sum_m w_{itjm} v_m$ on a support $v$ symmetric on $[-1, 1]$. This makes over-identified moment systems – in particular time-varying covariates such as household income – feasible, where the exact problem is unbounded.

Usage

markov_gce(
  data,
  id,
  time,
  states,
  covariates,
  v = NULL,
  nu = 0.5,
  p0 = NULL,
  w0 = NULL,
  control = list()
)

## S3 method for class 'markov_gce'
coef(object, ...)

## S3 method for class 'markov_gce'
fano_bounds(object, ...)

## S3 method for class 'markov_gce'
margins(
  object,
  average = TRUE,
  se = FALSE,
  se_method = c("sandwich", "bootstrap"),
  B = 500L,
  ...
)

## S3 method for class 'markov_gce'
fitted(object, ...)

## S3 method for class 'markov_gce'
residuals(object, ...)

## S3 method for class 'markov_gce'
print(x, digits = max(3L, getOption("digits") - 3L), ...)

## S3 method for class 'markov_gce'
summary(object, ...)

Arguments

data	A balanced long-format data.frame, one row per unit-period.
id, time	Names of the unit-identifier and time columns.
states	Character vector of state columns; each row over them is non-negative and sums to 1. K = length(states).
covariates	Character vector of covariate columns (S = length(covariates)); required. Time-varying covariates are supported.
v	Numeric error support, symmetric about 0 with at least two points; default c(-1, 0, 1) (Golan's symmetric $[-1,1]$ support).
nu	Signal/noise weight in $(0,1)$; default 0.5. Smaller nu penalizes the noise more strongly (driving $\hat\varepsilon \to 0$, toward the exact markov_ce); larger nu lets the noise absorb more of the moment discrepancy.
p0	Optional K-by-K prior transition matrix, strictly positive with rows summing to 1. Defaults to uniform.
w0	Optional error prior on the support (length length(v), positive, summing to 1). Defaults to uniform (mean-zero error prior).
control	A named list merged over the defaults and passed to optim (BFGS): maxit (default 1000) and reltol (default 1e-12). fnscale is forced to -1 (the dual is maximized).
object, x	A markov_gce object.
...	Unused.
average	Logical; if TRUE (default) return the origin-weighted average marginal effect (S x K); if FALSE the full K x K x S array of per-origin effects.
se	Logical (default FALSE); if TRUE (with average = TRUE) attach standard errors.
se_method	"sandwich" (default) or "bootstrap". The naive Hessian-inverse delta overstates ~4-5x here; the robust sandwich $(-H)^{-1}\hat V(-H)^{-1}$ matches the "bootstrap" (which resamples units and re-solves the dual).
B	Integer bootstrap resamples when se_method = "bootstrap" (default 500).
digits	Significant digits to print.

Details

The $\nu$-weighted concentrated dual maximized over the $S \times K$ multipliers $\lambda$ is

$\ell(\lambda) = \sum_{t\ge 2}\sum_j\sum_{i,s} y_{itj} z_{its}\lambda_{sj} - \nu \sum_k \log\Omega_k(\lambda) - (1-\nu)\sum_{i,t,j}\log\Psi_{itj}(\lambda),$

with $\hat p_{kj}\propto p0_{kj}\exp(\nu^{-1}\sum A_{ks}\lambda_{sj})$ and $\hat w_{itjm}\propto w0_{m}\exp((1-\nu)^{-1}\sum_s z_{its}v_m\lambda_{sj})$. The first-order conditions are the noisy moments $M - A'\hat P - \sum_{i,t} z\,\hat\varepsilon = 0$. The exact estimator markov_ce is recovered in two limits: as the support $v \to 0$ (the slack vanishes mechanically), and as $\nu \to 0$ (the noise-divergence term dominates, forcing $\hat\varepsilon \to 0$). Conversely, a wide support or $\nu \to 1$ lets the errors absorb the moments and pulls $\hat P$ toward the prior.

Standard errors of $\lambda$ (se_lambda) use the analytic Hessian of the dual (Golan 2008, eq. 4.7): the signal covariance $\nu^{-1}\sum_k A_{ks}A_{ks'}P_{kj}(\delta_{jj'}-P_{kj'})$ plus the noise covariance $(1-\nu)^{-1}\sum_{i,t} z_{its} z_{its'}\mathrm{Var}_w(v)_{itj}$ (block-diagonal in the destination $j$). The noise term makes the Hessian full rank – $\lambda$ is identified, so no reference normalization is needed (unlike markov_ce). The reported SE is the robust unit-clustered sandwich $H^{-1}\hat V H^{-1}$ (Golan 2008, Sec. 3.3), which matches a unit bootstrap; the naive $H^{-1}$ and the non-clustered sandwich are conservative because covariate scores are correlated within a unit.

Value

An object of class c("markov_gce", "infometrics"), which is a list containing the following components:

p, p_hat

K-by-K row-stochastic transition matrix, $\hat P_{kj} = k \to j$ (two names for the same object). Extracted by coef.

lambda, lambda_hat

S-by-K matrix of estimated Lagrange multipliers, reported on the original covariate scale (two names for the same object).

epsilon

Matrix of fitted per-transition errors $\hat\varepsilon_{itj}$, one row per transition and one column per destination state.

se

K-by-K matrix of large-sample multinomial standard errors of $\hat P$.

se_lambda

S-by-K matrix of analytic-Hessian, unit-clustered sandwich standard errors of $\lambda$. All entries are reported: the noise identifies $\lambda$, so no reference normalization is needed.

vcov_lambda

Covariance matrix of $\lambda$, or NULL if the sandwich could not be formed.

states

Character vector of length K: the state labels.

se_data

List of the cross moments, lead/lagged design blocks and scaling used by margins.markov_gce.

n_transitions

Integer: the number of unit-period transitions, $N(T-1)$.

n_from

Numeric vector of length K: the total weight of each origin state across transitions.

nu

Single number: the signal/noise weight used.

v, w0

The resolved error support and its prior.

H_p, H_signal, entropy

Single number: the mean row entropy of $\hat P$ (three names for the same value).

H_w

Single number: the mean noise entropy $H(\hat w)$.

S

Single number: the prior-relative normalized entropy, the mean row entropy of $\hat P$ divided by that of p0.

fitted.values

Matrix of fitted lead states, $X\hat P$, one row per transition. Extracted by fitted.

residuals

Matrix of the same dimension: $Y - X\hat P$. Extracted by residuals.

moment_residual

Single number: $\max|M - A'\hat P|$, the exact-moment gap. Expected to be nonzero – it is absorbed by the noise.

foc_residual

Single number: $\max|M - A'\hat P - \sum z\,\hat\varepsilon|$, the noisy first-order condition. Near 0 at a healthy optimum.

objective, value

Single number: the maximized dual objective (two names for the same value).

prior

The resolved K-by-K prior p0.

converged

Logical: TRUE when optim reported convergence.

convergence

Integer: the raw optim convergence code.

method

Character string naming the solver, "dual".

call

The matched call.

Choosing nu and v

There is no universally correct default. The noise enters per observation and the literal support $v = (-1,0,1)$ is permissive relative to share-valued moment contributions, so at a wide support or large nu the errors can absorb the signal and pull $\hat P$ toward the prior. In controlled recovery experiments the estimate approaches the exact markov_ce solution monotonically as the support narrows ($v \to 0$) or as nu decreases toward 0, though an ultra-narrow support can make the over-identified exact problem numerically unbounded (non-convergence). Users should report sensitivity to nu and v for their data rather than rely on a single setting.

References

Golan, A. (2008). Information and Entropy Econometrics - A Review and Synthesis. Foundations and Trends in Econometrics, 2(1-2), 1-145. Section 7.7.1, eqs. (7.25a), (7.28).

See Also

markov_ce, matrix_gce

Examples

set.seed(1)
P_true <- matrix(c(.6,.3,.1, .2,.5,.3, .1,.3,.6), nrow = 3, byrow = TRUE)
N <- 120L; Tn <- 5L; K <- 3L
Y <- matrix(rgamma(N * K, 1), N, K); Y <- Y / rowSums(Y); inc <- rnorm(N)
id <- rep(seq_len(N), each = Tn); tm <- rep(seq_len(Tn), N)
sh <- matrix(0, N * Tn, K); income <- numeric(N * Tn); infl <- numeric(N * Tn)
for (t in seq_len(Tn)) {
  idx <- which(tm == t); inc <- 0.8 * inc + rnorm(N, 0, 0.5)
  sh[idx, ] <- Y; income[idx] <- inc; infl[idx] <- 0.02 * t
  Y <- (Y %*% P_true) * exp(matrix(rnorm(N * K, 0, 0.03), N, K)); Y <- Y / rowSums(Y)
}
panel <- data.frame(i = id, t = tm, y1 = sh[, 1], y2 = sh[, 2], y3 = sh[, 3],
                    income = income, infl = infl)
fit <- markov_gce(panel, id = "i", time = "t", states = c("y1","y2","y3"),
                  covariates = c("income", "infl"), v = c(-0.2, 0, 0.2),
                  nu = 0.2)
coef(fit)

---