支持向量机的理论支持在此不细说，可以参考李航的《统计学习》，还有西瓜书。
简化版SMO算法处理小规模数据集
SMO算法是一种启发式算法。此简化版首先在数据集上遍历每一个alpha，然后在剩下的alpha集合中随机选择另一个alpha，从而建立alpha对。
# -*- coding: utf-8 -*-
from numpy import *
from time import sleep
# SMO算法的辅助函数
def loadDataSet(fileName): #加载并预处理数据集
dataMat = []; labelMat = []
fr = open(fileName,'r')
for line in fr.readlines():
lineArr = line.strip().split('\t') # 以制表符分割
dataMat.append([float(lineArr[0]), float(lineArr[1])]) #提取前两个元素存入data.Mat中
labelMat.append(float(lineArr[2])) # [].append(),最终的形式是矩阵
return dataMat,labelMat
def selectJrand(i,m): # 该辅助函数用于在某个区间范围内随机选择一个整数
j=i # m是所有alpha的数目，i是第一个alpha的下标
while (j==i):
j = int(random.uniform(0,m)) # random.uniform(0,m)用于生成指定范围内的随机浮点数
return j
def clipAlpha(aj,H,L): # 该辅助函数用于在数值太大时对其进行调整
if aj > H:
aj = H
if L > aj:
aj = L
return aj
# 简化版SMO算法
def smoSimple(dataMatIn, classLabels, C, toler, maxIter): # 参数：数据集，类别标签，常数c，容错率，循环次数
dataMatrix = mat(dataMatIn) # mat()转换成矩阵类型
labelMat = mat(classLabels).transpose() #转置之前是列表，转置后是一个列向量
b = 0; m,n = shape(dataMatrix) # 得到行，列数，m行，n列
alphas = mat(zeros((m,1))) # zeros(shape, dtype=float, order='C'),所以也可以写作zeros((10,1),)
iter = 0 # 该变量存储的是在没有任何alpha改变时遍历数据集的次数
while (iter < maxIter): # 限制循环迭代次数，也就是在数据集上遍历maxIter次，且不再发生任何alpha修改，则循环停止
alphaPairsChanged = 0 # 每次循环时先设为0，然后再对整个集合顺序遍历，该变量用于记录alpha是否已经进行优化
for i in range(m): # 遍历每行数据向量，m行
# 该公式是分离超平面，我们预测值
fXi = float(multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[i,:].T)) + b
#print 'fxi:',fxi
Ei = fXi - float(labelMat[i]) # 预测值和真实输出之差
# 如果误差很大就对该数据对应的alpha进行优化，正负间隔都会被测试，同时检查alpha值
if ((labelMat[i]*Ei < -toler) and (alphas[i] < C)) or ((labelMat[i]*Ei > toler) and (alphas[i] > 0)):
j = selectJrand(i,m) # 随机选择不等于i的0-m的第二个alpha值
fXj = float(multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[j,:].T)) + b
Ej = fXj - float(labelMat[j])
alphaIold = alphas[i].copy(); alphaJold = alphas[j].copy();
if (labelMat[i] != labelMat[j]): # 这里是对SMO最优化问题的子问题的约束条件的分析
L = max(0, alphas[j] - alphas[i]) # L和H分别是alpha所在的对角线端点的界
H = min(C, C + alphas[j] - alphas[i]) # 调整alphas[j]位于0到c之间
else:
L = max(0, alphas[j] + alphas[i] - C)
H = min(C, alphas[j] + alphas[i])
if L==H: print "L==H"; continue # L=H停止本次循环
# 是一个中间变量：eta=2xi*xi-xixi-xjxj，是alphas[j]的最优修改量
eta = 2.0 * dataMatrix[i,:]*dataMatrix[j,:].T - dataMatrix[i,:]*dataMatrix[i,:].T \
- dataMatrix[j,:]*dataMatrix[j,:].T
if eta >= 0: print "eta>=0"; continue # eta>=0停止本次循环，这里是简化计算
alphas[j] -= labelMat[j]*(Ei - Ej)/eta # 沿着约束方向未考虑不等式约束时的alpha[j]的解
alphas[j] = clipAlpha(alphas[j],H,L) # 此处是考虑不等式约束的alpha[j]解
if (abs(alphas[j] - alphaJold) < 0.00001):
print "j not moving enough"; continue # 如果该alpha值不再变化，就停止该alpha的优化
alphas[i] += labelMat[j]*labelMat[i]*(alphaJold - alphas[j]) # 更新alpha[i]
# 完成两个alpha变量的更新后，都要重新计算阈值b
b1 = b - Ei- labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i,:]*dataMatrix[i,:].T \
- labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[i,:]* dataMatrix[j,:].T #李航统计学习7.115式
b2 = b - Ej- labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i,:]*dataMatrix[j,:].T \
- labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[j,:]*dataMatrix[j,:].T #李航统计学习7.116式
if (0 < alphas[i]) and (C > alphas[i]): b = b1
elif (0 < alphas[j]) and (C > alphas[j]): b = b2
else: b = (b1 + b2)/2.0 # alpha[i]和alpha[j]是0或者c,就取中点作为b
alphaPairsChanged += 1 # 到此的话说明已经成功改变了一对alpha
print "iter: %d i:%d, pairs changed %d" % (iter,i,alphaPairsChanged)
if (alphaPairsChanged == 0): iter += 1 # 如果alpha不再改变迭代次数就加1
else: iter = 0
print "iteration number: %d" % iter
return b,alphas
# 主函数
dataArr,labelArr=loadDataSet('testSet.txt') # 因为在同一文件夹下，就不用写绝对路径
#print dataArr
b,alphas=smoSimple(dataArr, labelArr, 0.6, 0.001, 40)
print 'b:',b
print 'alphas[alphas>0]:',alphas[alphas>0] # 数组过滤
print shape(alphas[alphas>0]) # 得到支持向量的个数
for i in range(100): # 得到是支持向量的数据点
if alphas[i]>0.0: print dataArr[i],labelArr[i]其中要注意的python语法：
- 数组和矩阵的转换
In [17]:type(dataArr) # 转换之前是list的类型
Out[17]: list
In [18]:mat(dataArr) # 转换之后是矩阵类型
Out[18]:
matrix([[ 3.542485, 1.977398],
[ 3.018896, 2.556416],
[ 7.55151 , -1.58003 ],
...,
[ 2.912122, -0.202359],
[ 1.731786, 0.589096],
[ 2.387003, 1.573131]])- Python自带的copy()，deepcopy()，numpy的copy（）之间的区别
>>> import copy
>>> origin = [1, 2, [3, 4]]
>>> cop1 = copy.copy(origin) # 浅拷贝
>>> cop2 = copy.deepcopy(origin) # 深拷贝
>>> cop1 == cop2
True
>>> cop1 is cop2
False #cop1 和 cop2 看上去相同，但已不再是同一个object
>>> origin[2][0] = "hey!"
>>> origin
[1, 2, ['hey!', 4]]
>>> cop1
[1, 2, ['hey!', 4]]
>>> cop2
[1, 2, [3, 4]] #把origin内的子list [3, 4] 改掉了一个元素，观察 cop1 和 cop2可以看到 cop1，也就是 shallow copy 跟着 origin 改变了。而 cop2 ，也就是 deep copy 并没有变。
似乎 deep copy 更加符合我们对「复制」的直觉定义: 一旦复制出来了，就应该是独立的了。如果我们想要的是一个字面意义的「copy」，那就直接用 deep_copy 即可。
注意：这里指的是python自带的copy( )函数，copy.copy(object)是浅拷贝，而numpy的copy( )函数，即object.copy( )是深拷贝的效果。
from numpy import *
In [26]:origin
Out[26]:
matrix([[1, 2],
[3, 4]])
In [28]:origin[0][0,1]
Out[28]: 2
In [29]:origin[0][0,0]
Out[29]: 1
In [31]:old=origin[0][0,1].copy()
In [32]:origin[0][0,1]=5
In [33]:origin[0][0,1]
Out[33]: 5
In [34]:old
Out[34]: 2分类图示：
在园圈中的就是支持向量
代码中的公式含义：
（1）
fXi = float(multiply(alphas,labelMat).T*(dataMatrix*dataMatrix[i,:].T)) + bfxi 对应的是《统计学习》中7.104的 g(x) ：
g(x)=∑i=1Nα∗iyixi⋅x+b∗
其代表的是分离超平面
（2）
Ei = fXi - float(labelMat[i])
Ei=g(xi)−yi
表示预测值和真实输出之差
（3）
if (labelMat[i] != labelMat[j]):
L = max(0, alphas[j] - alphas[i])
H = min(C, C + alphas[j] - alphas[i])
else:
L = max(0, alphas[j] + alphas[i] - C)
H = min(C, alphas[j] + alphas[i])L和H是alphas的所在对角线端点的界：
如果
y1≠y2
则
L=max(0,αoldj−αoldi),H=min(c,c+αoldj−αoldi)
如果
y1=y2
则
L=max(0,αoldj+αoldi−c),H=min(c,αoldj+αoldi)
（4）
eta = 2.0 * dataMatrix[i,:]*dataMatrix[j,:].T - dataMatrix[i,:]*dataMatrix[i,:].T \
- dataMatrix[j,:]*dataMatrix[j,:].T对应
η=2xj∗xi−xi∗xi−xj∗xj
表示的是alphas[j]的最优修改量，是一个中间变量
（5）
alphas[j] -= labelMat[j]*(Ei - Ej)/eta 对应
αnew,uncj=αoldj−yj(Ei−Ej)η
是沿着约束方向未考虑不等式约束时的alpha[j]的解
（6）
alphas[j] = clipAlpha(alphas[j],H,L) 此处是考虑不等式约束的alpha[j]解：
αnewj=⎧⎩⎨⎪⎪⎪⎪Hαnew,uncjLαnew,uncj>HL⩽αnew,uncj⩽Hαnew,uncj<L
得到了第一个变量的值
（7）
alphas[i] += labelMat[j]*labelMat[i]*(alphaJold - alphas[j])
αnewi=αoldi+y1y2(αoldj−αnewj)
（8）
b1 = b - Ei- labelMat[i]*(alphas[i]-alphaIold)*dataMatrix[i,:]*dataMatrix[i,:].T \
- labelMat[j]*(alphas[j]-alphaJold)*dataMatrix[i,:]* dataMatrix[j,:].T得到的是b：
bnew1=bold−Ei−yixi∗xi(αnewi−αoldi)−y2xj∗xi(αnewj−αoldj)+bold
代码中的主要公式大概就这些。
运行结果：
j not moving enough
j not moving enough
iteration number: 1
j not moving enough
j not moving enough
...,
j not moving enough
iteration number: 38
j not moving enough
j not moving enough
iteration number: 39
j not moving enough
j not moving enough
iteration number: 40
b: [[-3.83495394]]
alphas[alphas>0]: [[ 0.15601002 0.14181599 0.06893826 0.36676427]]
(1L, 3L)
[4.658191, 3.507396] -1.0
[3.457096, -0.082216] -1.0
[6.080573, 0.418886] 1.0“不忘初心，方得始终”搞了那么多代码和公式，不要沉浸进去而不知道最后要求的是什么？SMO算法的目标是求出一系列的
α
和
b
，而此时运行得到的就是那些支持向量和相对应的参数
∑i=1Nα∗iyixi⋅x+b∗=0
分类决策函数可以写成：
f(x)=sign(∑i=1Nα∗iyixi⋅x+b∗)
由此验证了在决定分离超平面时只有支持向量起作用。如果移动支持向量将改变所求的解；但如果在间隔边界以外移动其他实例点，甚至去掉这些点，则解不会改变。
由于SMO算法的随机性，每次的运行结果可能不同。
还有一点要注意：就是优化结束的同时必须确保合适的时机结束循环，如果程序执行到for循环的最后一行都不执行continue语句，那么就成功地改变了一对 α 值，同时可以增加alphaPairsChanged的值，在for循环之外，需要检查 α 值是否有更新，如果有更新则将iter设为0后继续运行程序。只有在所有的数据集上遍历maxIter次，且不再发生任何 α 的修改后，程序才停止并退出while循环。
利用完整版platt SMO算法加速优化
在完整版的SMO算法中，实现alpha的更改和代数运算的优化环节一模一样，在优化过程中，唯一不同的是选择alpha的方式，完整版的应用了一些能够提速的启发式方法。
# -*- coding: utf-8 -*-
"""
Created on Wed Oct 18 20:53:40 2017
@author: LiLong
"""
from numpy import *
from time import sleep
# SMO算法的辅助函数
def loadDataSet(fileName): #加载并预处理数据集
dataMat = []; labelMat = []
fr = open(fileName,'r')
for line in fr.readlines():
lineArr = line.strip().split('\t') # 以制表符分割
dataMat.append([float(lineArr[0]), float(lineArr[1])]) #提取前两个元素存入data.Mat中
labelMat.append(float(lineArr[2])) # [].append(),最终的形式是矩阵
return dataMat,labelMat
# 完整版platt SMO 算法的支持函数
# 建立一个数据结构来保存所有的重要值
class optStruct:
def __init__(self,dataMatIn, classLabels, C, toler):
self.X = dataMatIn
self.labelMat = classLabels
self.C = C
self.tol = toler
self.m = shape(dataMatIn)[0] # 有多少行数据
self.alphas = mat(zeros((self.m,1)))
self.b = 0
self.eCache = mat(zeros((self.m,2))) # 误差缓存，第一列是ecache是否有效的标志位，第二列是实际的E值
def clipAlpha(aj,H,L):
if aj > H:
aj = H
if L > aj:
aj = L
return aj
def selectJrand(i,m): # 该辅助函数用于在某个区间范围内随机选择一个整数
j=i # m是所有alpha的数目，i是第一个alpha的下标
while (j==i):
j = int(random.uniform(0,m)) # random.uniform(0,m)用于生成指定范围内的随机浮点数
return j
# 计算E值并返回,E值是函数对输入xi的预测值与真实输出的差
def calcEk(oS, k):
fXk = float(multiply(oS.alphas,oS.labelMat).T*(oS.X*oS.X[k,:].T))+ oS.b
Ek = fXk - float(oS.labelMat[k])
return Ek
# 用于选择合适的第二个alpha值以保证每次优化中采用最大步长，是内循环的启发式方法
def selectJ(i, oS, Ei): # 该函数的误差值与第一个alpha值Ei和下标i有关
maxK = -1; maxDeltaE = 0; Ej = 0
oS.eCache[i] = [1,Ei] # 设置有效，有效意味着它已经计算好了
validEcacheList = nonzero(oS.eCache[:,0].A)[0] # 构建出一个非零表，返回的列表中包含以输入列表为目录的列表值
if (len(validEcacheList)) > 1:
for k in validEcacheList:
if k == i: continue # 跳出本次循环
Ek = calcEk(oS, k) # 传递对象和k，计算误差值
deltaE = abs(Ei - Ek)
if (deltaE > maxDeltaE): # 选择具有最大步长的j
maxK = k; maxDeltaE = deltaE; Ej = Ek # 会在所有的值上循环，并选择其中使得改变最大的那个值
return maxK, Ej
else: # 在这种情况下（第一次，我们没有任何有效的eCache值 ），随机选择一个alpha值
j = selectJrand(i, oS.m)
Ej = calcEk(oS, j)
return j, Ej
def updateEk(oS, k): # alpha改变时更新缓存中的值
Ek = calcEk(oS, k)
oS.eCache[k] = [1,Ek]
# 完整platt SMO算法中的优化例程
def innerL(i, oS):
Ei = calcEk(oS, i) # 计算误差值
if ((oS.labelMat[i]*Ei < -oS.tol) and (oS.alphas[i] < oS.C)) or \
((oS.labelMat[i]*Ei > oS.tol) and (oS.alphas[i] > 0)):
j,Ej = selectJ(i, oS, Ei) # 第二个alpha选择中的启发式方法
alphaIold = oS.alphas[i].copy(); alphaJold = oS.alphas[j].copy();
if (oS.labelMat[i] != oS.labelMat[j]):
L = max(0, oS.alphas[j] - oS.alphas[i])
H = min(oS.C, oS.C + oS.alphas[j] - oS.alphas[i])
else:
L = max(0, oS.alphas[j] + oS.alphas[i] - oS.C)
H = min(oS.C, oS.alphas[j] + oS.alphas[i])
if L==H: print "L==H"; return 0
eta = 2.0 * oS.X[i,:]*oS.X[j,:].T - oS.X[i,:]*oS.X[i,:].T - oS.X[j,:]*oS.X[j,:].T
if eta >= 0: print "eta>=0"; return 0
oS.alphas[j] -= oS.labelMat[j]*(Ei - Ej)/eta
oS.alphas[j] = clipAlpha(oS.alphas[j],H,L)
updateEk(oS, j) # 更新误差缓存
if (abs(oS.alphas[j] - alphaJold) < 0.00001): print "j not moving enough"; return 0
oS.alphas[i] += oS.labelMat[j]*oS.labelMat[i]*(alphaJold - oS.alphas[j])
updateEk(oS, i)
b1 = oS.b - Ei- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.X[i,:]*oS.X[i,:].T \
- oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.X[i,:]*oS.X[j,:].T
b2 = oS.b - Ej- oS.labelMat[i]*(oS.alphas[i]-alphaIold)*oS.X[i,:]*oS.X[j,:].T \
- oS.labelMat[j]*(oS.alphas[j]-alphaJold)*oS.X[j,:]*oS.X[j,:].T
if (0 < oS.alphas[i]) and (oS.C > oS.alphas[i]): oS.b = b1
elif (0 < oS.alphas[j]) and (oS.C > oS.alphas[j]): oS.b = b2
else: oS.b = (b1 + b2)/2.0
return 1 # 如果有任意一对alpha发生改变，那么就会返回1，其他返回0
else: return 0
# 完整版platt SMO的外循环代码
def smoP(dataMatIn, classLabels, C, toler, maxIter):
# 建立一个数据结构来容纳所有的数据
oS = optStruct(mat(dataMatIn),mat(classLabels).transpose(),C,toler)
iter = 0
entireSet = True; alphaPairsChanged = 0 # 退出循环的变量的一些初始化
# 迭代次数超过指定的最大值或者遍历整个集合都未对任意的alpha对进行修改时就退出循环
while (iter < maxIter) and ((alphaPairsChanged > 0) or (entireSet)):
alphaPairsChanged = 0
if entireSet:
for i in range(oS.m): # 一开始在数据集上遍历任意可能的alpha
# 选择第二个alpha，并在可能时对其进行优化处理,有任一一对alpha发生变化化了alphaPairsChanged+1
alphaPairsChanged += innerL(i,oS)
print "fullSet, iter: %d i:%d, pairs changed %d" % (iter,i,alphaPairsChanged)
iter += 1
else: # 遍历所有的非边界alpha值，也就是不在边界0或c上的值
nonBoundIs = nonzero((oS.alphas.A > 0) * (oS.alphas.A < C))[0]
for i in nonBoundIs:
alphaPairsChanged += innerL(i,oS)
print "non-bound, iter: %d i:%d, pairs changed %d" % (iter,i,alphaPairsChanged)
iter += 1
if entireSet: entireSet = False # 在非边界循环和完整遍历之间进行切换
elif (alphaPairsChanged == 0): entireSet = True
print "iteration number: %d" % iter
return oS.b,oS.alphas
# 分类超平面的w计算
def calcWs(alphas,dataArr,classLabels):
X = mat(dataArr); labelMat = mat(classLabels).transpose()
m,n = shape(X)
w = zeros((n,1))
for i in range(m):
w += multiply(alphas[i]*labelMat[i],X[i,:].T)
return w
# 主函数
dataArr,labelArr=loadDataSet('testSet.txt')
b,alphas=smoP(dataArr,labelArr,0.6,0.001,40)
print 'b:',b
print 'alphas:',alphas # 输出w和b
ws=calcWs(alphas,dataArr,labelArr)
print 'ws:',ws
datmat=mat(dataArr)
result=datmat[0]*mat(ws)+b # 进行分类
print 'result',result一些Python的技巧：
nonzero的用法
matrix.A用法
In [5]:a=mat(zeros((3,2)))
In [6]:a
Out[6]:
matrix([[ 0., 0.],
[ 0., 0.],
[ 0., 0.]])
In [7]:nonzero(a[:,0].A)[0]
Out[7]: array([], dtype=int64)
In [8]:a[:,0].A
Out[8]:
array([[ 0.],
[ 0.],
[ 0.]])
In [15]: a=mat([[0,1],[1,0],[0,1]])
In [16]: a[:,0].A
Out[16]:
array([[0],
[1],
[0]])
In [19]: nonzero(a[:,0].A)[0]
Out[19]: array([1], dtype=int64)
In [20]: a=mat([[0,1],[1,0],[1,1]])
In [21]: a[:,0].A
Out[21]:
array([[0],
[1],
[1]])
In [22]: nonzero(a[:,0].A)[0]
Out[22]: array([1, 2], dtype=int64) >>> from numpy import *
>>> aa=[0,1,2,3,4,10,6,7,8,0]
>>> b=mat(aa)
>>> b
matrix([[ 0, 1, 2, 3, 4, 10, 6, 7, 8, 0]])
>>> b.A
array([[ 0, 1, 2, 3, 4, 10, 6, 7, 8, 0]])
>>> b.A>3
array([[False, False, False, False, True, True, True, True, True,
False]], dtype=bool)
>>> b.A<4
array([[ True, True, True, True, False, False, False, False, False,
True]], dtype=bool)
>>> (b.A>3)*(b.A<4)
array([[False, False, False, False, False, False, False, False, False,
False]], dtype=bool)
>>> nonzero((b.A>3)*(b.A<4))
(array([], dtype=int64), array([], dtype=int64))
>>> nonzero((b.A>3)*(b.A<4))[0]
array([], dtype=int64)
>>> SMO中拉格朗日乘子的启发式选择方法：
所谓的启发式选择方法主要思想是每次选择拉格朗日乘子的时候，优先选择样本前面系数
0<ai<c
的
ai
作优化（论文中称为无界样例），因为在界上（
ai
为0或C）的样例对应的系数
ai
一般不会更改。
这条启发式搜索方法是选择第一个拉格朗日乘子用的，那么这样选择的话，是否最后会收敛。可幸的是Osuna定理告诉我们只要选择出来的两个 ai 中有一个违背了KKT条件，那么目标函数在一步迭代后值会减小。违背KKT条件不代表 0<ai<c ，在界上也有可能会违背。是的，因此在给定初始值 ai =0后，先对所有样例进行循环，循环中碰到违背KKT条件的（不管界上还是界内）都进行迭代更新。等这轮过后，如果没有收敛，第二轮就只针对 0<ai<c 的样例进行迭代更新。
在第一个乘子选择后，第二个乘子也使用启发式方法选择，第二个乘子的迭代步长大致正比于 |E1−E2| ，选择第二个乘子能够最大化 |E1−E2| 。即当 E1 为正时选择负的绝对值最大的 E2 ，反之，选择正值最大的 E2 。
最后的收敛条件是在界内（ 0<ai<c ）的样例都能够遵循KKT条件，且其对应的 ai 只在极小的范围内变动。
参考：
http://www.cnblogs.com/jerrylead/archive/2011/03/18/1988419.html#3793878
运行结果：
L==H
fullSet, iter: 0 i:0, pairs changed 0
L==H
fullSet, iter: 0 i:1, pairs changed 0
fullSet, iter: 0 i:2, pairs changed 1
L==H
...,
j not moving enough
fullSet, iter: 2 i:97, pairs changed 0
fullSet, iter: 2 i:98, pairs changed 0
fullSet, iter: 2 i:99, pairs changed 0
iteration number: 3
b: [[-2.89901748]]
alphas: [[ 0.06961952]
[ 0. ]
[ 0. ]
...,
[ 0. ]
[ 0. ]
[ 0. ]]
ws: [[ 0.65307162]
[-0.17196128]]
result [[-0.92555695]]得到
α
和
b
后就可以进行分类了：
完整版的SMO算法和简化版的几点区别：
- while的循环退出条件更多一些
- maxiter变量的简化版的作用不同，简化版的是当没有任何alpha发生变化时会将整个集合的一次遍历过程计成一次迭代，而完整版的的一次迭代计成定义为一次循环过程，而不管该循环具体做了什么事
- while循环内部和简化版的有所不同
- 速度有很大提升
支持向量机有些复杂，以后在应用中再进一步学习！！
3644
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
