1% Generated by roxygen2: do not edit by hand 2% Please edit documentation in R/MCMCoprobitChange.R 3\name{MCMCoprobitChange} 4\alias{MCMCoprobitChange} 5\title{Markov Chain Monte Carlo for Ordered Probit Changepoint Regression Model} 6\usage{ 7MCMCoprobitChange( 8 formula, 9 data = parent.frame(), 10 m = 1, 11 burnin = 1000, 12 mcmc = 1000, 13 thin = 1, 14 tune = NA, 15 verbose = 0, 16 seed = NA, 17 beta.start = NA, 18 gamma.start = NA, 19 P.start = NA, 20 b0 = NULL, 21 B0 = NULL, 22 a = NULL, 23 b = NULL, 24 marginal.likelihood = c("none", "Chib95"), 25 gamma.fixed = 0, 26 ... 27) 28} 29\arguments{ 30\item{formula}{Model formula.} 31 32\item{data}{Data frame.} 33 34\item{m}{The number of changepoints.} 35 36\item{burnin}{The number of burn-in iterations for the sampler.} 37 38\item{mcmc}{The number of MCMC iterations after burnin.} 39 40\item{thin}{The thinning interval used in the simulation. The number of 41MCMC iterations must be divisible by this value.} 42 43\item{tune}{The tuning parameter for the Metropolis-Hastings step. Default 44of NA corresponds to a choice of 0.05 divided by the number of categories in 45the response variable.} 46 47\item{verbose}{A switch which determines whether or not the progress of the 48sampler is printed to the screen. If \code{verbose} is greater than 0 the 49iteration number, the \eqn{\beta} vector, and the error variance are 50printed to the screen every \code{verbose}th iteration.} 51 52\item{seed}{The seed for the random number generator. If NA, the Mersenne 53Twister generator is used with default seed 12345; if an integer is passed 54it is used to seed the Mersenne twister. The user can also pass a list of 55length two to use the L'Ecuyer random number generator, which is suitable 56for parallel computation. The first element of the list is the L'Ecuyer 57seed, which is a vector of length six or NA (if NA a default seed of 58\code{rep(12345,6)} is used). The second element of list is a positive 59substream number. See the MCMCpack specification for more details.} 60 61\item{beta.start}{The starting values for the \eqn{\beta} vector. 62This can either be a scalar or a column vector with dimension equal to the 63number of betas. The default value of of NA will use the MLE estimate of 64\eqn{\beta} as the starting value. If this is a scalar, that value 65will serve as the starting value mean for all of the betas.} 66 67\item{gamma.start}{The starting values for the \eqn{\gamma} vector. 68This can either be a scalar or a column vector with dimension equal to the 69number of gammas. The default value of of NA will use the MLE estimate of 70\eqn{\gamma} as the starting value. If this is a scalar, that value 71will serve as the starting value mean for all of the gammas.} 72 73\item{P.start}{The starting values for the transition matrix. A user should 74provide a square matrix with dimension equal to the number of states. By 75default, draws from the \code{Beta(0.9, 0.1)} are used to construct a proper 76transition matrix for each raw except the last raw.} 77 78\item{b0}{The prior mean of \eqn{\beta}. This can either be a scalar 79or a column vector with dimension equal to the number of betas. If this 80takes a scalar value, then that value will serve as the prior mean for all 81of the betas.} 82 83\item{B0}{The prior precision of \eqn{\beta}. This can either be a 84scalar or a square matrix with dimensions equal to the number of betas. If 85this takes a scalar value, then that value times an identity matrix serves 86as the prior precision of beta. Default value of 0 is equivalent to an 87improper uniform prior for beta.} 88 89\item{a}{\eqn{a} is the shape1 beta prior for transition probabilities. 90By default, the expected duration is computed and corresponding a and b 91values are assigned. The expected duration is the sample period divided by 92the number of states.} 93 94\item{b}{\eqn{b} is the shape2 beta prior for transition probabilities. 95By default, the expected duration is computed and corresponding a and b 96values are assigned. The expected duration is the sample period divided by 97the number of states.} 98 99\item{marginal.likelihood}{How should the marginal likelihood be calculated? 100Options are: \code{none} in which case the marginal likelihood will not be 101calculated, and \code{Chib95} in which case the method of Chib (1995) is 102used.} 103 104\item{gamma.fixed}{1 if users want to constrain \eqn{\gamma} values 105to be constant. By default, \eqn{\gamma} values are allowed to vary 106across regimes.} 107 108\item{...}{further arguments to be passed} 109} 110\value{ 111An mcmc object that contains the posterior sample. This object can 112be summarized by functions provided by the coda package. The object 113contains an attribute \code{prob.state} storage matrix that contains the 114probability of \eqn{state_i} for each period, the log-likelihood of 115the model (\code{loglike}), and the log-marginal likelihood of the model 116(\code{logmarglike}). 117} 118\description{ 119This function generates a sample from the posterior distribution of an 120ordered probit regression model with multiple parameter breaks. The function 121uses the Markov chain Monte Carlo method of Chib (1998). The user supplies 122data and priors, and a sample from the posterior distribution is returned as 123an mcmc object, which can be subsequently analyzed with functions provided 124in the coda package. 125} 126\details{ 127\code{MCMCoprobitChange} simulates from the posterior distribution of an 128ordinal probit regression model with multiple parameter breaks. The 129simulation of latent states is based on the linear approximation method 130discussed in Park (2011). 131 132The model takes the following form: 133 134\deqn{\Pr(y_t = 1) = \Phi(\gamma_{c, m} - x_i'\beta_m) - \Phi(\gamma_{c-1, m} - x_i'\beta_m)\;\; m = 1, \ldots, M} 135 136Where \eqn{M} is the number of states, and \eqn{\gamma_{c, m}} and 137\eqn{\beta_m} are paramters when a state is \eqn{m} at \eqn{t}. 138 139We assume Gaussian distribution for prior of \eqn{\beta}: 140 141\deqn{\beta_m \sim \mathcal{N}(b_0,B_0^{-1}),\;\; m = 1, \ldots, M} 142 143And: 144 145\deqn{p_{mm} \sim \mathcal{B}eta(a, b),\;\; m = 1, \ldots, M} 146 147Where \eqn{M} is the number of states. 148 149Note that when the fitted changepoint model has very few observations in any 150of states, the marginal likelihood outcome can be ``nan," which indicates 151that too many breaks are assumed given the model and data. 152} 153\examples{ 154 155set.seed(1909) 156N <- 200 157x1 <- rnorm(N, 1, .5); 158 159## set a true break at 100 160z1 <- 1 + x1[1:100] + rnorm(100); 161z2 <- 1 -0.2*x1[101:200] + rnorm(100); 162z <- c(z1, z2); 163y <- z 164 165## generate y 166y[z < 1] <- 1; 167y[z >= 1 & z < 2] <- 2; 168y[z >= 2] <- 3; 169 170## inputs 171formula <- y ~ x1 172 173## fit multiple models with a varying number of breaks 174out1 <- MCMCoprobitChange(formula, m=1, 175 mcmc=100, burnin=100, thin=1, tune=c(.5, .5), verbose=100, 176 b0=0, B0=0.1, marginal.likelihood = "Chib95") 177out2 <- MCMCoprobitChange(formula, m=2, 178 mcmc=100, burnin=100, thin=1, tune=c(.5, .5, .5), verbose=100, 179 b0=0, B0=0.1, marginal.likelihood = "Chib95") 180 181## Do model comparison 182## NOTE: the chain should be run longer than this example! 183BayesFactor(out1, out2) 184 185## draw plots using the "right" model 186plotState(out1) 187plotChangepoint(out1) 188 189} 190\references{ 191Jong Hee Park. 2011. ``Changepoint Analysis of Binary and 192Ordinal Probit Models: An Application to Bank Rate Policy Under the Interwar 193Gold Standard." \emph{Political Analysis}. 19: 188-204. <doi:10.1093/pan/mpr007> 194 195Andrew D. Martin, Kevin M. Quinn, and Jong Hee Park. 2011. ``MCMCpack: 196Markov Chain Monte Carlo in R.'', \emph{Journal of Statistical Software}. 19742(9): 1-21. \doi{10.18637/jss.v042.i09}. 198 199Siddhartha Chib. 1998. ``Estimation and comparison of multiple change-point 200models.'' \emph{Journal of Econometrics}. 86: 221-241. 201} 202\seealso{ 203\code{\link{plotState}}, \code{\link{plotChangepoint}} 204} 205\keyword{models} 206