### 10.1 Simple Seemingly Unrelated Regressions (SUR)

Mathematical representation

$\begin{array}{ccccc}\hfill {y}_{1i}\hfill & \hfill =\hfill & \hfill {\mathbf{x}}_{1i}^{\prime }{\beta }_{1}\hfill & \hfill +\hfill & \hfill {𝜀}_{1i}\hfill \\ \hfill {y}_{2i}\hfill & \hfill =\hfill & \hfill {\mathbf{x}}_{2i}^{\prime }{\beta }_{2}\hfill & \hfill +\hfill & \hfill {𝜀}_{2i}\hfill \\ \hfill ⋮\hfill & \hfill \hfill & \hfill ⋮\hfill & \hfill \hfill & \hfill ⋮\hfill \\ \hfill {y}_{Mi}\hfill & \hfill =\hfill & \hfill {\mathbf{x}}_{Mi}^{\prime }{\beta }_{M}\hfill & \hfill +\hfill & \hfill {𝜀}_{Mi}\hfill \\ \hfill \hfill \end{array}$

• the model consists of $M$ equations
• the model is estimated using $N$ observations ($i=1,2,\dots ,N$)
• ${y}_{mi}$ is the value of equation $m$’s dependent variable for observation $i$
• ${\mathbf{x}}_{mi}$ is a ${K}_{m}\phantom{\rule{0.3em}{0ex}}×\phantom{\rule{0.3em}{0ex}}1$ vector that stores the values of the ${K}_{m}$ independent variables for observation $i$, as they appear in equation $m$
• the same independent variable can appear in multiple equations, associated with diﬀerent coeﬃcients
• ${\beta }_{m}$ is a ${K}_{m}\phantom{\rule{0.3em}{0ex}}×\phantom{\rule{0.3em}{0ex}}1$ vector of parameters associated with equation $m$’s independent variables
• in total, there are $K\phantom{\rule{0.3em}{0ex}}=\phantom{\rule{0.3em}{0ex}}\sum _{m=1}^{M}{K}_{m}$ $\beta$ slope parameters to be estimated
• the $M$ error terms jointly follow a multivariate Normal distribution with mean $\mathbf{0}$ and precision matrix $\Omega$

An equivalent and more compact representation of the model is:

${\mathbf{y}}_{i}={\mathbf{X}}_{i}\beta +{𝜀}_{i},\phantom{\rule{2em}{0ex}}\phantom{\rule{1em}{0ex}}{𝜀}_{i}\sim \mathrm{N}\left(\mathbf{0},{\Omega }^{-1}\right)$

where:

$\underset{M\phantom{\rule{0.3em}{0ex}}×\phantom{\rule{0.3em}{0ex}}1}{{\mathbf{y}}_{i}}=\left[\begin{array}{c}\hfill {y}_{1i}\hfill \\ \hfill {y}_{2i}\hfill \\ \hfill ⋮\hfill \\ \hfill {y}_{Mi}\hfill \end{array}\right],\phantom{\rule{2em}{0ex}}\underset{M\phantom{\rule{0.3em}{0ex}}×\phantom{\rule{0.3em}{0ex}}K}{{\mathbf{X}}_{i}}=\left[\begin{array}{cccc}\hfill {\mathbf{x}}_{1i}^{\prime }\hfill & \hfill \mathbf{0}\hfill & \hfill \dots \hfill & \hfill \mathbf{0}\hfill \\ \hfill \mathbf{0}\hfill & \hfill {\mathbf{x}}_{2i}^{\prime }\hfill & \hfill \dots \hfill & \hfill \mathbf{0}\hfill \\ \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill \ddots \hfill & \hfill ⋮\hfill \\ \hfill \mathbf{0}\hfill & \hfill \mathbf{0}\hfill & \hfill \dots \hfill & \hfill {\mathbf{x}}_{Mi}^{\prime }\hfill \end{array}\right],\phantom{\rule{2em}{0ex}}\underset{K\phantom{\rule{0.3em}{0ex}}×\phantom{\rule{0.3em}{0ex}}1}{\beta }=\left[\begin{array}{c}\hfill {\beta }_{1}\hfill \\ \hfill {\beta }_{2}\hfill \\ \hfill ⋮\hfill \\ \hfill {\beta }_{M}\hfill \end{array}\right],\phantom{\rule{2em}{0ex}}\underset{M\phantom{\rule{0.3em}{0ex}}×\phantom{\rule{0.3em}{0ex}}1}{{𝜀}_{i}}=\left[\begin{array}{c}\hfill {𝜀}_{1i}\hfill \\ \hfill {𝜀}_{2i}\hfill \\ \hfill ⋮\hfill \\ \hfill {𝜀}_{Mi}\hfill \end{array}\right]\phantom{\rule{2.22198pt}{0ex}}\phantom{\rule{2.22198pt}{0ex}}$

Priors

 Parameter Probability density function Default hyperparameters $\beta$ $p\left(\beta \right)=\frac{|\mathbf{P}{|}^{1∕2}}{{\left(2\pi \right)}^{K∕2}}exp\left\{-\frac{1}{2}{\left(\beta -\mathbf{m}\right)}^{\prime }\mathbf{P}\left(\beta -\mathbf{m}\right)\right\}$ $\mathbf{m}={\mathbf{0}}_{K}$, $\mathbf{P}=0.001\cdot {\mathbf{I}}_{K}$ $\Omega$ $p\left(\Omega \right)=\frac{|\Omega {|}^{\frac{n-M-1}{2}}|{\mathbf{V}}^{-1}{|}^{n∕2}}{{2}^{nM∕2}{\Gamma }_{M}\left(\frac{n}{2}\right)}exp\left\{-\frac{1}{2}tr\left({\mathbf{V}}^{-1}\Omega \right)\right\}$ $n={M}^{2}$, $\mathbf{V}=\frac{100}{M}\cdot {\mathbf{I}}_{M}$

Syntax

$\left[$<model name> = $\right]$ sur( {
y1 ~ x11 x12  x1K${}_{1}$,
y2 ~ x21 x22  x2K${}_{2}$,
,
yM ~ xM1 xM2  xMK${}_{M}$ }
$\left[$, <options> $\right]$
);

where:

• y1, y2, …, yM are the dependent variable names, as they appear in the dataset used for estimation
• xm1 xm2 $\dots$xmK${}_{m}$ is a list of the ${K}_{m}$ independent variable names for equation $m=1,2,\dots ,M$, as they appear in the dataset used for estimation; when a constant term is to be included in an equation, this must be requested explicitly; $M$ such lists must be provided

The optional arguments for the Seemingly Unrelated Regressions model are:1

 Gibbs parameters "chains" number of chains to run in parallel (positive integer); the default value is 1 "burnin" number of burn-in draws per chain (positive integer); the default value is 10000 "draws" number of retained draws per chain (positive integer); the default value is 20000 "thin" value of the thinning parameter (positive integer); the default value is 1 "seed" value of the seed for the random-number generator (positive integer); the default value is 42 Hyperparameters "m" mean vector of the prior for $\beta$ ($K\phantom{\rule{0.3em}{0ex}}×\phantom{\rule{0.3em}{0ex}}1$ vector); the default value is ${\mathbf{0}}_{K}$ "P" precision matrix of the prior for $\beta$ ($K\phantom{\rule{0.3em}{0ex}}×\phantom{\rule{0.3em}{0ex}}K$ symmetric and positive-deﬁnite matrix); the default value is $0.001\phantom{\rule{0.3em}{0ex}}\cdot \phantom{\rule{0.3em}{0ex}}{\mathbf{I}}_{K}$ "V" scale matrix of the prior for $\Omega$ ($M\phantom{\rule{0.3em}{0ex}}×\phantom{\rule{0.3em}{0ex}}M$ symmetric and positive-deﬁnite matrix); the default value is $\frac{100}{M}\phantom{\rule{0.3em}{0ex}}\cdot \phantom{\rule{0.3em}{0ex}}{\mathbf{I}}_{M}$ "n" degrees-of-freedom parameter of the prior for $\Omega$ (real number greater than or equal to $M$); the default value is ${M}^{2}$ Dataset and log-marginal likelihood "dataset" the id value of the dataset that will be used for estimation; the default value is the ﬁrst dataset in memory (in alphabetical order) "logML_CJ" boolean indicating whether the Chib (1995)/Chib & Jeliazkov (2001) approximation to the log-marginal likelihood should be calculated (true$|$false); the default value is false

Reported Parameters

 $\beta$ variable_name vector of parameters associated with the independent variables; these are broken into groups according to the equation in which the independent variables appear

Stored values and post-estimation analysis
If a left-hand-side id value is provided when a SUR model is created, then the following results are saved in the model item and are accessible via the ‘.’ operator:

 Samples a matrix containing the draws from the posterior of $\beta$ (across all equations, starting from the ﬁrst equation) and the unique elements of $\Omega$ ym$xm1,$\dots$, ym$xmK${}_{m}$ vectors containing the draws from the posterior of the parameters associated with variables xm1,$\dots$,xmK${}_{m}$, for $m=1,2,\dots ,M$ (the names of these vectors are the names of the variables that were included in the right-hand side of equation $m$, prepended by ym$, where ym is the name of the dependent variable in equation $m$; this is done so that the samples on the parameters associated with a variable that appears in more than one equations can be distinguished) Omega_i_j vectors containing the draws from the posterior of the unique elements of $\Omega$; because $\Omega$ is symmetric, only $\frac{\left(M-1\right)M}{2}\phantom{\rule{0.3em}{0ex}}+\phantom{\rule{0.3em}{0ex}}M$ of its elements are stored (instead of all ${M}^{2}$ elements); i and j index the row and column of $\Omega$, respectively, at which the corresponding element is located Omega $M\phantom{\rule{0.3em}{0ex}}×\phantom{\rule{0.3em}{0ex}}M$ matrix that stores the posterior mean of $\Omega$ logML the Lewis & Raftery (1997) approximation of the log-marginal likelihood logML_CJ the Chib (1995)/Chib & Jeliazkov (2001) approximation to the log-marginal likelihood; this is available only if the model was estimated with the "logML_CJ"=true option nchains the number of chains that were used to estimate the model nburnin the number of burn-in draws per chain that were used when estimating the model ndraws the total number of retained draws from the posterior ($=$chains $\cdot$ draws) nthin value of the thinning parameter that was used when estimating the model nseed value of the seed for the random-number generator that was used when estimating the model Additionally, the following functions are available for post-estimation analysis (see section B.14): • diagnostics() • test() • pmp() Examples Example 1 myData = import("$BayESHOME/Datasets/dataset5.csv");
myData.constant = ones(rows(myData), 1);

model1 = sur( {
y1 ~ constant x1 x2 x3 x4 x5 x6 x7 x8 x9 x10,
y2 ~ constant x1 x2 x3,
y3 ~ constant x1 x2 x3
} );

Example 2

myData = import("$BayESHOME/Datasets/dataset5.csv"); myData.constant = ones(rows(myData), 1); model1 = sur( { y1 ~ constant x1 x2 x3 x4 x5 x6 x7 x8 x9 x10, y2 ~ constant x1 x2 x3, y3 ~ constant x1 x2 x3 }, "logML_CJ"=true ); print(mean([model1.y1$x1-model1.y2$constant, model1.y1$x2-model1.y3$constant])); model2 = sur( { y1 ~ constant x1 x2 x3 x4 x5 x6 x7 x8 x9 x10, y2 ~ constant x1 x2 x3, y3 ~ constant x1 x2 x3 }, "constraints" = { y1$x1-y2$constant=0, y1$x5-0.5*y2$x1=0, y1$x6-y2$x2=0, y1$x7-y2$x3=0, y1$x2-y3$constant=0, y1$x6-y3$x1=0, y1$x8-0.5*y3$x2=0, y1$x9-y3$x3=0, }, "Xi" = 1e7*eye(8,8), "logML_CJ"=true ); print(mean([model2.y1$x1-model2.y2$constant, model2.y1$x2-model2.y3\$constant]));

pmp( {model1, model2} );
pmp( {model2, model2}, "logML_CJ" = true);

1Optional arguments are always given in option-value pairs (eg. "chains"=3).