马尔可夫链模型

用手机看条目

出自 MBA智库百科(https://wiki.mbalib.com/)

马尔可夫链模型(Markov Chain Model)

目录

马尔可夫链模型概述

  马尔可夫链安德烈·马尔可夫Andrey Markov,1856-1922)得名,是数学中具有马尔可夫性质的离散时间随机过程。该过程中,在给定当前知识或信息的情况下,过去(即当期以前的历史状态)对于预测将来(即当期以后的未来状态)是无关的。

  时间和状态都是离散的马尔可夫过程称为马尔可夫链, 简记为{X_n=X(n),n=0,1,2,\cdots}

  马尔可夫链是随机变量X_1,X_2,X_3\cdots的一个数列。这些变量的范围,即他们所有可能取值的集合,被称为“状态空间”,而Xn的值则是在时间n的状态。如果Xn + 1对于过去状态的条件概率分布仅是Xn的一个函数,则

  P(X_{n+1}=x|X_0,X_1,X_2,\cdots,X_n)=P(X_{n+1}=x|X_n)

  这里x为过程中的某个状态。上面这个恒等式可以被看作是马尔可夫性质。

  马尔可夫在1906年首先做出了这类过程 。而将此一般化到可数无限状态空间是由柯尔莫果洛夫在1936年给出的。

  马尔可夫链与布朗运动以及遍历假说这两个二十世纪初期物理学重要课题是相联系的,但马尔可夫寻求的似乎不仅于数学动机,名义上是对于纵属事件大数法则的扩张。

  马尔可夫链是满足下面两个假设的一种随机过程:

  1、t+l时刻系统状态的概率分布只与t时刻的状态有关,与t时刻以前的状态无关;

  2、从t时刻到t+l时刻的状态转移与t的值无关。一个马尔可夫链模型可表示为=(S,P,Q),其中各元的含义如下:

  1)S是系统所有可能的状态所组成的非空的状态集,有时也称之为系统的状态空间,它可以是有限的、可列的集合或任意非空集。本文中假定S是可数集(即有限或可列)。用小写字母i,j(或Si,Sj)等来表示状态。

  2)P=[P_{ij}]_{n\times n}是系统的状态转移概率矩阵,其中Pij表示系统在时刻t处于状态i,在下一时刻t+l处于状态i的概率,N是系统所有可能的状态的个数。对于任意i∈s,有\sum_{j=1}^NP_{ij}=1

  3)Q=[q_1,q_2\cdots q_n]是系统的初始概率分布,qi是系统在初始时刻处于状态i的概率,满足\sum_{i=1}^Nq_i=1

马尔可夫链模型的性质

  马尔可夫链是由一个条件分布来表示的

  P(Xn + 1 | Xn)

  这被称为是随机过程中的“转移概率”。这有时也被称作是“一步转移概率”。二、三,以及更多步的转移概率可以导自一步转移概率和马尔可夫性质:

  P(X_{n+2}|X_n) = \int P(X_{n+2},X_{n+1}|X_n)dX_{n+1} = \int P(X_{n+2}|X_{n+1})P(X_{n+1}|X_n)dX_{n+1}

  同样:

  P(X_{n+3}|X_n) = \int P(X_{n+3}|X_{n+2}) \int P(X_{n+2}|X_{n+1})P(X_{n+1}|X_n)dX_{n+1}dX_{n+2}

  这些式子可以通过乘以转移概率并求k−1次积分来一般化到任意的将来时间n+k。

  边际分布P(Xn)是在时间为n时的状态的分布。初始分布为P(X0)。该过程的变化可以用以下的一个时间步幅来描述:

  P(X_{n+1}) = \int P(X_{n+1}|X_n)P(X_n)dX_n

  这是Frobenius-Perron equation的一个版本。这时可能存在一个或多个状态分布π满足:

  \pi(X) = \int P(X|Y)\pi(Y)dY

  其中Y只是为了便于对变量积分的一个名义。这样的分布π被称作是“平稳分布”(Stationary Distribution)或者“稳态分布”(Steady-state Distribution)。一个平稳分布是一个对应于特征根为1的条件分布函数的特征方程。

  平稳分布是否存在,以及如果存在是否唯一,这是由过程的特定性质决定的。“不可约”是指每一个状态都可来自任意的其它状态。当存在至少一个状态经过一个固定的时间段后连续返回,则这个过程被称为是“周期的”。

离散状态空间中的马尔可夫链模型

  如果状态空间是有限的,则转移概率分布可以表示为一个具有(i,j)元素的矩阵,称之为“转移矩阵”:

  Pij = P(Xn + 1 = i | Xn = j)

  对于一个离散状态空间,k步转移概率的积分即为求和,可以对转移矩阵求k次幂来求得。就是说,如果\mathbf{P}是一步转移矩阵,\mathbf{P}^k就是k步转移后的转移矩阵。

  平稳分布是一个满足以下方程的向量:

  \mathbf{P}\pi^* = \pi^*

  在此情况下,稳态分布π * 是一个对应于特征根为1的、该转移矩阵的特征向量。

  如果转移矩阵\mathbf{P}不可约,并且是非周期的,则\mathbf{P}^k收敛到一个每一列都是不同的平稳分布π * ,并且,

  \lim_{k\rightarrow\infty}\mathbf{P}^k\pi=\pi^*

  独立于初始分布π。这是由Perron-Frobenius theorem所指出的。

  正的转移矩阵(即矩阵的每一个元素都是正的)是不可约和非周期的。矩阵被称为是一个随机矩阵,当且仅当这是某个马尔可夫链中转移概率的矩阵。

  注意:在上面的定式化中,元素(i,j)是由j转移到i的概率。有时候一个由元素(i,j)给出的等价的定式化等于由i转移到j的概率。在此情况下,转移矩阵仅是这里所给出的转移矩阵的转置。另外,一个系统的平稳分布是由该转移矩阵的左特征向量给出的,而不是右特征向量。

  转移概率独立于过去的特殊况为熟知的Bernoulli scheme。仅有两个可能状态的Bernoulli scheme被熟知为贝努利过程

马尔可夫链模型的应用

科学中的应用

  马尔可夫链通常用来建模排队理论统计学中的建模,还可作为信号模型用于熵编码技术,如算法编码。马尔可夫链也有众多的生物学应用,特别是人口过程,可以帮助模拟生物人口过程的建模。隐蔽马尔可夫模型还被用于生物信息学,用以编码区域或基因预测。

  马尔可夫链最近的应用是在地理统计学(geostatistics)中。其中,马尔可夫链用在基于观察数据的二到三维离散变量的随机模拟。这一应用类似于“克里金”地理统计学(Kriging geostatistics),被称为是“马尔可夫链地理统计学”。这一马尔可夫链地理统计学方法仍在发展过程中。

人力资源中的应用

  马尔可夫链模型主要是分析一个人在某一阶段内由一个职位调到另一个职位的可能性,即调动的概率。该模型的一个基本假设就是,过去的内部人事变动的模式和概率与未来的趋势大体相一致。实际上,这种方法是要分析企业内部人力资源的流动趋势和概率,如升迁、转职、调配或离职等方面的情况,以便为内部的人力资源的调配提供依据。

  它的基本思想是:通过发现过去组织人事变动的规律,以推测组织在未来人员的供给情况。马尔可夫链模型通常是分几个时期收集数据,然后再得出平均值,用这些数据代表每一种职位中人员变动的频率,就可以推测出人员变动情况。

  具体做法是:将计划初期每一种工作的人数量与每一种工作的人员变动概率相乘,然后纵向相加,即得到组织内部未来劳动力的净供给量。其基本表达式为:

  N_i(t)=\sum_{j=1}^k*P_{ij}+V_i(t)

Ni(t):t时间内I类人员数量;
Pji:人员从j类向I类转移的转移率;
Vi(t):在时间(t-1,t)I类所补充的人员数。

  企业人员的变动有调出、调入、平调晋升与降级五种。表3 假设一家零售公司在1999至2000年间各类人员的变动情况。年初商店经理有12人,在当年期间平均90%的商店经理仍在商店内,10%的商店经理离职,期初36位经理助理有 11%晋升到经理,83%留在原来的职务,6%离职;如果人员的变动频率是相对稳定的,那么在2000年留在经理职位上有11人(12×90%),另外,经理助理中有4人(36×83%)晋升到经理职位,最后经理的总数是15人(11+4)。可以根据这一矩阵得到其他人员的供给情况,也可以计算出其后各个时期的预测结果。假设的零售公司的马尔可夫分析,见下表:

1999~2000商店经理经理助理区域经理部门经理销售员离职
商店经理
(n=12)
90%
11
10%
1
经理助理
(n=36)
11%
4
83%
30
6%
2
区域经理
(n=96)
11%
11
66%
63
8%
8
15%
14
部门经理
(=288)
10%
29
72%
207
2%
6
16%
46
销售员
(=1440)
6%
86
74%
1066
25%
228
供给预测1541923011072351

马尔可夫模型案例分析[1]

  案例:在信用卡账户行为变化预测中的应用

  信用卡业务是商业银行的零售业务,信用卡的消费金额是银行的应收账款.在此,我们可以借鉴零售行业应收账款状态变化的预测方法对信用卡账户的行为变化进行描述和预测。

  对信用卡账户的马尔可夫过程进行研究,主要解决新增贷款发生周期性变化的情况下利用马尔可夫过程预测不同时刻的信用卡账户各状态下的金额、已偿付态和坏帐态的金额、全部应收款的现值及它们的方差计算等内容,以为商业银行信用卡账户的行为风险管理提供方法依据。

马尔可夫模型的建立

  马尔可夫状态转移模型是在满足“马氏性”和“平稳性”的基础上建立的.假定银行的信用卡账户中每期处于不同期限的逾期贷款数量只与上期逾期贷款的数量与结构有关,而与前期的状态无关,这就满足了“马氏性”。同时,在外部经济环境稳定、人口特征比较稳定、银行的信用卡管理技术和方法没有发生重大变化的情况下,可以认为逾期贷款由一种状态转移到另一种状态的概率在各期是保持不变的,即每年的转移概率矩阵基本保持稳定,满足了马氏链的“平稳性”要求.这样,银行就可以通过往年的数据资料模拟出比较精确的转移概率矩阵,对信用卡账户的行为状态做出预测和评估,下面给出具体分析。

  假设某一银行在时间i有一定的信用卡应收账款,当前或者随后的时间内这些余额都可以划分为n个时间段(即状态。对于这批在时间i的应收账款而言,有:

  B0=逾期为0期的应收账款余额(也就是当前期);

  B1 = 逾期为1期的应收账款余额;

  …

  Bj = 逾期为j期的应收账款余额;

  …

  Bn − 1 = 逾期为n-1期的应收账款余额;

  Bn = 逾期为n期的应收账款余额。

  实践中,时间段的数目将视情况而定,最后一个时间段主要依赖于银行应收账款的“冲销”原则,美国的信用卡贷款一般拖欠180天以上即成为呆账予以“冲销”.虽然拖欠账款最终也可能得到偿还,但是将超过规定还款期限的应收账款归入坏帐种类中是很自然的会计程序

  一般而言,我们可以让Bjk表示从i时刻处于j状态转移到i+1时刻处于k状态的账户的金额.用这种方法,我们可以对处于i时刻的所有应收账款做出在i+1时刻的一步转移账户.需要注意的是,还应该有一个“时间”状态应该加入到先前所描述的分类中,这一状态就是已付款状态,用\overline{0}表示.在i时刻任何一种分类状态从0到n的账户在i+1时刻都可以转移到状态\overline{0}.这样,i时刻的应收账款账户可以用一个n+2维矩阵来表示,矩阵中的每一项Bjk表示i时刻j状态转移为i+1时刻k状态的金额,如下所示:

  B=\begin{bmatrix}B_{\overline{00}}&\cdots & B_{\overline{0}k}& \cdots & B_{\overline{0}n} \\ \bullet& &\bullet& &\bullet \\B_{j\overline{0}} &\cdots & B_{jk} & \cdots & B_{jn} \\ \bullet & &\bullet & &\bullet \\ B_{n\overline{0}} &\cdots & B_{nk} &\cdots & B_{nn}\end{bmatrix}

  对信用卡账户而言,需要注意的是,当状态Bjk中的j<i时,应理解为i时刻处于状态j的账户,在随后的i+1时刻(一般为30天后)偿还了部分的利息,使得应收账款(贷款)又转变为k状态。

  从n+2维应收账款矩阵B可以导出n+2维转移概率矩阵P.转移概率矩阵P中的每一项目表示在特定时间内某一账户由一种状态转移到另一状态的可能性.这样的话,一个隐含假设是,转移概率矩阵的考察周期和应收账款分类的考察周期是相同的.一般情况下,转移概率Pjk表示的是i时刻j状态的账款转移到i+1时刻k状态账款的可能性.根据应收账款矩阵B及Bjk,转移概率Pjk可被定义为:

  P_{jk}=\frac{B_{jk}}{{\sum^n_{\overline{0}}B_{jk}}'}  (k=\overline{0},0,1,\cdots,n)  (1)

  在应用转移概率矩阵时需要注意两点。一是\overline{0}状态的账款不可能转移到其它的状态,它只能停留在已付款状态,\overline{0}状态账户的转移概率依次为:p_{\overline{00}}=1.00,p_{\overline{0}0}=0,p_{\overline{0}1}=0,…,p_{\overline{0}k}=0,…,p_{\overline{0}n}=0。二是呆账类账户的状态,虽然有时候坏呆账类账款仍能收回现金,但在我们的模型里边假设呆账类账款只能停留在呆账类的状态,即:p_{n\overline{0}}=0,pn0 = 0,pn1 = 0,…,pnn = 1.00

  上面描述的模型可以被看作一个有n+2个状态的马尔可夫链过程,其转移概率矩阵为P.而且,它有两个吸收态(偿付态0和呆账态n),从其他任何一个暂态(非吸收态)都可以到达这两个吸收态,因此它是一个具有两个吸收态的马尔可夫链.我们将在充分利用马尔可夫理论和已有研究的基础上,研究如何利用马尔可夫链方法预测和估计信用卡账户行为的变化。

马尔可夫模型的应用

  在此,采用Kemeny和Snell的部分研究成果.为便于计算,将n+2维转移概率方阵重新排列,将吸收态的偿付态和呆账态放在一起,将另外的暂态0,1,2,…,n-1放在一起.这样矩阵P就可以被分割为:

  P=\begin{bmatrix} I & O \\ R & Q \end{bmatrix}

  其中I是一个2×2阶单位矩阵,O是一个2×n阶0矩阵,R是一个n×2阶矩阵,Q是一个n×n阶矩阵.其中,我们定义矩阵:

  N=(1-Q)^{-1}=I+Q+Q^2+Q^3+\cdots +Q^k+\cdots

  一定存在,并将其称为吸收态马尔可夫链的基本矩阵

  对于n×2阶矩阵的所有分项,N R给出了每一状态转移到吸收态\overline{0}和n的吸收概率.NR中的第一列给出了每一个状态转移到已偿付状态的概率,第二列给出了每一个状态下转移到呆账的概率。

  1.无新增贷款的情况

  假设在时刻i,具有n个分项向量的B_i=(B_{i0},B_{i1},\cdots,B_{in-1})给出来每一状态下应收账款的余额.让b等于所有这些余额之和,则向量\pi=(\frac{1}{b})B是一个没有非负分量且全部之和为1的概率向量,向量的分量代表了每一状态下应收账款的比例.如果我们假设上述状态中的余额的移动是独立的,那么我们就可以认定向量π为马尔可夫链的初始向量.另外,还假定:如果A是任一矩阵,那么我们让Asq表示A中每一项平方后的结果;让Art表示A中每一项取平方根后的结果.则有如下结论:

  结论1

  二维向量BNR中的分量可以给出来自应收账款向量B的期望还款和坏帐金额;分量给出来偿还态和呆帐态的方差,Art给出了这两种状态的标准差。

  A=b\begin{bmatrix}\pi N R-(\pi N R)_{sq}\end{bmatrix}(2)

  证明 如上所述,矩阵NR中第一列的分量给出来应收账款从每一暂态转移到吸收态(偿付态)的概率.向量\pi=(\frac{1}{b})B的分量给出了每次过程开始时账款转移到每一暂态的初始概率.因此,账款在最终时偿付态的概率可以由向量πNR的第一列分量给出.如果这一过程开始了b次,那么在最终时偿付态的平均数就是向量bπNR = BNR的第一列分量.向量πNR的第一分量是函数f的平均值,其中f表示在最终结束时偿付态的价值为全部价值,其它状态的价值为零.这一函数的方差可以由下式的第一分量给出:

  V(f)=M(f^2)-\begin{bmatrix} M(f)\end{bmatrix}^2

  因为f2 = f,所以M(f2) = M(f),因此f的方差可以由πNR − (πNR)sq的第一分量给出.如果过程开始了b次,那么偿付态的全部金额的方差可以由A=b\begin{bmatrix}\pi N R-(\pi N R)_{sq}\end{bmatrix}的第一分量给出.有关呆帐态的分析与偿付态的分析类似。

  此外,还可以对应收账款现值的计算进行了研究. 如果 r是利率,则\beta=\frac{1}{(1+r)}就表示了贴现率,应收账款现值的计算就可以由下面的计算给出。

  假定B是应收账款向量,R1是矩阵R的第一列分量,则BR1表示当前时期的收现额;从下一期的BQR1的价值就只有BBQR1;依此类推,在(k+1)周期时BQkR1的价值就只有ΒkBQkR1.将这些折现价值加在一起就可以得到应收账款的当前现值:

  B R_1+\beta B Q R_1+\cdots+\beta^k B Q^k R_1+\cdots =B\begin{bmatrix} I+\beta Q+\cdots+\beta^k Q^k+\cdots\end{bmatrix}R_1=B N_\beta R_1,其中的Nβ表示I+\beta Q+\cdots+\beta^k Q^k+\cdots

  在实践当中,银行一般都要对信用卡客户收取一定的年费,假定银行对客户收取b的费率,则β = 1 + b,那么完全可以利用上述公式来计算应收账款的现值.当然,如果考虑利率和年费率两种因素的话,将会有一个净折扣率或者一个费用率。

  2.新增贷款固定不变的情况

  假设每期又发生了金额为c的新应收款,这些新应收款被分不在不同的状态下,构成了向量C的各分量组成,即:C=(C_0,C_1,\cdots,C_{n-1}).定义向量\eta=(\frac{1}{c})C,则η为概率向量并且被认为是马尔可夫链的初始向量.假设,马尔可夫过程每期以初始概率η开始了c次.那么应收账款的稳定态分布会怎么样,这些账户的方差又是多少?每期期望付款和呆账的数量以及它们的期望方差又怎么样?

  结论2

  如果马尔可夫过程每期以初始概率η开始了c次,则向量CN的分量给出来所有时刻下稳定的应收账款金额,数值CNξ给出了稳定态的全部应收账款金额,其中ξ是各项为1的n维列向量.二维向量CNR给出来每期偿付款和呆账的稳定态的金额。

  证明 如果上述马尔可夫过程进行了许多个周期,则各状态的金额由当前η一个月前的ηQ、二个月前的ηQ2,等等组成.那么这些数量之和为:

  \eta+\eta Q^2+\eta Q^3+\cdots=\eta(I+Q+Q^2+Q^3+\cdots)=\eta N

  如果这个过程每周期开始了c次,每一状态下的应收账款可以由向量cηN = CN表示.如果ξ是一个各项为1的列向量,则CNξ是向量CN的分量之和,代表了应收账款的全部账户余额.

  如果上述过程进行了很多周期,将会有ηR的账款从第一期的新收款中转移到吸收态,将有ηQR的账款从接下来的一期的新收款中转移到吸收态,将有ηQ2R的账款从过期两个月的新收款中转移到吸收态,依此类推,那么所有这些之和为:

  \eta R+\eta Q^2R+\eta Q^3 R+\cdots=\eta(I+Q+Q^2+Q^3+\cdots)R=\eta NR

  如果这一过程开始了c次,每期稳定态的偿付款和呆账将有cηNR = CNR给出。证明完毕。

  综合定理1和定理2,我们能够得出一下推论.让t = CNξ,\pi=(\frac{1}{t})CN;那么CN2Rt\begin{bmatrix}\pi NR-(\pi NR)_{sq}\end{bmatrix}是偿付款和呆账的预测均值和方差.而且,可以根据对应收款的利率和费率来计算应收账款的现值。

  3.新增贷款发生周期性变化的情况

  上述讨论都没有考虑应收账款发生变化的情况,然而,在现实情况下,银行的信用卡消费呈现出一定的周期性,例如在春节、国庆节和秋季开学的时候消费比较高. 除此之外,商业银行每年的消费贷款也可能因为经济增长或萧条等原因而扩张或收缩. 因此,我们需要考虑这些因素对模型的一些影响.

  具体来讲,让Ci是给定月份i的新应收款的向量; ci是全部应收款的金额; η = (1 / ci)Ci是第 i时刻的初始向量,假设:

  ηiT = (ηi)   (3)

  CiT = αCi   (4)

  其中α是增长系数的倒数,例如某一贷款机构的信用卡业务以2%的年增长率扩张则α = 1 / (1 + 0.02) = 1 / 1.02T 为循环周期的长度,一般情况下周期T = 12. 从上面的两个式子里边我们可以推出ciT = αci

  结论3

  让Nα = (I − αQT) − 1,那么下列式子:

  A_i= \left[\sum^{r-1}_{k=0}C_{i-k}Q^k\right]N_{\alpha}  (5)

  \alpha_ii=\left[\sum^{t-1}_{k=0}C_{i-k}Q^k\right]N_{\alpha} \xi  (6)

  D_i=\left[\sum^{t-1}_{k=0}C_{i-k}Q^k\right]N_{\alpha}R  (7)

  给出了 i时刻不同状态下的金额、全部应收账款、以及吸收态的金额.

  证明 让C_i,C_{i-1},\cdots,C_{i-T+1}是第i月份及其之前T-1月的真实新收款. 在知道增长率的情况下,根据(4)式能够推出以前月份的所有应收款,其中第i月份不同状态的应收款是

  Ci;第(i-1)月份的是Ci − 1Q;第(i-2) 月份的是Ci − 2Q2,等等; 第(i-T+1) 月份的是

  CiT + 1QT − 1;第(i-T )月份的是CiQT(CiT = αCi),等等. 将这些向量加总后如下:

  A_i= C_i+ C_{i-1}Q+C_{i-2}Q^2+\cdots+C_{i-T+1}Q^{T-1}+\alpha C_iQ^T+\alpha C_{i-1}Q^{T+1}+\cdots+\alpha C_{i-T+1}Q^{2T-1}+\alpha^{2}C_iQ^{2T}+\alpha C_{i-1}Q^{2T+1}+\cdots + \alpha C_{i-T+1}Q^{3T-1}+\cdots

  =C_i(I+\alpha Q^T+\alpha^2Q^{2T}+\cdots) + C_{i-1}Q(I+\alpha Q^T+\alpha^2Q^{2T}+\cdots)+\cdots + C_{i-T+1}Q^{T+1}(I+\alpha Q^T+\alpha^2Q^{2T}+\cdots)

  = [C_i+C_{i-1}Q+C_{i-2}Q^2+\cdots + C_{i-T+1}Q^{T-1}]N_{\alpha}

  这就是AiiDi的证明与Ai类似.

  当然,对于 i时刻的这些估计依赖于第 i月及其前T - 1月的新增应收款,上面给出的估计结果比结论2给出的结果更准确一些. 当然,如果Qn快速趋于0,则用过去几个月的应收账款来估计一个合理的结果也是可以的.

  根据结论1和结论3的结论,我们可以用AiNR\alpha_i\left[\tau_iN R- (\tau_iN R)_{sq}\right],其中αi = Aiξτi = (1 / alphai)Ai来估计 i时刻偿付款和呆帐的均值和方差,而且也可以用AiNβR1用来估计 i时刻应收账款的现值.

参考文献

  1. 任金政 陈宝峰 庄传礼.马尔可夫链模型在信用卡账户行为变化预测中的应用.数学的实践与认识.2008年5月第9期第38卷
本条目对我有帮助134
MBA智库APP

扫一扫,下载MBA智库APP

分享到:
  如果您认为本条目还有待完善,需要补充新内容或修改错误内容,请编辑条目

本条目由以下用户参与贡献

Zfj3000,Cabbage,Vulture,Dan,Angle Roh,Way,Yixi,鲈鱼,HEHE林,4Fires哥哥.

评论(共4条)

提示:评论内容为网友针对条目"马尔可夫链模型"展开的讨论,与本站观点立场无关。
59.55.38.* 在 2010年8月4日 22:46 发表

好!

回复评论
218.26.221.* 在 2010年9月14日 10:02 发表

回复评论
155.69.2.* 在 2010年10月6日 17:29 发表

如果转移矩阵不可约,并且是非周期的,则收敛到一个每一列都是不同的平稳分布π *,怎么证明?

回复评论
14.109.116.* 在 2011年12月11日 21:09 发表

不错

回复评论

发表评论请文明上网,理性发言并遵守有关规定。

打开APP

以上内容根据网友推荐自动排序生成

闽公网安备 35020302032707号