★ FEATURED ARTICLE 简介「钻石价格分析-基于R语言」是一份面向R语言初学者与数据分析入门者的实操型文档围绕ggplot2内置的diamonds数据集展开尝试解决从数据清洗到回归建模预测价格这一完整分析链路中的学习需求可用作统计学课程作业、课程设计或数据分析练手项目的参考范例。资源包内仅含1个docx文档压缩后约2.4MB内容集中于单一主题便于一次性通读并对照复现。文档从数据探索切入依次覆盖缺失值与重复值检查、用箱线图识别异常值并剔除、借助价格直方图与核密度曲线判断右偏分布再延伸到把cut、color、clarity等定性变量编码为factor、构造新变量、做简单线性回归与多元回归的变量筛选同时配有各属性与价格关系的小提琴图、透明度分组对比等可视化结果与解读最后给出价格走势预测的思路。已有284人学习适合希望在真实数据集上跑通探索、建模、解释全流程的读者。1. 54000 颗钻石摆在一起先看清这份 R 语言分析的真实难点钻石价格分析这个题目几乎每本 R 语言教材都拿 ggplot2 自带的 diamonds 数据集练一遍——54000 多条记录、10 个变量看着像个玩具数据集。真把回归跑完就会发现不对劲carat 和 price 的关系明显不是直线cut 的系数跑出来是负的R² 顶到 0.88 就上不去了。这份资源走的是完整链路从缺失值、重复行、异常值清洗到描述性统计、分组可视化再到简单线性回归、带交互项的多元回归、非参数显著性检验和回归诊断。它适合已经能写几行 R、但没完整跑过一次「清洗—可视化—建模—诊断」的人也适合想拿一个真实分布右偏、长尾、混着定性变量练回归的人。下面按我实际复现的顺序拆代码能直接抄。2. 数据装载与清洗哪些行不该进回归模型这一步决定后面所有结论的上限。diamonds 数据集本身是干净的坑不在缺失值而在因子类型和物理上不可能存在的尺寸记录。前者会让你的回归系数含义整个跑偏后者会把 x、y、z 的系数拉到一个没法解释的位置。2.1 装载与类型确认ordered factor 是第一个陷阱# diamonds 挂在 ggplot2 包里先确认包在不在不要无脑重装 if (!requireNamespace(ggplot2, quietly TRUE)) { install.packages(ggplot2) } library(ggplot2) # 显式指定 package避免工作区里有同名对象时拿到旧数据 data(diamonds, package ggplot2) dim(diamonds) # 行数、列数 str(diamonds) # 变量类型 因子水平 sapply(diamonds, class) # 只看 class一眼扫完 10 个变量逻辑说明data()的第二个参数package是关键。R 的data()会先在工作区找、再去已加载的包里找如果之前load()过一个叫 diamonds 的 rda直接data(diamonds)拿到的东西可能不是你以为的那份。str()要看两件事——数值列的 num/int 类型对不对因子列的水平数和顺序对不对。参数说明quietly TRUE让requireNamespace不打印加载信息脚本source()的时候输出干净install.packages()只在包缺失时触发不会每次都跑。跑完str()你会看到cut、color、clarity都是Ord.factor。这是第一个必须处理的地方。R 对有序因子ordered factor默认使用多项式对比contr.poly也就是把 5 个切工等级拆成线性、二次、三次、四次趋势项回归结果里出现的是cut.L、cut.Q、cut.C这种名字。你想要的「Ideal 相对 Fair 贵多少」根本不在输出里。所以建模前必须显式重建dm - diamonds # 重建为普通因子levels 按业务从差到好手工排 dm$cut - factor(dm$cut, levels c(Fair, Good, Very Good, Premium, Ideal), ordered FALSE) dm$color - factor(dm$color, levels c(J, I, H, G, F, E, D), ordered FALSE) dm$clarity - factor(dm$clarity, levels c(I1, SI2, SI1, VS2, VS1, VVS2, VVS1, IF), ordered FALSE) # 确认对比方式变成 treatment 对比 contrasts(dm$cut)逻辑说明levels的第一个元素就是回归里的基准组。cut 的基准是 Fair最差切工color 的基准是 J最差颜色clarity 的基准是 I1最差净度——三个都选最差档做基准好处是所有系数都是正数读起来就是「比最差档贵多少」方向和业务直觉一致。参数说明ordered FALSE不能省。哪怕你levels排得再对只要保留 ordered 属性lm()照样走多项式对比。我见过不少人的模型里cut.L显著、cut.Q不显著然后对着输出发半天呆就是这么来的。2.2 缺失值与重复行先量化再决定删不删library(mice) # 缺失值模式矩阵最后一列是每行缺失个数最后一行是每列缺失个数 md.pattern(dm, rotate.names TRUE) # 不装 mice 的话用这两行也够 colSums(is.na(dm)) sum(duplicated(dm))逻辑说明md.pattern()返回的矩阵里1 表示有值、0 表示缺失最后一行统计每个变量的缺失条数最后一列统计每条记录的缺失变量数。diamonds 这个数据集跑出来全是 1说明没有缺失值。duplicated()返回逻辑向量sum()得到重复行数。参数说明rotate.names TRUE把变量名转 90 度10 个变量的表横着排屏幕塞不下转过来才读得清。重复行这里有个判断要做。diamonds 里确实存在属性完全一致的行但它们可能是两颗独立打磨的钻石——克拉、切工、颜色、净度、尺寸都撞上并不是不可能。这份资源的做法是直接去重我的做法是先记数量再决定n_before - nrow(dm) dm - dm[!duplicated(dm), ] cat(去重前:, n_before, 去重后:, nrow(dm), 差异占比:, round((n_before - nrow(dm)) / n_before * 100, 3), %\n)逻辑说明:-号在[ 里配合逻辑向量是标准的「保留非重复行」写法。打出去重前后的行数和差异占比用来判断这个操作要不要写进最终流程。参数说明差异占比在千分之几的量级删掉对结论没影响如果超过 1%就要考虑这些「重复」是不是清洗阶段的副产物比如你把连续变量四舍五入过那就不能删。2.3 异常值判定箱线图只能告诉你「有」不能告诉你「删哪些」# 第一步物理规则过滤三维尺寸必须为正 dm - subset(dm, x 0 y 0 z 0) # 第二步尺寸上限y 最大 58.9mm 这类是典型的小数点录入错误 dm - subset(dm, x 20 y 20 z 20) # 第三步对 price 只做标记先不删 q - quantile(dm$price, c(0.25, 0.75)) iqr - q[2] - q[1] dm$flag_price_out - dm$price (q[1] - 1.5 * iqr) | dm$price (q[2] 1.5 * iqr) table(dm$flag_price_out)逻辑说明这里刻意拆成三步而不是一行subset()。第一步去掉尺寸为 0 的记录——0 毫米的钻石不存在这类行通常是缺测被填成 0 的结果。第二步卡上限原始数据里 y 的最大值到了 58.9 毫米而 carat 最大才 5.01一颗 5 克拉的钻长边通常在 10 毫米量级58.9 是数量级错误。第三步对价格只打标记不删除因为 price 是右偏长尾大克拉高价钻恰恰是这个数据集最有信息量的部分一刀切会把模型往中间压缩。参数说明1.5 * iqr是箱线图触须的标准倍率对应「离群但不算离谱」的边界想更严格就换3 * iqr。flag_price_out是一个新列后面做敏感性分析时用它分层重跑一遍模型看系数有没有大变化比直接删掉稳妥得多。对比一下原文用三倍标准差剔除异常值的做法。在正态假设下 3σ 是合理的但 price 明显不服成正态分布下一章会看到用正态的尺子去量长尾分布等于系统性地砍掉高价样本最后模型会低估大克拉钻石的价格。所以我一般对价格类变量用 IQR 标记 分层验证对物理尺寸用硬规则剔除。2.4 衍生变量与建模底表把清洗流程固化成一次# 体积三个维度的乘积和 carat 高度相关但不等价 dm$volume - dm$x * dm$y * dm$z # 每克拉单价用来观察「越大越贵」是不是线性 dm$price_per_carat - dm$price / dm$carat # 把清洗结果落盘后面所有分析都从这份快照读 saveRDS(dm, diamonds_clean.rds)逻辑说明volume和x、y、z之间存在极强的线性关系建模时这两组变量二选一同时放进去会让系数标准误爆炸。price_per_carat是个诊断变量不是建模变量——把它按 carat 分箱画出来能直接看到每克拉单价在 1 克拉、1.5 克拉附近的变化比看回归系数直观。参数说明saveRDS()存单个对象读回来用readRDS()比write.csv()快且不会丢因子水平。清洗脚本从头到尾跑一遍只要几秒但如果你每次分析都重跑清洗、每次都微调一下过滤阈值前后结论就对不上了——固化快照是保证可复现的最低成本手段。这一步别用attach()。attach(diamonds)之后工作区里凭空多出 10 个变量名price到底是原始的diamonds$price还是你改过的dm$price你自己三天后也说不清这是黑匣子级的事故源。3. 描述性统计与可视化价格右偏、1000 美元峰值和四张必看图清洗完先别急着lm()。这个数据集的形状决定了后面模型怎么设价格右偏到什么程度、峰值在哪、分组之后分布有没有重叠、carat 和 price 到底是直线还是分段。这四件事没搞清楚就上回归系数解释一定出事。3.1 价格分布三件套直方图、核密度与正态曲线叠加library(ggplot2) ggplot(dm, aes(x price)) geom_histogram(aes(y after_stat(density)), bins 60, fill grey75, colour white) geom_density(colour #c0392b, linewidth 0.8) stat_function(fun dnorm, colour steelblue, linetype dashed, args list(mean mean(dm$price), sd sd(dm$price))) labs(x 价格美元, y 密度) theme_minimal()逻辑说明after_stat(density)把直方图的纵轴从计数换成密度这样核密度曲线和正态曲线才能画在同一根轴上比较。三条线叠在一起的目的很明确——灰柱是真实分布红线是平滑后的形状蓝虚线是「如果它服从正态会是什么样」。三者一比右偏有多严重一目了然。参数说明bins 60是对 5 万多条记录比较合适的粒度太少看不出断层太多全是毛刺。linewidth在新版 ggplot2 里替代了旧的size写size会报警告。stat_function的args必须传list()直接写mean ..., sd ...会被当成未知参数。跑完能看到三件事。第一峰值在 1000 美元附近0 到 2000 美元区间堆了绝大多数样本。第二2000、4000、8000 美元附近有明显断层——这是定价习惯造成的整数关口附近的样本密度会突然掉下去。第三正态曲线和核密度曲线在 8000 美元之前几乎完全不重合之后才逐渐贴近横轴。这意味着任何假设 price 正态的检验都不成立后面做显著性检验只能走非参数路线这也是原文选非参数检验的原因。3.2 分组比较小提琴图为什么比箱线图更适合 cutggplot(dm, aes(x cut, y price, fill cut)) geom_violin(trim FALSE, scale width, alpha 0.7) geom_boxplot(width 0.12, outlier.shape NA, fill white) scale_y_log10() theme_minimal() theme(legend.position none)逻辑说明小提琴图给的是完整的密度形状箱线图给的是中位数和四分位。两张图叠在一起既能看出「Ideal 组是不是整体更贵」也能看出「Ideal 组内部是不是有两坨」。这个数据集里 Ideal 组确实存在一批低价样本只看中位数会被平均值带跑。参数说明trim FALSE保留尾部不然小提琴的两端会被截断看不出长尾scale width让每组小提琴等宽默认是按样本量缩放Fair 组样本少会缩成一条缝视觉上造成误导scale_y_log10()是必须的price 右偏线性轴上所有箱子都压在底部箱线图里outlier.shape NA是因为小提琴已经表达了尾部再画一遍离群点会糊成一团。同样套路换color和clarity再跑两遍。颜色从 J 到 D价格中位数基本单调上升这个方向符合直觉。净度那边就有意思了——IF、VVS1、VVS2 价格最高但 I1 的中位数并不是最低的SI2 反而比 I1 低。原文对这个现象的解释是「净度到达优质状态后区分钻石价格的维度更多」说白了就是净度不是唯一变量低净度区间里 carat 的权重更大。3.3 克拉与价格线性坐标会骗人换对数轴再看ggplot(dm, aes(x carat, y price)) geom_point(alpha 0.15, size 0.6, colour grey30) geom_smooth(method lm, formula y ~ x, colour red, linewidth 0.8) scale_x_log10() scale_y_log10() labs(x 克拉对数轴, y 价格对数轴) theme_minimal()逻辑说明线性坐标下 5 万多个点全挤在左下角什么都看不出来。双对数轴上能看到三层结构——整体趋势是直线说明 carat 和 price 之间接近幂律关系在 0.5、1.0、1.5 克拉处有垂直的「墙」说明整数克拉有溢价定价不是连续的大克拉区域点的离散度明显变大说明高价钻的定价更依赖其他因素。参数说明alpha 0.15是 5 万个点的合理透明度再高就成实心块了size 0.6让单点足够小避免重叠成色块formula y ~ x显式写给geom_smooth新版 ggplot2 不写会警告。这张图直接决定了回归怎么设。如果只用price ~ carat残差图一定有明显的弧形因为真实关系是对数级的。后面第 4 章会看到单变量模型的 R² 已经不低但残差结构不对这就是原因。3.4 相关性矩阵与共线性初判num_vars - c(carat, depth, table, x, y, z, price) round(cor(dm[num_vars], use complete.obs), 3) # 只看尺寸和 carat 的共线性 library(car) vif(lm(price ~ carat depth table x y z, data dm))逻辑说明cor()给的是两两线性相关vif()给的是「某个变量被其他变量解释的程度」后者对多元回归更有意义。x、y、z 三者之间的相关系数普遍在 0.95 以上VIF 会直接冲到几十。参数说明use complete.obs是保险写法虽然这个数据集没缺失但脚本换个数据源就跑加一行省事。VIF 的经验阈值是 10超过就说明共线严重系数符号可能反掉。结果是 x、y、z 的 VIF 全部超标。这不奇怪一颗钻的长宽深本来就同步变化。处理方式有两种要么只留volume一个综合维度要么留carat然后把 x、y、z 全砍掉。我一般选后者因为 carat 是行业里真正的交易单位系数解释起来有业务含义而「体积每增加 1 立方毫米贵多少钱」没人听得懂。4. 回归建模从 carat 单变量系数 7760 到带交互项的多元模型到这一步才开始写lm()。顺序不能乱先单变量看方向和量级再全变量看整体解释力最后加交互项看边际效应是不是随分组变化。每一步都要留残差图看结构有没有被吃掉。4.1 简单线性回归7760 美元/克拉这个数怎么读fit1 - lm(price ~ carat, data dm) summary(fit1) # 系数置信区间比只看 p 值信息量大 confint(fit1, level 0.95) # 四张诊断图一次出 par(mfrow c(2, 2)) plot(fit1) par(mfrow c(1, 1))逻辑说明summary()里重点看三行——carat 的 Estimate、它的 Pr(|t|)、还有最下面的 R²。这份资源里跑出来 carat 的系数约 7760也就是「克拉数每增加 1价格平均增加约 7760 美元」。截距是负的但在业务上没有意义因为不存在 0 克拉的钻石别硬解释。参数说明confint()默认 95% 置信区间可以改成level 0.99。par(mfrow c(2,2))把四张图排成 2×2跑完记得用par(mfrow c(1,1))复位不然后面所有图都会被挤到一个小格子里——这是新手最常见的「图怎么画不出来」问题。四张诊断图里先看 Residuals vs Fitted。如果点云呈喇叭口左窄右宽说明存在异方差残差的方差随拟合值增大这时候系数的标准误会算小p 值会虚高。再看 Normal Q-Q点偏离对角线两端翘起说明残差不正态——price 右偏残差不可能正态。这两条都在提示同一件事直接对 price 建线性模型假设层面就不成立后面才会去考虑对数变换和交互项。4.2 哑变量处理基准组选谁系数含义就变了# 确认三个因子的对比方式 contrasts(dm$cut) contrasts(dm$color) contrasts(dm$clarity) # 显式指定基准组写成代码里比靠 levels 顺序更稳 dm$cut - relevel(dm$cut, ref Fair) dm$color - relevel(dm$color, ref J) dm$clarity - relevel(dm$clarity, ref I1)逻辑说明contrasts()打印出来的是对比矩阵行是水平名列是哑变量名第一行全 0 的那个水平就是基准组。relevel()的作用是把指定水平挪到第一位也就是改基准。参数说明ref参数只接受已经存在于 levels 里的字符串写错会直接报错。这一步我建议在建模前单独跑一次relevel()并写死而不是依赖上游factor()里 levels 的排列顺序——上游改一行下游所有系数含义全变这种坑排查起来极其费时间。为什么基准组这么重要因为系数是「相对基准组的差值」。如果 color 的基准是 J最差那 D 的系数就是正的一大截「颜色越好越贵」一眼就懂。如果按默认字母序基准变成 D最好所有系数都变成负数输出看起来像「颜色越差越贵」实际上只是参照系反了。这个符号问题在下一章的避坑清单里还会再出现一次——它是这份数据集中最容易翻车的地方。4.3 无交互项多元回归R²0.8798 说明什么fit2 - lm(price ~ carat cut color clarity depth table x y z, data dm) summary(fit2)逻辑说明这是标准的加性模型每个变量对 price 的边际效应是常数不随其他变量变化。这份资源里 R² ≈ 0.8798调整 R² 略低一点。也就是说还有 12% 的价格波动没被解释。参数说明lm()里公式写就是加性不写交互。注意这里同时放了carat和x、y、z共线性严重你会看到x的系数可能是负的——这不是「钻石越长越便宜」而是「在 carat 和其他尺寸不变的前提下x 单独变大反而降价」业务上没意义。正确做法是砍掉 x、y、z或者用上一章算的 volume 替代。depth和table这两个变量在这个模型里通常不显著符号还可能反。depth 是深度百分比定义是2z/(xy)它和 x、y、z 是函数关系放进来纯属添乱。我的处理是先在模型里保留它们看 p 值不显著就从最终模型里去掉别舍不得——每多一个不显著的变量调整 R² 就掉一点。R² 0.8798 这个数怎么评价单看绝对值不算低但结合残差图看就不够。残差在低价区和高价区呈现不同的离散度说明模型对便宜钻和贵钻的拟合能力不一样。这就是这份资源往下走、加交互项的原因。4.4 交互项与 AIC 逐步回归三个交互项全进模型之后# carat 与三个定性变量的交互主效应会被自动保留 fit3 - lm(price ~ carat * (cut color clarity) depth table x y z, data dm) summary(fit3)$r.squared AIC(fit3) # 向后逐步按 AIC 删变量 fit4 - step(fit3, direction backward, trace FALSE) summary(fit4) AIC(fit3, fit4)逻辑说明carat * (cut color clarity)是展开写法等价于 carat、cut、color、clarity 四个主效应加上 carat:cut、carat:color、carat:clarity 三个交互项。交互项的含义是「克拉数对价格的影响在不同切工/颜色/净度下不一样」——这正是加性模型缺的东西。参数说明step()的direction可选backward、forward、both向后是从全模型出发逐步删trace FALSE关掉每一步的打印不然几万行输出刷屏向前选择时必须给scope参数指定候选变量范围。5 万多行数据、几十个参数step()会反复拟合几十次模型总耗时也就一两秒行数再上一个量级就要考虑换方案了。原文的结论是「所有变量都被选入模型」这个结果本身值得琢磨——AIC 逐步回归在 5 万样本上非常保守因为 AIC 的惩罚项是2k而样本量一大任何一个微小改进的似然增益都超过这个惩罚。所以别把「全都留下了」理解成「所有变量都重要」它只说明在这么大样本下删掉任何一个都会让 AIC 变差一点点。交互项进来之后 R² 会往上走一点但更值得看的是交互项的系数。如果carat:cutIdeal是负的意思就是「同样增加一克拉Ideal 切工的钻石涨价比 Fair 少」——这在业务上是说得通的因为小克拉钻更强调切工大克拉钻更强调保重。这类结论才是交互项真正的价值光看 R² 涨了零点零几没什么感觉。5. 踩坑排查R 语言钻石回归里最容易翻车的五个地方下面五条都是我在这份数据上真实撞过的。每条按「跑出来什么样 → 为什么 → 怎么改」写遇到了直接对号入座。5.1 cut 的系数是负的是不是数据错了现象模型里cutIdeal的系数是负数而且显著可你单独按 cut 分组算平均价格Ideal 组一点也不便宜。两个结论打架第一反应是数据有问题。原因有两层。第一层是参照系——如果 cut 还是 ordered factor 没改过来你看到的根本不是cutIdeal而是cut.L、cut.Q这些趋势项符号含义完全不同。第二层是混杂Ideal 切工大量出现在小克拉钻上切工好能让小钻显大而 Fair 切工往往出现在大克拉钻上保重优先。当你控制住 carat 之后Ideal 在同等克拉下反而便宜这是典型的辛普森悖论式现象不是数据错。解决先contrasts(dm$cut)确认对比方式再tapply(dm$price, dm$cut, median)看分组中位数两个都正常说明是混杂不是错误。解释系数时必须带限定语——「在控制克拉重量、颜色、净度等其他变量不变的前提下」。想进一步验证就加carat:cut交互项看不同克拉区间里 cut 的效应是不是反过来的。5.2 对 price 取对数后回归系数不能直接读成美元现象为了处理右偏把模型改成lm(log(price) ~ log(carat) ...)跑出来 carat 的系数是 1.8然后有人直接写「克拉数每增加 1价格增加 1.8 美元」。原因双对数模型里系数是弹性1.8 的意思是「carat 增加 1%price 平均增加 1.8%」跟美元没有半点关系。如果是半对数模型只对 y 取 log系数 β 对应的是「x 增加 1 单位y 变化(exp(β)-1)×100%」β 小的时候近似等于 β×100%β 大的时候这个近似误差很难看。解决写个转换函数固定下来别每次心算。# 半对数模型系数转百分比 pct_change - function(beta) (exp(beta) - 1) * 100 # 双对数模型的弹性直接读但要写成百分比表述 fit_log - lm(log(price) ~ log(carat) cut color clarity, data dm) summary(fit_log)$coefficients[log(carat), Estimate]参数说明pct_change()传入的是summary()里的 Estimate 列不是 t 值。双对数模型的弹性解释可以直接用系数值乘以 100 变成百分比不需要再套exp()。5.3 逐步回归结果每次跑都不一样现象同一份数据同事那边step()选出来的变量集跟你不一样或者你自己隔一天重跑变量少了一个。原因step()本身不吃随机数结果不稳定通常来自三个地方。一是起始模型不同——你从全模型出发他从只含主效应的模型出发贪心搜索的路径当然不同。二是样本量变了——上游某一步na.omit()或者过滤条件改过参与拟合的行数不一样AIC 没法比。三是工作区里有残留对象lm()实际用的是旧的dm。解决把清洗结果saveRDS()落盘建模脚本第一行dm - readRDS(diamonds_clean.rds)彻底切断对工作区状态的依赖。起始模型写成固定公式别用~ .这种简写——~ .会把数据框里所有列都拉进来包括你后面加的flag_price_out、price_per_carat这些诊断列直接造成数据泄漏。每次跑完把AIC(fit3, fit4)和变量列表一起记在脚本注释里。5.4 交互项全加进去之后报错或出现 NA 系数现象summary()里冒出几个NA的系数或者直接报对比相关的错误提示某个因子水平有问题。原因某个因子水平的组合样本量为 0。比如你按某个条件过滤之后cut Fair且color D的记录一条都不剩交互项在这个格子上就没有数据设计矩阵秩亏R 只能给 NA。另一个常见原因是把depth又按定义2*z/(xy)算了一遍塞进模型和已有的 x、y、z 构成精确共线。解决建模前先查空格子。# 因子组合的样本量出现 0 的地方就是问题所在 table(dm$cut, dm$color) # 过滤后重新降水平防止空水平残留 dm$cut - droplevels(dm$cut) dm$color - droplevels(dm$color)参数说明droplevels()只删掉样本量为 0 的水平不影响其他水平。table()输出里出现 0 的格子要么合并水平要么在模型里去掉这个交互项没有第三种办法。5.5 拿训练集 R² 当预测精度现象模型 R² 0.98看着很美拿几颗新钻石的报价去预测误差两三千美元。原因R² 是样本内指标变量越多它越高跟预测能力不是一回事。尤其是这个数据集里 carat 和 x、y、z 高度共线模型相当于用四个变量反复描述同一件事样本内拟合自然好。另外异常值对 R² 的影响很大几个极端点就能把它拉上去或拽下来。解决划分训练集和测试集在测试集上算 RMSE 和 MAE这两个才是能拿出去说的数。同时看残差的分布如果残差在不同价格区间上系统性偏正或偏负说明模型有结构性偏差加变量也救不回来。具体做法放下一章。6. 模型验证与一个可复用技巧留出集 RMSE 和 log 反变换的修正模型跑完到能交付中间还差两步在没见过的数据上量一次误差以及把对数模型的预测值正确换回美元。第二步最容易被忽略也最容易造成系统性低估。6.1 留出集验证RMSE 怎么算才有意义set.seed(2024) n - nrow(dm) idx - sample.int(n, size floor(0.7 * n)) train - dm[idx, ] test - dm[-idx, ] fit_tr - lm(log(price) ~ log(carat) cut color clarity, data train) pred_log - predict(fit_tr, newdata test) rmse - function(actual, pred) sqrt(mean((actual - pred)^2)) mae - function(actual, pred) mean(abs(actual - pred))逻辑说明sample.int()抽行号而不是直接抽数据框后面用-idx取补集时不会出错。70/30 是个常规比例5 万多条记录里拿 1.6 万条做测试足够稳定。predict()的newdata必须包含模型里用到的所有变量缺一个就报错——这是好事能提前发现测试集列名不一致的问题。参数说明set.seed(2024)固定划分换机器重跑得到同一个测试集否则两次结果没法比。RMSE 的单位和响应变量一致对 log 模型来说是「对数价格」的单位不能直接读成美元得先反变换。6.2 log 反变换直接 exp 会系统性低估朴素做法是exp(pred_log)但这一步有问题。对数变换是非线性的E[exp(Y)] ≠ exp(E[Y])直接取指数得到的是中位数预测而不是均值预测整体偏低。修正办法叫 smearing estimate用训练集残差构造一个缩放因子# 朴素反变换 pred_naive - exp(pred_log) # smearing 修正因子训练集残差取 exp 后的均值 smear - mean(exp(residuals(fit_tr))) pred_smear - pred_naive * smear # 两种反变换在测试集上的误差对比 rbind( naive c(RMSE rmse(test$price, pred_naive), MAE mae(test$price, pred_naive)), smear c(RMSE rmse(test$price, pred_smear), MAE mae(test$price, pred_smear)) )逻辑说明residuals(fit_tr)是训练集上的对数残差exp()之后取均值就是那个缩放因子。它衡量的是「对数尺度上的误差在原始尺度上会被放大多少」价格分布越分散这个因子离 1 越远。参数说明smear必须用训练集残差算用测试集算等于偷看了答案。这个因子通常在 1.05 到 1.3 之间价格右偏越厉害越大。如果算出来接近 1说明对数变换没什么必要直接建线性模型更省事。对比表我一般这么看指标含义判断标准训练集 R²样本内解释力只用于横向比模型不能当精度测试集 RMSE预测误差的绝对量级和 price 中位数比超过 30% 就要警惕测试集 MAE平均绝对误差比 RMSE 抗极端值两个一起看naive 与 smear 的差值反变换偏差差值大说明分布右偏严重必须修正6.3 从那以后我每次都先跑一遍反变换对比吃了两次亏之后我养成一个习惯任何带 log 变换的模型在写结论之前先把 naive 和 smear 两列误差并排打出来。第一次踩坑是给一个采购场景做报价预估模型 R² 很好看交付之后对方反馈「你们的本文还有配套的精品资源点击获取
阅读完成 · 觉得有帮助？
