SemiParBIVProbit(SemiParBIVProbit)
SemiParBIVProbit()所属R语言包：SemiParBIVProbit
Semiparametric Bivariate Probit Modelling
半参数的二元概率模型
译者：生物统计家园网 机器人LoveR
描述----------Description----------
SemiParBIVProbit can be used to fit bivariate probit models where the linear predictors can be flexibly specified using parametric and regression spline components, and nonparametric random effects. During the model fitting process, the possible presence of correlated error equations, endogeneity or sample selection is accounted for. Regression spline bases are extracted from the package mgcv. Multi-dimensional smooths are available via the use of penalized thin plate regression splines (isotropic). The current implementation does not support scale invariant tensor product smooths.
SemiParBIVProbit可以用来适合二元probit模型的线性预测，可以灵活使用参数和回归样条组件，和非参数随机效应指定。在模型拟合过程中，可能存在的相关误差方程，内生性或样本选择。回归样条基中提取的包mgcv。多维平滑可通过处罚薄板样条回归（各向同性）使用。目前的实现并不支持规模不变的张量积平滑。
用法----------Usage----------
SemiParBIVProbit(formula.eq1, formula.eq2, data=list(),
selection=FALSE, H=FALSE, iterlimSP=50, pr.tol=1e-6,
gamma=1, aut.sp=TRUE, fp=FALSE, start.v=NULL,
rinit=1, rmax=100, fterm=sqrt(.Machine$double.eps),
mterm=sqrt(.Machine$double.eps),
control=list(maxit=50,tol=1e-6,step.half=25,
rank.tol=sqrt(.Machine$double.eps)),
npRE=FALSE, K=3, id=NULL, e.npRE=TRUE, e.npREsp=TRUE)
参数----------Arguments----------
参数：formula.eq1
A GAM formula for equation 1. s terms are used to specify smooth smooth functions of predictors. SemiParBIVProbit supports the use shrinkage smoothers for variable selection purposes. See the examples below and the documentation of mgcv for further details on GAM formula specifications. Note that if selection=TRUE then the formula MUST refer to the selection equation.
一个GAM公式为式（1）。 s条款用于指定平滑顺畅的预测功能。 SemiParBIVProbit支持使用收缩平滑的变量选择的目的。请参阅下面的例子和文档mgcv GAM公式规格的进一步详情。需要注意的是，如果selection=TRUE那么该公式必须参考的选择公式。
参数：formula.eq2
A GAM formula for equation 2.
一个GAM公式，式（2）。
参数：data
An optional data frame, list or environment containing the variables in the model. If not found in data, the variables are taken from environment(formula), typically the environment from which SemiParBIVProbit is called.
一个可选的数据框，列表或包含在模型中的变量的环境。如果没有找到data，变量environment(formula)，通常是SemiParBIVProbit被称为环境。
参数：selection
If TRUE, then the numerical routine for bivariate probit modelling in the presence of sample selection is employed.
如果TRUE，那么数值的常规二元概率模型中存在的样本选择采用。
参数：H
If TRUE, then the hessian (rather than the Fisher) information matrix is used for fully parametric bivariate probit modelling.
如果TRUE，然后麻（而不是费舍尔）信息矩阵是用于完全参数化的二元概率模型。
参数：iterlimSP
A positive integer specifying the maximum number of loops to be performed before the smoothing parameter estimation step is terminated.
一个正整数，指定之前必须执行的平滑参数估计步骤中的最大数量的循环终止。
参数：pr.tol
Tolerance to use in judging convergence of the algorithm when automatic smoothing parameter selection is used.
公差使用，在使用时自动平滑参数的选择，判断算法的收敛性。
参数：gamma
It is is an inflation factor for the model degrees of freedom in the UBRE score. Smoother models can be obtained setting this parameter to a value greater than 1. Typically gamma=1.4 achieves this.
它是通货膨胀因素的模型程度的自由的UBRE得分。将该参数设置为一个大于1的值，可以得到更流畅的模型。通常gamma=1.4达致这个目标。
参数：aut.sp
If TRUE, then automatic multiple smoothing parameter selection is carried out. If FALSE, then smoothing parameters are set to the values obtained from the univariate fits.
如果TRUE，然后自动多平滑参数的选择进行。如果FALSE，然后平滑化参数被设置为从单变量拟合得到的值。
参数：fp
If TRUE, then a fully parametric model with regression splines if fitted. This only makes sense if used jointly with aut.sp=FALSE. See the example below.
如果TRUE，如果安装了全参数化模型回归样条曲线。这仅是有道理的，如果共同使用aut.sp=FALSE。请看下面的例子。
参数：start.v
Although not recommended, starting values for the parameters of the two equations and correlation coefficient can be provided here.
虽然不建议，可以提供这两个方程的参数和相关系数的初始值，用于这里。
参数：rinit
Starting trust region radius. The trust region radius is adjusted as the algorithm proceeds. See the documentation of trust for further details.
开始信赖域半径。信赖域半径的算法进行调整。有关进一步详情，请参阅文档的trust。
参数：rmax
Maximum allowed trust region radius. This may be set very large. If set small, the algorithm traces a steepest descent path.
允许的最大信赖域半径。这可以被设置非常大的。如果设置小，算法的痕迹，一个最速下降路径。
参数：fterm
Positive scalar giving the tolerance at which the difference in objective function values in a step is considered close enough to zero to terminate the algorithm.
正标量给予的差异被认为是在一个步骤中的目标函数的值在足够接近为零，以终止该算法的公差。
参数：mterm
Positive scalar giving the tolerance at which the two-term Taylor-series approximation to the difference in objective function values in a step is considered close enough to zero to terminate the algorithm.
正标量的公差在两个术语的泰勒级数近似给出在一个步骤中的目标函数的值的差在被认为是足够接近为零，以终止该算法。
参数：control
It is a list containing iteration control constants with the following elements: maxit: maximum number of iterations of the magic algorithm; tol: tolerance to use in judging convergence; step.half: if a trial step fails then the method tries halving it up to a maximum of step.half times; rank.tol: constant used to test for numerical rank deficiency of the problem. See the documentation of magic in mgcv for further details.
这是一个列表，其中包含迭代控制参数包含下列元素：maxit：magic算法的迭代的最大数量，“tol：宽容使用在判断收敛，”step.half如果一个试验步骤失败，则该方法试图减半它最大的step.half倍rank.tol：恒用于测试数值秩不足的问题。请参阅文档magicmgcv进一步的细节。
参数：npRE
If TRUE, then the numerical routine for bivariate probit modelling with nonparametric random effects is employed.
如果TRUE，然后采用常规的二元概率与非参数随机效应模型的数值。
参数：K
If npRE=TRUE, then K represents the number of bivariate mass points and must be supplied. Default is 3.
如果npRE=TRUE，那么K代表二元质点的数量，必须提供。默认值是3。
参数：id
If npRE=TRUE, then id represents the individual random effect identifier and must be supplied.
如果npRE=TRUE，那么id代表的个人随机效应标识符和必须提供。
参数：e.npRE
If npRE=TRUE, then e.npRE indicates whether some extra iterations in the starting value procedure have to be carried out.
如果npRE=TRUE，然后e.npRE表示是否必须进行一些额外的迭代初始值的过程。
参数：e.npREsp
If npRE=TRUE, then e.npREsp indicates whether some extra iterations in the smoothing parameter selection loop have to be carried out.
如果npRE=TRUE，然后e.npREsp表示是否必须进行一些额外的平滑参数选择循环迭代。
Details
详细信息----------Details----------
The bivariate probit model has a probit link for each of the two equations, and models the association between the responses by the correlation parameter ρ of a standardised bivariate normal distribution. In a semiparametric bivariate probit model the linear predictors are flexibly specified using parametric components and smooth functions of covariates. Replacing the smooth components with their regression spline expressions yields a fully parametric bivariate probit model. In principle, classic maximum likelihood estimation can be employed. However, to avoid overfitting, penalized likelihood maximization has to be employed instead. Here the use of penalty matrices allows for the suppression of that part of smooth term complexity which has no support from the data. The trade-off between smoothness and fitness is controlled by smoothing parameters associated with the penalty matrices. Smoothing parameters are chosen to minimize the approximate Un-Biased Risk Estimator (UBRE).
二元概率模型有两个方程，概率链接和模型的相关参数ρ一个标准的二元正态分布的反应之间的关联。在一个半参数二元probit模型的线性预测灵活使用的协变量的参数化元件和平滑的功能。更换的顺利回归样条曲线的表达得到了全参数化的二元概率模型。原则上，经典的最大似然估计可以采用。然而，为了避免过拟合，处罚的可能性最大化，而不是被雇用。这里使用的刑罚矩阵可以抑制的那部分光滑术语复杂的数据不支持从。贸易之间的平滑度和健身的处罚矩阵相关的平滑参数控制。平滑参数的选择，以尽量减少近似无偏风险估计（UBRE）。
The optimization problem is solved by Fisher scoring's method. Automatic smoothing parameter selection is integrated using a performance-oriented iteration approach (Gu, 1992; Wood, 2004). At each iteration, (i) the penalized weighted least squares problem is solved, then (ii) the smoothing parameters of that problem estimated by approximate UBRE. Steps (i) and (ii) are iterated until convergence. Details of the underlying fitting methods are given in Marra and Radice (2011).
由Fisher评分的方法解决该优化问题。自动平滑参数的选择是综合使用以业绩为导向的迭代方法（顾，1992年，木，2004年）。在每次迭代中，（i）本惩罚加权最小二乘问题就解决了，然后（ii）平滑参数这个问题的的近似UBRE估计的。步骤（i）及（ii）进行迭代直到收敛为止。详细的基本拟合方法是给定的在马拉和雷迪斯（2011）。
The use of nonparametric random effects is allowed for via the option npRE and is illustrated in Example 4 below. See Marra and Papageorgiou (submitted) for full details.
允许使用非参数随机效应为通过选项“npRE下面的示例4所示。马拉和Papageorgiou（提交）的全部细节。
值----------Value----------
The function returns an object of class SemiParBIVProbit as described in SemiParBIVProbitObject.
该函数返回一个类的对象SemiParBIVProbit所描述的SemiParBIVProbitObject。
警告----------WARNINGS----------
Any automatic smoothing parameter selection procedure is not likely to work well when the data have low information content. In bivariate probit models, this issue is especially relevant when ρ is high and the number of observations low. Here, convergence failure is typically associated with an infinite cycling between the two steps detailed above. If this occurs, as some practical solutions, one might either (i) lower the total number of parameters to estimate by reducing the dimension of the regression spline bases, (ii) set the smoothing parameters to the values obtained from the univariate fits (aut.sp=FALSE), or (iii) set the smoothing parameters to the values obtained from the non-converged algorithm. The default option is (iii).
任何自动平滑参数选择程序是不可能很好地工作的数据时，具有较低的信息内容。这个问题是在二元probit模型，特别是有关当ρ是高和低的观测数。在这里，收敛失败通常是与上文详述的在两个步骤之间的无限循环。如果发生这种情况，因为一些实际的解决方案，其中一个或可能（ⅰ）降低的总数量的参数来估计通过减少回归样条基的尺寸，（ⅱ）设置的平滑化参数的值从单变量拟合得到（ aut.sp=FALSE），或（iii）设置从非融合算法得到的值的平滑化参数。默认的选项是（III）。
Fully parametric modelling is allowed for. However, it is not possible to specify one linear predictor as a function of parametric and smooth components, and the other as a function of parametric terms only. If continuous covariates are available, then we should let the data determine which effects are linear or non-linear and for which equations.
全参数化建模是允许的。然而，这是不可能指定一个线性预测参数和平滑分量作为一个功能，和其他的作为的函数的参数术语只。连续协变量的情况下，那么我们就应该让数据确定它的效果是线性或非线性的方程。
（作者）----------Author(s)----------
Maintainer: Giampiero Marra <a href="mailto:giampiero@stats.ucl.ac.uk">giampiero@stats.ucl.ac.uk</a>
参考文献----------References----------
Journal of Statistics, 39(2), 259-279.
参见----------See Also----------
AT, InfCr, LM.bpm, plot.SemiParBIVProbit, SemiParBIVProbit-package, SemiParBIVProbitObject, summary.SemiParBIVProbit, residuals.SemiParBIVProbit
AT，InfCr，LM.bpm，plot.SemiParBIVProbit，SemiParBIVProbit-package，SemiParBIVProbitObject，summary.SemiParBIVProbit，residuals.SemiParBIVProbit
实例----------Examples----------
library(SemiParBIVProbit)
############[＃＃＃＃＃＃＃＃＃＃＃]
## EXAMPLE 1[＃示例1]
############[＃＃＃＃＃＃＃＃＃＃＃]
## Generate data[＃生成数据]
## Correlation between the two equations 0.5 - Sample size 400 [两个方程的相关性0.5 - 样品大小400]
set.seed(0)
n <- 400
Sigma <- matrix(c(1,0.5,0.5,1),2,2)
u <- rmvnorm(n,rep(0,2),Sigma)
x1 <- round(runif(n))
x2 <- runif(n)
x3 <- runif(n)
f1 <- function(x) (cos(pi*2*x)) + sin(pi*x)
f2 <- function(x) (x+exp(-30*(x-0.5)^2))
y1 <- ifelse(-1.55 + 2*x1 + f1(x2) + u[,1] > 0, 1, 0)
y2 <- ifelse(-0.25 - 1.25*x1 + f2(x2) + u[,2] > 0, 1, 0)
dataSim <- data.frame(y1,y2,x1,x2,x3)
#[]
#[]
## CLASSIC BIVARIATE PROBIT[＃经典二元PROBIT]
out <- SemiParBIVProbit(y1 ~ x1 + x2 + x3,
y2 ~ x1 + x2 + x3,
data=dataSim)
summary(out)
InfCr(out)
InfCr(out,cr="BIC")
## SEMIPARAMETRIC BIVARIATE PROBIT[＃半参数二元PROBIT]
## "cr" cubic regression spline basis - "cs" shrinkage version of "cr"[＃“CR”三次回归样条曲线的基础 - “CS”缩水版“CR”]
## "tp" thin plate regression spline basis - "ts" shrinkage version of "tp"[＃“TP”薄板回归样条曲线的基础 - “TS”缩水版“TP”]
## for smooths of one variable, "cr/cs" and "tp/ts" achieve similar results [＃平滑的一个变量，“CR / CS”和“TP / TS”达到类似的效果]
## k is the basis dimension - default is 10[＃k为基础的维度 - 默认值是10]
## m is the order of the penalty for the specific term - default is 2[＃m是为了处罚的具体期限 - 默认为2]
out <- SemiParBIVProbit(y1 ~ x1 + s(x2,bs="cr",k=10,m=2) + s(x3,bs="cr",k=10),
y2 ~ x1 + s(x2,bs="cr",k=10) + s(x3,bs="cr",k=10),
data=dataSim)
summary(out)
InfCr(out)
## estimated smooth function plots - red lines are true curves[＃估计光滑函数图 - 红色的线是真正的曲线]
x2 <- sort(x2)
f1.x2 <- f1(x2)[order(x2)]-mean(f1(x2))
f2.x2 <- f2(x2)[order(x2)]-mean(f2(x2))
f3.x3 <- rep(0,length(x3))
par(mfrow=c(2,2),mar=c(4.5,4.5,2,2))
plot(out, eq=1, select=1); lines(x2, f1.x2, col="red")
plot(out, eq=1, select=2); lines(x3, f3.x3, col="red")
plot(out, eq=2, select=1); lines(x2, f2.x2, col="red")
plot(out, eq=2, select=2); lines(x3, f3.x3, col="red")
#[]
## same plots but CIs 'with intercept' [＃一样的图，但证明书“拦截”]
par(mfrow=c(2,2),mar=c(4.5,4.5,2,2))
plot(out, eq=1, select=1, seWithMean=TRUE); lines(x2, f1.x2, col="red")
plot(out, eq=1, select=2, seWithMean=TRUE); lines(x3, f3.x3, col="red")
plot(out, eq=2, select=1, seWithMean=TRUE); lines(x2, f2.x2, col="red")
plot(out, eq=2, select=2, seWithMean=TRUE); lines(x3, f3.x3, col="red")
## p-values suggest to drop x3 from both equations, with a stronger [＃p值建议下降X3的两个方程，具备了更强]
## evidence for eq. 2. This can be also achieved using shrinkage smoothers[＃的证据当量。 2。这也可通过使用收缩平滑]
outSS <- SemiParBIVProbit(y1 ~ x1 + s(x2,bs="cs") + s(x3,bs="cs"),
y2 ~ x1 + s(x2,bs="cs") + s(x3,bs="cs"),
data=dataSim)
par(mfrow=c(2,2),mar=c(4.5,4.5,2,2))
plot(outSS, eq=1, select=1)
plot(outSS, eq=1, select=2, ylim=c(-0.1,0.1))
plot(outSS, eq=2, select=1)
plot(outSS, eq=2, select=2, ylim=c(-0.1,0.1))
#[]
#[]
############[＃＃＃＃＃＃＃＃＃＃＃]
## EXAMPLE 2[＃示例2]
############[＃＃＃＃＃＃＃＃＃＃＃]
## Generate data with one endogenous variable and exclusion restriction[＃生成数据的一个内生变量和约束]
set.seed(0)
n <- 400
Sigma <- matrix(c(1,0.5,0.5,1),2,2)
u <- rmvnorm(n,rep(0,2),Sigma)
x1 <- round(runif(n))
x2 <- runif(n)
f1 <- function(x) (cos(pi*2*x)) + sin(pi*x)
f2 <- function(x) (x+exp(-30*(x-0.5)^2))
y1 <- ifelse(-1.55 + 2*x1 + f1(x2) + u[,1] > 0, 1, 0)
y2 <- ifelse(-0.25 - 1.25*y1 + f2(x2) + u[,2] > 0, 1, 0)
dataSim <- data.frame(y1,y2,x1,x2)
#[]
## Testing the hypothesis of absence of endogeneity... [＃测试的假设情况下的内生性...]
LM.bpm(y1 ~ x1 + s(x2),y2 ~ y1 + s(x2),dataSim)
# p-value suggests presence of endogeneity, hence fit a bivariate model[p值显示内生性的存在，因此适合一个双变量模型]
## CLASSIC RECURSIVE BIVARIATE PROBIT[＃经典的递归二元PROBIT]
# Algorithm based on Fisher Scoring[一种基于Fisher评分]
out <- SemiParBIVProbit(y1 ~ x1 + x2,
y2 ~ y1 + x2,
data=dataSim)
summary(out)
# Algorithm based on Newton's method[一种基于牛顿的方法]
out <- SemiParBIVProbit(y1 ~ x1 + x2,
y2 ~ y1 + x2,
data=dataSim,H=TRUE)
summary(out)
## SEMIPARAMETRIC RECURSIVE BIVARIATE PROBIT[＃半参数递归二元PROBIT]
out <- SemiParBIVProbit(y1 ~ x1 + s(x2),
y2 ~ y1 + s(x2),
data=dataSim)
summary(out)
#[]
## average treatment effect on the treated (ATT) with CIs[＃平均处理后的治疗效果与独联体（ATT）]
AT(out,eq=2,nm.bin="y1",E=FALSE)
#[]
## recursive bivariate probit modelling with unpenalized splines [＃递归二元概率模型与unpenalized样条]
## can be achieved as follows[＃就可以实现如下]
outFP <- SemiParBIVProbit(y1 ~ x1 + s(x2,bs="cr",k=5),
y2 ~ y1 + s(x2,bs="cr",k=6), aut.sp=FALSE, fp=TRUE,
data=dataSim)
summary(outFP)
## example using a 2D smooth [例如，使用一个2D光滑]
f.biv <- function(x,z) -0.7*exp(((3*x+3) + 0.7*(3*z-3)^2)/5)
y1 <- ifelse(-1.55 + 2*x1 + f1(x2) + u[,1] > 0, 1, 0)
y2 <- ifelse(3.5 - 1.25*y1 + f.biv(x2,x3) + u[,2] > 0, 1, 0)
dataSim <- data.frame(y1,y2,x1,x2,x3)
outSb <- SemiParBIVProbit(y1 ~ x1 + s(x2,bs="cr"),
y2 ~ y1 + s(x2,x3,bs="tp"),
data=dataSim)
summary(outSb)
par(mfrow=c(1,2),mar=c(4.5,4.5,2,2))
x2 <- x3 <- seq(0, 1, length=40)
z <- outer(x2, x3, f.biv)
persp(x2, x3, z, theta = 30, phi = 25,
xlab="x2", ylab="x3", zlab="f.biv")
plot(outSb, eq=2, select=1,theta = 30, phi = 25)
#[]
#[]
############[＃＃＃＃＃＃＃＃＃＃＃]
## EXAMPLE 3[＃示例3]
############[＃＃＃＃＃＃＃＃＃＃＃]
## Generate data with sample selection mechanism and exclusion restriction[＃生成数据与样本选择机制和排除限制]
set.seed(0)
n <- 2000
Sigma <- matrix(c(1,0.5,0.5,1),2,2)
u <- rmvnorm(n,rep(0,2),Sigma)
SigmaC <- matrix( c(1,0.5,0.5,0.5,1,0.5,0.5,0.5,1), 3 , 3)
cov <- rmvnorm(n,rep(0,3),SigmaC, method="svd")
cov <- pnorm(cov)
bi <- round(cov[,1]); x1 <- cov[,2]; x2 <- cov[,3]
f11 <- function(x) -0.7*(4*x + 2.5*x^2 + 0.7*sin(5*x) + cos(7.5*x))
f12 <- function(x) -0.4*( -0.3 - 1.6*x + sin(5*x))
f21 <- function(x) 0.6*(exp(x) + sin(2.9*x))
ys <- 0.58 + 2.5*bi + f11(x1) + f12(x2) + u[, 1] > 0
y <- -0.68 - 1.5*bi + f21(x1) + + u[, 2] > 0
yo <- y*(ys > 0)
dataSim <- data.frame(y,ys,yo,bi,x1,x2)
## Testing the hypothesis of absence of sample selection... [＃测试假设的情况下，样本选择...]
LM.bpm(ys ~ bi + s(x1) + s(x2),yo ~ bi + s(x1),dataSim,selection=TRUE)
# p-value suggests presence of sample selection, hence fit a bivariate model[p值表明存在样本选择，因此适合一个双变量模型]
#[]
## SEMIPARAMETRIC SAMPLE SELECTION BIVARIATE PROBIT[＃半参数的SAMPLE选择二元PROBIT]
## the first equation MUST be the selection equation[＃第的方程必需的选择方程]
out <- SemiParBIVProbit(ys ~ bi + s(x1) + s(x2),
yo ~ bi + s(x1),
data=dataSim, selection=TRUE)
AT(out,eq=2,nm.bin="bi",E=FALSE) ## ATT[＃ATT]
## compare the two summary outputs[＃比较两个摘要输出，]
## the second output produces a summary of the results obtained when only [＃第二个输出产生一个总结时得到的结果只]
## the outcome equation is fitted, i.e. selection bias is not accounted for[＃结果方程拟合，即选择偏倚不占]
summary(out)
summary(out$gam2)
## mean predicted probabilities that 'yo' is equal to 1[＃平均预测“哟”的概率为1]
## the second predicted probability does not account for selection bias[＃第二个预测的概率不考虑选择偏倚]
mean(pnorm(out$eta2[out$eta2!=0]))
mean(pnorm(predict(out$gam2)))
## estimated smooth function plots[＃估计光滑的函数曲线]
## the red line is the true curve[＃红线是真正曲线]
## the blue line is the naive curve not accounting for selection bias[＃蓝线是天真的曲线，不占选择偏倚]
x1 <- sort(x1)
f21.x1 <- f21(x1)[order(x1)]-mean(f21(x1))
plot(out, eq=2, select=1,ylim=c(-1.2,1)); lines(x1, f21.x1, col="red")
par(new=TRUE)
plot(out$gam2, select=1, se=FALSE, col="blue",ylim=c(-1.2,1),ylab="",rug=FALSE)
#[]
#[]
############[＃＃＃＃＃＃＃＃＃＃＃]
## EXAMPLE 4[＃示例4]
############[＃＃＃＃＃＃＃＃＃＃＃]
## Generate data with one endogenous variable, exclusion restriction[＃生成数据的一个内生变量，排除限制。]
## and subject specific features [并受特定的功能]
set.seed(0)
n.cl <- 100
n.ob.cl <- 13
sre <- round( runif(n.cl,n.ob.cl,n.ob.cl+2) )
t.n <- sum(sre)
id <- rep(seq(1,n.cl),sre)
sigma1 <- 1
sigma2 <- 2
rho <- rho.re <- 0.5
SigmaRE <- matrix(c(sigma1^2, sigma1*sigma2*rho.re,
sigma1*sigma2*rho.re,sigma2^2 ),2,2)
u1 <- rmvnorm(n.cl,rep(0,2),SigmaRE,method="svd")
u <- cbind(rep(u1[,1],sre),rep(u1[,2],sre))
Sigma <- matrix(c(1,rho,rho,1),2,2)
eps <- rmvnorm(t.n,rep(0,2),Sigma)
SigmaC <- matrix(c(1,0.5,0.5,
0.5,1,0.5,
0.5,0.5,1),3,3)
cov <- rmvnorm(t.n,rep(0,3),SigmaC, method="svd")
x1 <- round(pnorm(cov[,1]))
x2 <- pnorm(cov[,2])
x3 <- cov[,3]
f1 <- function(x) 0.5*(exp(x) + sin(2.9*x))
f2 <- function(x) 1.8*(0.25*exp(x) - x^3)
f1.x2 <- f1(x2) - mean(f1(x2))
f2.x2 <- f2(x2) - mean(f2(x2))
y1 <- ifelse(u[,1] + -0.75*x1 + f1.x2 + x3 + eps[,1] > 0, 1, 0)
y2 <- ifelse(u[,2] + -1.5*y1 + 1*x1 + f2.x2 + eps[,2] > 0, 1, 0)
dataSim <- data.frame(id,y1,y2,x1,x2,x3)
out <- SemiParBIVProbit(y1 ~ x1 + s(x2) + x3,
y2 ~ y1 + x1 + s(x2),
data=dataSim, aut.sp=TRUE,
npRE=TRUE, K=3, id=dataSim$id)
summary(out)
AT(out,eq=2,nm.bin="y1",E=FALSE)
x2s <- sort(x2)
f1.x2 <- f1.x2[order(x2)]
f2.x2 <- f2.x2[order(x2)]
par(mfrow=c(1,2))
plot(out, eq=1, select=1); lines(x2s, f1.x2, col="red")
plot(out, eq=2, select=1); lines(x2s, f2.x2, col="red")
#[]
#[]
转载请注明:出自 生物统计家园网(http://www.biostatistic.net)。
注：
注1：为了方便大家学习，本文档为生物统计家园网机器人LoveR翻译而成，仅供个人R语言学习参考使用，生物统计家园保留版权。
注2：由于是机器人自动翻译，难免有不准确之处，使用时仔细对照中、英文内容进行反复理解，可以帮助R语言的学习。
注3：如遇到不准确之处，请在本贴的后面进行回帖，我们会逐渐进行修订。
|
发表于 2012-9-30 00:35:24
