你对乐高（LEGO）品牌的玩具了解吗？乐高公司生产拼装类玩具，由很多大小不同的塑料插块组成。这些塑料插块的设计非常出色，不需要任何粘合剂就可以随意拼装起来。除了简单玩具之外，乐高玩具在一些成人中也很流行。一般来说，这些插块都成套出售，它们可以拼装成很多不同的东西，如船、城堡、一些著名建筑，等等。乐高公司每个套装包含的部件数目从10件到5000件不等。
一种乐高套装基本上在几年后就会停产，但乐高的收藏者之间仍会在停产后彼此交易。Dangler喜欢为乐高套装估价，下面将用本章的回归技术帮助他建立一个预测模型。
示例：用回归法预测乐高套装的价格
- 收集数据：用Google Shopping的API收集数据。
- 准备数据：从返回的JSON数据中抽取价格。
- 分析数据：可视化并观察数据。
- 训练算法：构建不同的模型，采用逐步线性回归和直接的线性回归模型。
- 测试算法：使用交叉验证来测试不同的模型，分析哪个效果最好。
- 使用算法：这次练习的目标就是生成数据模型。
在这个例子中，我们将从不同的数据集上获取价格，然后对这些数据建立回归模型。需要做的第一件事就是如何获取数据。
8.6.1 收集数据：使用Google购物的API
Google已经为我们提供了一套购物的API来抓取价格。在使用API之前，需要注册一个Google账号，然后访问Google API的控制台来确保购物API可用。完成之后就可以发送HTTP请求，API将以JSON格式返回所需的产品信息。Python提供了JSON解析模块，我们可以从返回的JSON格式里整理出所需数据。详细的API介绍可以参见：http://code.google.com/apis/shopping/search/v1/getting_started.html。
打开regression.py文件并加入如下代码。
程序清单8-5 购物信息的获取函数
from time import sleepimport jsonimport urllib2def searchForSet(retX, retY, setNum, yr, numPce, origPrc): sleep(10) myAPIstr = /'get from code.google.com/' searchURL = /'https://www.googleapis.com/shopping/search/v1/public/products?key=%s&country=US&q=lego+%d&alt=json/' % (myAPIstr, setNum) pg = urllib2.urlopen(searchURL) retDict = json.loads(pg.read) for i in range(len(retDict[/'items/'])): try: currItem = retDict[/'items/'][i] if currItem[/'product/'][/'condition/'] == /'new/': newFlag = 1 else: newFlag = 0 listOfInv = currItem[/'product/'][/'inventories/'] for item in listOfInv: sellingPrice = item[/'price/'] if sellingPrice > origPrc * 0.5: #❶ 过滤掉不完整的套装 print /"%dt%dt%dt%ft%f/" % (yr,numPce,newFlag,origPrc, sellingPrice) retX.append([yr, numPce, newFlag, origPrc]) retY.append(sellingPrice) except: print /'problem with item %d/' % idef setDataCollect(retX, retY): searchForSet(retX, retY, 8288, 2006, 800, 49.99) searchForSet(retX, retY, 10030, 2002, 3096, 269.99) searchForSet(retX, retY, 10179, 2007, 5195, 499.99) searchForSet(retX, retY, 10181, 2007, 3428, 199.99) searchForSet(retX, retY, 10189, 2008, 5922, 299.99) searchForSet(retX, retY, 10196, 2009, 3263, 249.99) 上述程序清单中的第一个函数是searchForSet，它调用Google购物API并保证数据抽取的正确性。这里需要导入新的模块：time.sleep、json和urllib2。但是一开始要休眠10秒钟，这是为了防止短时间内有过多的API调用。接下来，我们拼接查询的URL字符串，添加API的key和待查询的套装信息，打开和解析操作通过json.loads方法实现。完成后我们将得到一部字典，下面需要做的是从中找出价格和其他信息。
部分返回结果的是一个产品的数组，我们将在这些产品上循环迭代，判断该产品是否是新产品并抽取它的价格。我们知道，乐高套装由很多小插件组成，有的二手套装很可能会缺失其中一两件。也就是说，卖家只出售套装的若干部件（不完整）。因为这种不完整的套装也会通过检索结果返回，所以我们需要将这些信息过滤掉（可以统计描述中的关键词或者是用贝叶斯方法来判断）。我在这里仅使用了一个简单的启发式方法：如果一个套装的价格比原始价格低一半以上，则认为该套装不完整。程序清单8-5在代码❶处过滤掉了这些套装，解析成功后的套装将在屏幕上显示出来并保存在list对象retX和retY中。
程序清单8-5的最后一个函数是 setDataCollect，它负责多次调用searchForSet。函数searchForSet的其他参数是从www.brickset.com收集来的，它们也一并输出到文件中。
下面看一下执行结果，添加程序清单8-5中的代码之后保存regression.py，在Python提示符下输入如下命令：
>>> lgX = ; lgY = >>> regression.setDataCollect(lgX, lgY)2006 800 1 49.990000 549.9900002006 800 1 49.990000 759.0500002006 800 1 49.990000 316.9900002002 3096 1 269.990000 499.9900002002 3096 1 269.990000 289.990000 . .2009 3263 0 249.990000 524.9900002009 3263 1 249.990000 672.0000002009 3263 1 249.990000 580.000000 检查一下lgX和lgY以确认一下它们非空。下节我们将使用这些数据来构建回归方程并预测乐高玩具套装的售价。
8.6.2 训练算法：建立模型
上一节从网上收集到了一些真实的数据，下面将为这些数据构建一个模型。构建出的模型可以对售价做出预测，并帮助我们理解现有数据。看一下Python是如何完成这些工作的。
首先需要添加对应常数项的特征X0（X0=1），为此创建一个全1的矩阵：
>>> shape(lgX)(58, 4)>>> lgX1=mat(ones((58,5)))
接下来，将原数据矩阵lgX复制到新数据矩阵lgX1的第1到第5列：
>>> lgX1[:,1:5]=mat(lgX) 确认一下数据复制的正确性：
>>> lgX[0][2006.0, 800.0, 0.0, 49.990000000000002]>>> lgX1[0]matrix([[ 1.00000000e+00, 2.00600000e+03, 8.00000000e+02, 0.00000000e+00, 4.99900000e+01]]) 很显然，后者除了在第0列加入1之外其他数据都一样。最后在这个新数据集上进行回归处理：
>>> ws=regression.standRegres(lgX1,lgY)>>> wsmatrix([[ 5.53199701e+04], [ -2.75928219e+01], [ -2.68392234e-02], [ -1.12208481e+01], [ 2.57604055e+00]]) 检查一下结果，看看模型是否有效：
>>> lgX1[0]*wsmatrix([[ 76.07418853]])>>> lgX1[-1]*wsmatrix([[ 431.17797672]])>>> lgX1[43]*wsmatrix([[ 516.20733105]]) 可以看到模型有效。下面看看具体的模型。该模型认为套装的售价应该采用如下公式计算：
$55319.97-27.59*Year-0.00268*NumPieces-11.22*NewOrUsed+2.57*original price 这个模型的预测效果非常好，但模型本身并不能令人满意。它对于数据拟合得很好，但看上去没有什么道理。从公式看，套装里零部件越多售价反而会越低。另外，该公式对新套装也有一定的惩罚。
下面使用缩减法中一种，即岭回归再进行一次实验。前面讨论过如何对系数进行缩减，但这次将会看到如何用缩减法确定最佳回归系数。打开regression.py并输入下面的代码。
程序清单8-6 交叉验证测试岭回归
def crossValidation(xArr,yArr,numVal=10): m = len(yArr) indexList = range(m) errorMat = zeros((numVal,30)) for i in range(numVal): #❶（以下两行）创建训练集和测试集容器 trainX=; trainY= testX = ; testY = random.shuffle(indexList) for j in range(m): if j < m*0.9: #❷（以下五行）数据分为训练集和测试集 trainX.append(xArr[indexList[j]]) trainY.append(yArr[indexList[j]]) else: testX.append(xArr[indexList[j]]) testY.append(yArr[indexList[j]]) wMat = ridgeTest(trainX,trainY) for k in range(30): #❸（以下三行）用训练时的参数将测试数据标准化 matTestX = mat(testX); matTrainX=mat(trainX) meanTrain = mean(matTrainX,0) varTrain = var(matTrainX,0) matTestX = (matTestX-meanTrain)/varTrain yEst = matTestX * mat(wMat[k,:]).T + mean(trainY) errorMat[i,k]=rssError(yEst.T.A,array(testY))meanErrors = mean(errorMat,0)minMean = float(min(meanErrors))bestWeights = wMat[nonzero(meanErrors==minMean)]xMat = mat(xArr); yMat=mat(yArr).TmeanX = mean(xMat,0); varX = var(xMat,0)#❹（以下三行）数据还原unReg = bestWeights/varXprint /"the best model from Ridge Regression is:n/",unRegprint /"with constant term: /",-1*sum(multiply(meanX,unReg)) + mean(yMat) 上述程序清单中的函数crossValidation有三个参数，前两个参数lgX和lgY存有数据集中的X和Y值的list对象，默认lgX和lgY具有相同的长度。第三个参数numVal是算法中交叉验证的次数，如果该值没有指定，就取默认值10。函数crossValidation首先计算数据点的个数 m。创建好了训练集和测试集容器❶，之后创建了一个list并使用Numpy提供的random.shuffle函数对其中的元素进行混洗（shuffle），因此可以实现训练集或测试集数据点的随机选取。❷处将数据集的90%分割成训练集，其余10%为测试集，并将二者分别放入对应容器中。
一旦对数据点进行混洗之后，就建立一个新的矩阵wMat来保存岭回归中的所有回归系数。我们还记得在8.4.1节中，函数ridgeTest使用30个不同的λ值创建了30组不同的回归系数。接下来我们也在上述测试集上用30组回归系数来循环测试回归效果。岭回归需要使用标准化后的数据，因此测试数据也需要用与测试集相同的参数来执行标准化。在❸处用函数rssError计算误差，并将结果保存在errorMat中。
在所有交叉验证完成后，errorMat保存了ridgeTest里每个λ对应的多个误差值。为了将得出的回归系数与standRegres作对比，需要计算这些误差估计值的均值1。有一点值得注意：岭回归使用了数据标准化，而standRegres则没有，因此为了将上述比较可视化还需将数据还原。在❹处对数据做了还原并将最终结果展示。
1. 此处为10折，所以每个λ应该对应10个误差。应选取使误差的均值最低的回归系数。----译者注
来看一下整体的运行效果，在regression.py中输入程序清单8-6中的代码并保存，然后执行如下命令：
>>> regression.crossValidation(lgX,lgY,10)The best model from Ridge Regression is:[[ -2.96472902e+01 -1.34476433e-03 -3.38454756e+01 2.44420117e+00]]with constant term: 59389.2069537 为了便于与常规的最小二乘法进行比较，下面给出当前的价格公式：
$59389.21-29.64*Year-0.00134*NumPieces-33.85*NewOrUsed+2.44*original price. 可以看到，该结果与最小二乘法没有太大差异。我们本期望找到一个更易于理解的模型，显然没有达到预期效果。为了达到这一点，我们来看一下在缩减过程中回归系数是如何变化的，输入下面的命令：
>>> regression.ridgeTest(lgX,lgY)array([[ -1.45288906e+02, -8.39360442e+03, -3.28682450e+00, 4.42362406e+04],[ -1.46649725e+02, -1.89952152e+03, -2.80638599e+00,4.27891633e+04], . .[ -4.91045279e-06, 5.01149871e-08, 2.40728171e-05,8.14042912e-07]]) 这些系数是经过不同程度的缩减得到的。首先看第1行，第4项比第2项的系数大5倍，比第1项大57倍。这样看来，如果只能选择一个特征来做预测的话，我们应该选择第4个特征，也就是原始价格。如果可以选择2个特征的话，应该选择第4个和第2个特征。
这种分析方法使得我们可以挖掘大量数据的内在规律。在仅有4个特征时，该方法的效果也许并不明显；但如果有100个以上的特征，该方法就会变得十分有效：它可以指出哪些特征是关键的，而哪些特征是不重要的。