The four-parameter bivariate logistic distribution
has a density that can be written as
$$f(y_1,y_2;l_1,s_1,l_2,s_2) = 2 \frac{\exp[-(y_1-l_1)/s_1 -
(y_2-l_2)/s_2]}{
s_1 s_2 \left( 1 + \exp[-(y_1-l_1)/s_1] + \exp[-(y_2-l_2)/s_2]
\right)^3}$$
where \(s_1>0\) and \(s_2>0\) are the scale parameters,
and \(l_1\) and \(l_2\) are the location parameters.
Each of the two responses are unbounded, i.e.,
\(-\infty<y_j<\infty\).
The mean of \(Y_1\) is \(l_1\) etc.
The fitted values are returned in a 2-column matrix.
The cumulative distribution function is
$$F(y_1,y_2;l_1,s_1,l_2,s_2) =
\left( 1 + \exp[-(y_1-l_1)/s_1] + \exp[-(y_2-l_2)/s_2]
\right)^{-1}$$
The marginal distribution of \(Y_1\) is
$$P(Y_1 \leq y_1) = F(y_1;l_1,s_1) =
\left( 1 + \exp[-(y_1-l_1)/s_1] \right)^{-1} .$$
By default, \(\eta_1=l_1\),
\(\eta_2=\log(s_1)\),
\(\eta_3=l_2\),
\(\eta_4=\log(s_2)\) are the linear/additive predictors.