1、结构方程模型的约束最小二乘解 1)童恒庆 邓艳芳 刘天桢武汉理工大学数学系,武汉,430070摘 要本文首先改进了结构方程模型(SEM)的偏最小二乘(PLS)算法,在单位向量约束下找到了最小二乘迭代初值,提高了收敛速度,说明了收敛性;然后我们将这种模型和算法扩展到多层结构方程模型;最后我们研究多个对象的结构方程模型,借助于作者提出的广义凸约束线性模型将这多个 SEM 统一成一个模型并给出了算法。本文成果都通过了编程计算并收入 DASC 软件。关键词:结构方程模型,偏最小二乘算法,迭代初值,收敛性,约束解MR(2000)主题分类:62H12, 62J05THE LEST SQUARE SOLUT
2、ION IN STRUCTRUAL EQUATION MODELS WITH CONSTRAINTTong Hengqing, Deng Yanfang, Liu TianzhenDepartment of Mathematics, Wuhan University of Technology, Wuhan, 430070AbstractThis paper first improves the partial least square (PLS) algorithm in structural equation models(SEM), and finds a suitable iterat
3、ive initial value with the constraint of unit vector in the meaning of LSE. The algorithm enhances the convergence rate greatly and its convergence is illuminated. Then we build a multi-level structural equation models and find our algorithm is suitable too. At last we build a uniform structural equ
4、ation model with many objections using the generalized linear model with convex constraint and give the algorithm of the model. The results of this paper have been received in software DASC.Keywords: Structural equation models, Partial least square, Iterative initial value, Convergence, Constraint s
5、olution.2000 Mathematics Subject Classification:62H12, 62J051) 国家科技部技术创新基金资助项目(02C26214200218).国家自然科学基金资助项目(35070611)1. 结构方程模型与偏最小二乘结构方程模型(SEM)是应用统计领域近来发展迅速的一个分支,主要包括因果模式、路径图、偏最小二乘(PLS)算法等。SEM 广泛应用于心理学、社会学等领域,尤其是顾客满意指数(CSI)分析模型(Claes Fornel、Michael D. Johnson, et al. 1996)。由于 ISO9000 系列标准和卓越绩效国家标准要求
6、顾客满意指数分析,SEM 的计算就显得非常重要,许多国际知名的软件公司不断推出升级版的 SEM 应用软件(如 LISREL, AOMS, EQS,SEPATH, MPLUS和 SAS 中 CALIS 模块等)。关于 SEM 的算法已有许多进展,包括采用 ML 算法、EM 算法(Liang,et al,2004)等,实际应用中 PLS 算法仍然是主要的和有效的(Michel et al, 2005) 。但是过去都是采用任给的迭代初值,因此不能确保收敛或者收敛速度太慢。关于 SEM 的研究已经扩展到非线性模型(Lee Sik-Yum Lee, Lu Bin, 2003)和多层模型( Lee Sik
7、-Yum; Zhu Hong-Tu, 2002),但是还没有文献研究多对象的 SEM 并将其统一成一个模型。关于 SEM 研究的最新进展已经有了很好的综述报告(Andrew J,2005) 。SEM 包括两个方程组,一个是结构变量之间的关系方程组,称为结构方程组;一个是结构变量与观测变量之间的关系方程组,称为观测方程组。在典型的中国顾客满意指数模型里,结构方程组包含 6 个结构变量(隐含变量) 、 与 11 个关系(自变量作用15的关系为 ,因变量之间的作用关系为 ) ,如(1)所示。41ji(1) 543211432155434213542 0000在一般情形,结构变量不一定是 5 个,结构
8、方程系数形式除了要求对角线是 0 且是下三角矩阵以外也可以不同于(1) ,自变量的个数也可以多于一个。我们采用向量与矩阵记法进行一般描述。设因变量有 个,将 排成列向量,记为 ;自变量有 个,将mm,1 k排成列向量,记为 ; 的系数矩阵为 阶方阵记为 , 的系数矩阵为k,1 B阶矩阵记为 ,残差向量为 ,则结构方程组(1 )可以扩展为:m(2)BSEM 的结构变量是隐含的,不能直接观测,而是对应若干个观测变量。设一共有 个M观测变量,对每一个观测变量有 个观测,在顾客满意指数分析中就是有 个顾客的测NN评,这样我们手里数据是一个 矩阵。M结构变量与观测变量之间的作用关系也可以用方程表示出来,
9、按作用的因果路径有两种表示方式。设与自变量 对应的 个观测变量为 ,与因变量 对t)(tS)(,1,tSjxti应的 个观测变量为 ,于是从观测变量到结构变量的观测方程组可以)(iL,1,iLjyi表达为:, (3)(1tSj tjttxk,1, (4)iiLjjiiy)(1 m,反之,从结构变量到观测变量的观测方程可以表达为:, (5)(1)(1)(1txSxtttSttStx k,, (6)(1)(1)(1iLyiyiLiiiLiy m,其中 、 为载荷项。jtji我们将(2) 、 (3) 、 (4)合在一起(或者(2) 、 (5) 、 (6)合在一起)称为结构方程模型,有时也称为路径分析
10、模型。目前国际通行的路径分析模型求解是采用偏最小二乘(PLS)迭代算法。以顾客满意度模型为例,PLS 方法是先在(3) (4)给 任意初值,于是可以根据已知的观测变jit,量 的值,利用公式(3) (4)求出结构变量 的值;结构方程的变量 有了jityx, itit,数值,就可以对方程组(2)求最小二乘解,于是有了系数 的估计值,进而有了变jit,量 的估计值 。结构变量 与观测变量 都有了数值,就可以对方程组it,it,it, jityx,(3) (4)求解,此时的 就有了估计值,而不再是任给的初值。 有了估计jit, jit,值,就可以回到迭代起点,开始新一轮迭代。迭代收敛控制可以选取迭代
11、过程中所有对应元素的差值小于指定的误差精度。上述迭代过程可以表示为: )0(,(jit )4(3)0()(,外 生外 生 it )2( )0(,jit )2( )0()(,内 生内 生 it )4(3)1(jijt文献1给出的计算方法,迭代初值是任给的,一般取(1,0,0, ,0) ,没有给出收敛性的证明。文献2 介绍他们的 PLS 迭代,250 个样本,18 个指标,大约需要迭代 45 分钟,这对于 21 世纪的专用计算机来说显然收敛速度太慢。2. 单位向量约束下结构方程的最小二乘解本文作者发现,在合理的一定约束条件下,结构方程求解可以根据最小二乘原理计算出偏最小二乘迭代初值,数值计算表明迭
12、代过程收敛很快,我们也可以在一定条件下说明迭代的收敛性。先叙述结构方程解的一些基本性质。(1).方程(2) 、 (3) 、 (4)组成的模型的解并不唯一,可以相差一个常数倍,即若是一组解,则 也是一组解,这里 是任一常数。因此我们可以在 为单it, itc, cit,位向量的条件下求解。(2).方程(3) 、 (4)存在当然 0 解,但是方程(5) (6)并不存在当然 0 解。(3)方程(4)和(6)在 的条件下等价,因而(6)的最小二乘解也就1)(1iLjji是(4)的最小二乘解。已往的 PLS 迭代算法是直接利用方程(3) 、 (4)任给初值进行的,没有利用(5) (6) ,而我们推导约束
13、条件下的最小二乘解将利用方程(5)或(6) 。下面我们推导结构方程模型在单位向量约束条件下的最小二乘迭代初值。将(6)记为,作乘积 ,如果取结构变量为单位向量,即iiiy iiiiiy,则有 。这是两个 的矩阵在最小二乘意义下的近似相等,1i ii )(L写详细一些就是:(7) 2)(2)(1)( )(22111)(2)(1)( )(22 1211 iLiLiiiL iiii LiLiiiLiLi iiii Lyyy 注意左边的元素是两个向量相乘得到的数,右边的元素是数与数相乘得到的数。取对角线的元素相等,即得:, (8)jijiy2)(,1i这样我们得到了观测变量与结构变量之间的系数在最小二
14、乘意义下的估计 。)(1,iLi我们再来估计结构变量 。设 ,我们要逐个估计它的分量。观i),(21Niii 察(6)式,对于 的第 个分量,有is, , (9)sijsjiy)(,1iLjNs,1写成矩阵形式就是:, (10)isLiisiLsiy)(1)(1,记向量 ,它是矩阵 的横截向量。还是根据最小二乘原理,我们又,()(1sisisyY i有了 的最小二乘解:si(11)sisiLsisiisi Yy )(11,(, (12)isisYN,这里的 是已经估计出来的值。类似我们可以估计出 与 ,这样我们得到了全部结构i jtt变量的最小二乘估计值。由本段性质(3) ,从式(6)获得的最
15、小二乘估计在一定条件下也是式(4)的最小二乘估计,因此它满足(13)min|)(1iLjjiiy其几何意义是求一个单位球面与一个线性子空间的距离,在向量 之间不存在线性相关jiy的条件下,解是唯一的。要保证这一点,必须要求这个单位球面与这个线性子空间没有公共点。若不然,即(13)左边可以为 0,也即(4)存在精确解或者线性关系 ,代入(6)将)(1iLjjiiy得出向量 之间存在线性相关的结论。jiy在根据最小二乘原则从观测方程组得到结构变量的最小二乘解以后,回到结构方程组,可以得到路径作用系数 和 的唯一解。将方程(2)移项得:B(14)(I对于顾客满意指数模型,其路径关系表达的方程中,矩阵
16、 总是主对角线为 0 的下三角矩B阵,所以 是下三角矩阵,其行列式为 1,故可逆。于是有BI(15)BBI 1这里 ,所以对于给定的 和 ,结构方程的最小二乘解存在。如果矩阵 列向1IB 量线性无关,这个最小二乘解还是唯一的。上述迭代过程可以表示为: ),(jityx 1|,|),6(5 )0()(,外 生外 生 it)0()(,外 生外 生 it )2( )0(jit )2( )()(内 生内 生 it )4(3)1(,jijt )4(3 )1()(,外 生外 生 it方程(15)将联立方程模型转化成了普通回归分析模型,有利于我们说明解的存在与唯一性。但是实际计算我们不能采用(15) ,因为
17、普通回归估计出的 ,并不一定满足B(1) (2)中矩阵 和 对 0 元素的要求,即对路径关系的要求。实际计算还是需要根据B(1) (2) ,采用一般联立方程模型的二阶段最小二乘,以改进估计的相合性。上述算法已经收入了 DASC 软件。我们要说明这里给出的迭代初值是最小二乘意义下最优的。SEM 的最后解要满足两个方程组,其中的结构方程组在开始时其变量与系数都是未知的,不可能计算出迭代初值。我们根据观测方程组计算出来的迭代初值在最小二乘意义下已经满足了这个方程组,当然就是最小二乘意义下最优的。当结构方程组根据这个初值计算出系数时,结构方程组也就满足了。那么为什么还要继续 PLS 迭代呢?那是因为在
18、第一次迭代过程中,还没有体现结构变量之间的相互作用,所以还要继续迭代,直到两个方程组都被满足,特别是结构方程组内部平衡时为止。顺带说明的是,CSI 最后计算公式需要改进与研究:文献1的公式缺乏统计稳健性(数据的最大值最小值作为单项参与计算):, , (16)min()ax(CSIEniiX1)max()ax(niiX1)m()(我们建议 最终计算公式为:I(17)kt ii iLjjLjijtSjtijkt i ijjijijtjjtij yyx101)0(1)()( )(0)()(0I 这相当于计算出所有观测变量到顾客满意度变量的影响系数,分清了主从关系。与顾客满意度直接关联的观测变量较亲,
19、它们的影响系数就是 ;其余变量较疏,它们的影响系ji数是 与 的乘积,即它们要先通过 影响自己的结构变量,再通过 影响顾客满jii0jii0意度;顾客满意度以后的结构变量将不对顾客满意度发生影响。3. 多层结构方程模型及其算法单层结构方程模型其结构变量都是直接与观测变量相联系,实际工作提出需要多层路径分析模型。例如手机质量感知是一个结构变量,派生出若干个质量类别,如外观质量、通话质量、可靠性质量等类别,每个类别再与若干个观测变量相联系。我们可以将多层路径分析模型画图表达,仔细分析其结构,发现可以将派生的结构变量理解为自变量,于是模型的基本框架就能确定下来。在第一节叙述的中国顾客满意指数模型中,
20、假设其 、 、 不直接与观测变量联系,而是分别派生出 2 个、3 个、2 个低123层的结构变量,它们再与观测变量相联系。这派生的低层结构变量就可以理解为自变量,连同原来的自变量一起就有 8 个: ,因变量还是 5 个: 。于是基本结构方8151程如下: 54321872137265413241532154342135421 00000000 (18)自变量对应的观测变量以 表示,因变量对应的观测变量以 表示。自变量对应的itxjiy结构变量与观测变量关系方程为:, (19)(1)(1)(1txSxtttSttStx 8,在一个算例里 分别为 5,4,3,2,5,4,3,2;同时 分别与 4
21、个观测变量相联)(t 54,系,则直接因变量对应的结构变量与观测变量关系方程为:(20)854154421418500yyyyy间接因变量对应的结构变量与低层结构自变量关系方程为:(21) 876543232121321876543200其中 为载荷项, 为误差项。ijji,i在有了多层结构方程模型的正确描述以后,我们就可以作出它的基于单位向量约束下的最小二乘解。虽然路径系数的计算要复杂得多,变量编号查找与属性判定也非常复杂,但是根据前一节的论述,计算是切实可行的。4. 多对象评估的结构方程模型在凸约束下的最小二乘解我们现在考虑更为深入更为实际的结构方程模型。实际工作中的指标测评包括顾客满意指
22、数测评都是同时针对多个对象进行的。我们应该建立统一的测评模型,对每一个对象的汇总计算系数应该是一致的;同时模型系数应该是合理客观而不是事先指定的。我们考虑每一个对象,都满足第一节建立的结构方程模型,都是 个观测变量,对每M一个观测变量都有 个观测,是一个 矩阵。现在假设有 W 个对象,应该有 W 个NMN模型,问题是这 W 个模型的系数必须是一致的。如果将这 W 个模型统一起来,其数据结构应该是 行、 列的数据块。可以看到模型表面上理解并不复杂,困难在于模型的求解计算。如果对每一个模型分别套用第一节的方法,得到的系数彼此不能保证一致。为了求解多对象评估的结构方程模型,我们需要了解多对象评估的广
23、义配方模型,它是本文作者早期的研究成果(童恒庆,1993) 。我们这里简要叙述广义配方模型结构与主要定理。设 个指标是变量,分别以 表示。一张评估表是对某一对象的一次观测,p)()1(,px可取得数据 ,对 个对象各取得 次观测,就得 阵。pjxki,)(Wk,Ni,1X一张评估表是 阵的一行,一个对象的 次观测是 阵的一块。对每个变量的加权系数XNX待定,但需 (即 ) ; (即p,1jj,1,00cp21)。这是一种配方约束。对每个对象必须且只须给出一个评估分 ,c/ Wiy,1,它也是事先未知而待定的,这就是所谓广义。因此评估模型是:(22) 0min|/ ,2cXDyGPpy1这里 ,
24、 , , 将 块pyD1pWI )1,( ,),(1Wy数据块按列分别求平均,得到压缩的数据阵 ,这N )0(XXDNIp1就引入了广义配方模型。定理 1. 令 ,若 ,则在约束 下, DNIPWD1pPrk)( cp),(yQ有唯一解: ,这里 , 2|Xymin, y )0(Xy/1/pAX。如果 各分量非负,则也就是 模型的解。PND)0( GP若定理 1 叙述的解 有分量为负时,要考虑模型 的解的存在性、唯一性,有定理2。定理 2. 若 ,则 P模型有唯一解。若定理 1 中 有分量为负,则pWXrk)( GP模型的解 一定有分量为 0,并且 的零分量是 的分量之一。00要具体给出 模型
25、的解,需要引入凸集间的交互投影。求一点 到闭凸集 之间的0AB最短欧氏距离,若 , 则可以称 为 到 的投影。要求BAdB00),(),(两个闭凸集 之间的最短欧氏距离,可以使用交互投影法:任取 ,求 ,A, 00使 ;对于 ,求 ,使 ;对于 ,),()(00d01 ),(),(10ABdi求 Bi,使 ;对于 ,求 ,使 。),(BAii iAi),(1dii当 时,停止迭代,完成计算。这个迭代过程收敛的意思是: ),(iA lmii,我们有定理 3。),(,*d定理 3. 设 两个闭凸集之一有界,则其交互投影的迭代过程收敛。B根据定理 3,求两个闭凸集之间的距离可以化为累次求一点到闭凸集
26、间的距离。于是求解广义配方模型可以化为累次求解配方模型,求解凸约束广义配方模型可以化为累次求解一般凸约束模型,实际计算表明,收敛过程非常快。现在回到多对象评估的结构方程模型。每一个对象的观测数据为 矩阵,将MN个对象的观测数据按块纵向排列成 矩阵 。我们按如下步骤求解:WMWNX第一步:对 矩阵 ,按照本文第 1 节叙述的模型和第 2 节叙述的解法,求出NX每个观测变量的汇总系数:(3) (4)中的 。取 , , jit,)(1tSjjtckt, 。)(1iLjjiCmi,第二步:注意 ,我们实际上已经将矩阵 按列划分成 块,MiLtSikt 11)()( Xmk每一块对应一个结构变量。对于每
27、一个结构变量对应的数据块,套用多对象评估的广义配方模型,它的约束条件是 或者 ,当然还ttSc)(21 iiLC)(21有 ;这里的 或 就是广义配方模型的指标数 ;这样得到的广义配方模型0j)(tSiLp评估值 是一个 维的向量,我们记作 或 。我们一共得到了 个 维向量yWjtxjiymkW或 。jtxji第三步:将 个 维向量 或 带回到第 1 节叙述的结构方程模型,运用第 2mkjtxjiy节叙述的模型算法,先计算出(5) (6)中的 与 ;再计算出(3) (4)中的中jit,it,的 ,注意这里的 与第一步的并不相同;然后在(2)中计算出路径影响系jit, jit,数 ,这样开始了偏
28、最小二乘迭代并结束。jti第四步:取最后的计算结果 ,它是一个 维向量,它的第 个分量就是第 个对象0iWii获得的最后评估分。模型的系数还是反映了这些变量之间影响作用的大小。本文数值算例这里省略,可以在作者主页上查到:http:/ 。参 考 文 献1 Claes Fornel、Michael D. Johnson, et al. The American Customer Satisfaction Index: Nature, Popurse, and Findings, Journal of Marketing,Vol. No.60, Oct.1996,7-182 国家质检总局质量管理司,
29、清华大学中国企业研究中心,中国顾客满意指数指南,中国标准出版社,2003. 3 S.Y. Sohn, T.H. MoonStructural Equation Model for Predicting Technology:Commercialization Success Index (TCSI)J Technological Forecasting Zhu, Hong-Tu, Muximum likelihood estimation of nonlinear structural equation models, Psychometrika, Vol. 67 No.2, Jun. 200
30、2, 189-210.8 Liang J. and Bentler PM, An EM algorithm for fitting two-level structural equation models, Psychometrika, Vol. 69, 2004, 101-122.9 Andrew J. Tomarken and Niels G. Waller, Structural Equation Modeling: Strengths, Limitations and Misconceptions, Annu. Rev. Clin. Psychol. Vol. 1, 2005, 31-65.10 Tong Hengqing. Evaluation Model and Its Iterative Algorithm by Alternating Projection, Mathematical and Computer Modelling, Vol. 18, No. 8, 1993, 55-60.11 童恒庆,理论计量经济学,科学出版社,2005.12 童恒庆,数据分析与统计计算软件 DASC,科学出版社,2005.