# Bayesian Single-Index Regression with B-Spline Link and One-to-One Polar Transformation Fits a single-index model $Y_i \sim \mathcal{N}(f(X_i'\theta), \sigma^2), i = 1,\cdots,n$ where the link $f(\cdot)$ is represented by B-spline link function and the index vector $\theta$ is parameterized on the unit sphere via a one-to-one polar transformation. ```r bsPolar( formula, data, prior = NULL, init = NULL, monitors = NULL, niter = 10000, nburnin = 1000, thin = 1, nchain = 1, setSeed = FALSE ) bsPolar_setup( formula, data, prior = NULL, init = NULL, monitors = NULL, niter = 10000, nburnin = 1000, thin = 1, nchain = 1, setSeed = FALSE ) ``` ## Arguments - `formula`: an object of class formula . Interaction term is not allowed for single-index model. - `data`: an data frame. - `prior`: Optional named list of prior settings. Further descriptions are in Details section. - `init`: Optional named list of initial values. If the values are not assigned, they are randomly sampled from prior or designated value. Further descriptions are in Details section. - `monitors`: Optional character vector of additional monitor nodes. To check the names of the nodes, fit the `model_setup` function and then inspect the variable names stored in the model object using getVarMonitor . - `niter`: Integer. Total MCMC iterations (default `10000`). - `nburnin`: Integer. Burn-in iterations (default `1000`). - `thin`: Integer. Thinning for monitors (default `1`). - `nchain`: Integer. Number of MCMC chains (default `1`). - `setSeed`: Logical or numeric argument. Further details are provided in runMCMC `setSeed` argument. ## Returns A `list` typically containing: - **`coefficients`**: Mean estimates of index vector. Return list of `model_setup` does not include it. - **`ses`**: Mean standard error of index vector. Return list of `model_setup` does not include it. - **`residuals`**: Residuals with mean estimates of fitted values. Return list of `model_setup` does not include it. - **`fitted.values`**: Mean estimates of fitted values. Return list of `model_setup` does not include it. - **`linear.predictors`**: Mean estimates of single-index values. Return list of `model_setup` does not include it. - **`gof`**: Goodness-of-fit. Return list of `model_setup` function does not include it. - **`samples`**: Posterior draws of variables for computing fitted values of the model, including $\theta$, $\sigma^2$. Return list of `model_setup` does not include it. - **`input`**: List of data used in modeling, formula, input values for prior and initial values, and computation time without compiling. - **`defModel`**: Nimble model object. - **`defSampler`**: Nimble sampler object. - **`modelName`**: Name of the model. ## Details Model The single–index model uses a $m$-order polynomial spline with $k$ interior knots as follows: $$ f(t) = \sum_{j=1}^{m+k} B_j(t)\,\beta_j $$ on $[a, b]$ with $t_i = X_i' \theta, i = 1,\cdots, n$ and $\|\theta\|_2 = 1$. $\{\beta_j\}_{j=1}^{m+k}$ are spline coefficient and $a_\theta$, $ b_\theta$ are boundary knots where $a_\theta = min(t_i, i = 1, \cdots, n) - \delta$, and $b_\theta = max(t_i, i = 1,\cdots, n) + \delta$. $\theta$ lies on the unit sphere ($\|\theta\|_2=1$) to ensure identifiability and is parameterized via a one-to-one polar transformation with angle $\psi_1,...,\psi_{p-1}$, where $p$ is the dimension of predictor. Sampling is performed on the angular parameters $\\psi$ defining the index vector. The mapping is $$ \begin{aligned}\theta_1 &= \sin(\psi_1),\\\theta_i &= \Big(\prod_{j=1}^{i-1}\cos(\psi_j)\Big)\sin(\psi_i), \quad i=2,\dots,p-1,\\\theta_p &= \prod_{j=1}^{p-1}\cos(\psi_j).\end{aligned} $$ The vector is then scaled to unit length. Priors * $\psi$ is $p-1$ dimension of polar angle of index vector and re-scaled Beta distribution on $[0, \pi]$ is assigned. * Conditioned on $\theta$ and $\sigma^2$, the link coefficients $\beta = (\beta_1,\ldots,\beta_{m+k})^\top$ follow normal distribution with estimated mean vector $\hat{\beta}_{\theta} = (X_{\theta}'X_{\theta})^{-1}X_{\theta}'Y$ and covariance $\sigma^2 (X_{\theta}^\top X_{\theta})^{-1}$, where $X_{\theta}$ is the B-spline basis design matrix. * Inverse gamma prior on $\sigma^2$ with shape parameter $a_\sigma$ and rate parameter $b_\sigma$. Sampling Samplers are automatically assigned by nimble. Prior hyper-parameters These are the prior hyper-parameters set in the function. You can define new values for each parameter in prior_param . 1. Index vector: Only shape parameter `index_psi_alpha` of $p-1$ dimension vector is needed since rate parameters is computed to satisfy $\theta_{j, \text{MAP}}$. By default, the shape parameter for each element of the polar vector is set to `5000`. 2. Link function: B-spline basis and coefficient of B-spline setup. * `basis`: For the basis of B-spline, `link_basis_df` is the number of basis functions (default `21`), `link_basis_degree` is the spline degree (default `2`) and `link_basis_delta` is a small jitter for boundary-knot spacing control (default `0.001`). * `beta`: For the coefficient of B-spline, multivariate normal prior is assigned with mean `link_beta_mu`, and covariance `link_beta_cov`. By default, $\mathcal{N}_p\!\big(0, \mathrm{I}_p\big)$ is assigned. 3. Error variance (`sigma2`): An Inverse gamma prior is assigned to $\sigma^2$ where `sigma2_shape` is shape parameter and `sigma2_rate` is rate parameter of inverse gamma distribution. (default `sigma2_shape = 0.001, sigma2_rate = 100`) Initial values These are the initial values set in the function. You can define new values for each initial value in init_param . 1. Index vector: Initial vector of polar angle `index_psi` ($p-1$ dimension). Each element of angle is between 0 and $\pi$. By default, it is randomly draw from uniform distribution. 2. Link function: Initial spline coefficients(`link_beta`). By default, $\big(X_{\theta}^\top X_{\theta} + \rho\, \mathrm{I}\big)^{-1} X_{\theta}^\top Y$ is computed, where $X_{\theta}$ is the B-spline basis design matrix. 3. Error variance (`sigma2`): Initial scalar error variance (default `0.01`). ## Examples ```r set.seed(123) n <- 200; d <- 4 theta <- c(2, 1, 1, 1); theta <- theta / sqrt(sum(theta^2)) f <- function(u) u^2 * exp(u) sigma <- 0.5 X <- matrix(runif(n * d, -1, 1), nrow = n) index_vals <- as.vector(X %*% theta) y <- f(index_vals) + rnorm(n, 0, sigma) simdata <- data.frame(x = X, y = y) colnames(simdata) <- c(paste0("X", 1:4), "y") # One tool version fit1 <- bsPolar(y ~ ., data = simdata, niter = 5000, nburnin = 1000, nchain = 1) # Split version models <- bsPolar_setup(y ~ ., data = simdata) Ccompile <- compileModelAndMCMC(models) nimSampler <- get_sampler(Ccompile) initList <- getInit(models) mcmc.out <- runMCMC(nimSampler, niter = 5000, nburnin = 1000, thin = 1, nchains = 1, setSeed = TRUE, inits = initList, samplesAsCodaMCMC = TRUE) fit2 <- as_bsim(models, mcmc.out) summary(fit2) ``` ## References Antoniadis, A., Grégoire, G., & McKeague, I. W. (2004). Bayesian estimation in single-index models. **Statistica Sinica**, 1147-1164. Hornik, K., & Grün, B. (2014). movMF: an R package for fitting mixtures of von Mises-Fisher distributions. **Journal of Statistical Software**, 58, 1-31. Dhara, K., Lipsitz, S., Pati, D., & Sinha, D. (2019). A new Bayesian single index model with or without covariates missing at random. **Bayesian analysis**, 15(3), 759.

bsPolar function