hskkw | R Documentation |
Computes the analytic 4x4 Hessian matrix (matrix of second partial derivatives)
of the negative log-likelihood function for the Kumaraswamy-Kumaraswamy (kkw)
distribution with parameters alpha
(\alpha
), beta
(\beta
), delta
(\delta
), and lambda
(\lambda
).
This distribution is the special case of the Generalized Kumaraswamy (GKw)
distribution where \gamma = 1
. The Hessian is useful for estimating
standard errors and in optimization algorithms.
hskkw(par, data)
par |
A numeric vector of length 4 containing the distribution parameters
in the order: |
data |
A numeric vector of observations. All values must be strictly between 0 and 1 (exclusive). |
This function calculates the analytic second partial derivatives of the
negative log-likelihood function based on the kkw log-likelihood
(\gamma=1
case of GKw, see llkkw
):
\ell(\theta | \mathbf{x}) = n[\ln(\delta+1) + \ln(\lambda) + \ln(\alpha) + \ln(\beta)]
+ \sum_{i=1}^{n} [(\alpha-1)\ln(x_i) + (\beta-1)\ln(v_i) + (\lambda-1)\ln(w_i) + \delta\ln(z_i)]
where \theta = (\alpha, \beta, \delta, \lambda)
and intermediate terms are:
v_i = 1 - x_i^{\alpha}
w_i = 1 - v_i^{\beta} = 1 - (1-x_i^{\alpha})^{\beta}
z_i = 1 - w_i^{\lambda} = 1 - [1-(1-x_i^{\alpha})^{\beta}]^{\lambda}
The Hessian matrix returned contains the elements - \frac{\partial^2 \ell(\theta | \mathbf{x})}{\partial \theta_i \partial \theta_j}
for \theta_i, \theta_j \in \{\alpha, \beta, \delta, \lambda\}
.
Key properties of the returned matrix:
Dimensions: 4x4.
Symmetry: The matrix is symmetric.
Ordering: Rows and columns correspond to the parameters in the order
\alpha, \beta, \delta, \lambda
.
Content: Analytic second derivatives of the negative log-likelihood.
This corresponds to the relevant submatrix of the 5x5 GKw Hessian (hsgkw
)
evaluated at \gamma=1
. The exact analytical formulas are implemented directly.
Returns a 4x4 numeric matrix representing the Hessian matrix of the
negative log-likelihood function, -\partial^2 \ell / (\partial \theta_i \partial \theta_j)
,
where \theta = (\alpha, \beta, \delta, \lambda)
.
Returns a 4x4 matrix populated with NaN
if any parameter values are
invalid according to their constraints, or if any value in data
is
not in the interval (0, 1).
Lopes, J. E.
Cordeiro, G. M., & de Castro, M. (2011). A new family of generalized distributions. Journal of Statistical Computation and Simulation
Kumaraswamy, P. (1980). A generalized probability density function for double-bounded random processes. Journal of Hydrology, 46(1-2), 79-88.
hsgkw
(parent distribution Hessian),
llkkw
(negative log-likelihood for kkw),
grkkw
(gradient for kkw),
dkkw
(density for kkw),
optim
,
hessian
(for numerical Hessian comparison).
# Assuming existence of rkkw, llkkw, grkkw, hskkw functions for kkw
# Generate sample data
set.seed(123)
true_par_kkw <- c(alpha = 2, beta = 3, delta = 1.5, lambda = 0.5)
if (exists("rkkw")) {
sample_data_kkw <- rkkw(100, alpha = true_par_kkw[1], beta = true_par_kkw[2],
delta = true_par_kkw[3], lambda = true_par_kkw[4])
} else {
sample_data_kkw <- rgkw(100, alpha = true_par_kkw[1], beta = true_par_kkw[2],
gamma = 1, delta = true_par_kkw[3], lambda = true_par_kkw[4])
}
# --- Find MLE estimates ---
start_par_kkw <- c(1.5, 2.5, 1.0, 0.6)
mle_result_kkw <- stats::optim(par = start_par_kkw,
fn = llkkw,
gr = if (exists("grkkw")) grkkw else NULL,
method = "BFGS",
hessian = TRUE,
data = sample_data_kkw)
# --- Compare analytical Hessian to numerical Hessian ---
if (mle_result_kkw$convergence == 0 &&
requireNamespace("numDeriv", quietly = TRUE) &&
exists("hskkw")) {
mle_par_kkw <- mle_result_kkw$par
cat("\nComparing Hessians for kkw at MLE estimates:\n")
# Numerical Hessian of llkkw
num_hess_kkw <- numDeriv::hessian(func = llkkw, x = mle_par_kkw, data = sample_data_kkw)
# Analytical Hessian from hskkw
ana_hess_kkw <- hskkw(par = mle_par_kkw, data = sample_data_kkw)
cat("Numerical Hessian (kkw):\n")
print(round(num_hess_kkw, 4))
cat("Analytical Hessian (kkw):\n")
print(round(ana_hess_kkw, 4))
# Check differences
cat("Max absolute difference between kkw Hessians:\n")
print(max(abs(num_hess_kkw - ana_hess_kkw)))
# Optional: Use analytical Hessian for Standard Errors
# tryCatch({
# cov_matrix_kkw <- solve(ana_hess_kkw)
# std_errors_kkw <- sqrt(diag(cov_matrix_kkw))
# cat("Std. Errors from Analytical kkw Hessian:\n")
# print(std_errors_kkw)
# }, error = function(e) {
# warning("Could not invert analytical kkw Hessian: ", e$message)
# })
} else {
cat("\nSkipping kkw Hessian comparison.\n")
cat("Requires convergence, 'numDeriv' package, and function 'hskkw'.\n")
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.