- 在线时间
- 514 小时
- 最后登录
- 2023-12-1
- 注册时间
- 2018-7-17
- 听众数
- 15
- 收听数
- 0
- 能力
- 0 分
- 体力
- 40219 点
- 威望
- 0 点
- 阅读权限
- 255
- 积分
- 12777
- 相册
- 0
- 日志
- 0
- 记录
- 0
- 帖子
- 1419
- 主题
- 1178
- 精华
- 0
- 分享
- 0
- 好友
- 15
TA的每日心情 | 开心 2023-7-31 10:17 |
|---|
签到天数: 198 天 [LV.7]常住居民III
- 自我介绍
- 数学中国浅夏
 |
【R语言】回归分析案例:北京市商品房价格影响因素分析% n: z5 O9 M6 m+ _
这一案例是王汉生老师《应用商务统计分析》方差分析章节的案例,主要对离散型变量进行了处理。
4 G. t" N- @7 R5 }3 Z n1 V这里将连续型变量也加进来,进行协方差分析,建立完整的模型。 首先对房价进行对数变换,解决异方差问题: 0 r" f6 T1 c& X: i3 ?4 l; V
行描述性统计分析,各连续型变量之间的相关关系如下:
% x$ m7 t1 y! i/ m$ w" N' g1 ~![]()
+ z' _. U& s8 b: T名义变量的EDA一般做箱型图。 模型按照全模型-变量处理(分箱等)-变量选择-回归诊断等步骤建立。
: G0 w o" J% C. H# r![]()
8 q/ \$ x) v: n% a & c* F6 U9 Y8 ]# @4 t
最终模型残差图:
, X' \/ a9 Z7 f$ v & G" N& b7 k8 `4 b" w/ q
通过模型分析结果可知,影响北京市商品房平均销售价格的主要因素有:" E H% y P4 O( p* M+ a6 S6 e t
属性变量:所在辖区、所在环线、物业类别、装修状况、容积率大小(新引入);连续变量:绿化率、停车位住户比+ B2 P- a' l8 Z* o1 S
属性变量的具体影响在此处分析略去。
^/ b; D% ~) `( I% _' a1 d+ h连续型变量的影响主要为:
% s. Y4 q# F. G; O+ E" b9 S 绿化率:绿化率的影响十分显著,由系数估计值为正,说明对房价有正向影响,绿化率越高的楼盘房价越高;$ f# n2 G9 [; t% h# A6 ^( |
停车位住户比:有较显著的影响,停车位住户比越高,价格越高;
% A# @1 X0 a4 K0 Q同时,原本为连续型变量的容积率经过离散化变为属性变量后:$ R2 Z3 X# }9 n
容积率大小:容积率分组有较显著的影响,高容积率的小区商品房价格更贵;
6 V6 a J: V/ K2 z# P6 s7 j 容积率与环线之间存在着交互效应。' I: V. D# q2 R8 G! `$ I
rm(list=ls()) #清空当前工作空间- a, q0 ~% h+ ^' c
setwd("D:/回归分析")
, b; \9 @9 L F5 Ka=read.csv("real.csv",header=T) #读入csv格式的数据,赋值为a
8 O6 J; ?; s# C6 w3 {: X7 p2 l$ Y' mView(a)
% G" ^. h7 g- Vattach(a)1 R6 h) k% z3 |- f7 x
names(a)
' C" r& b$ z/ \; s8 W, T) A
5 z2 H# c+ n% }& x Z( D% {
S+ a$ r7 a2 D% }+ d2 i$ F, B##描述性统计
$ v0 J8 W9 y0 Q1 z, n
k* R$ I& c5 p/ d# n) o/ S: K* ^
#未做处理的响应变量分布情况
0 _ z) b, s3 P. Z% epar(mfrow=c(1,1))+ S( {4 @- Y1 Z0 _2 C* D: x( M
hist(price)
2 n- Z+ C. P8 l6 ?7 N2 x% c7 csummary(price) #查看响应变量的描述统计量
' m0 {* w& g P j f: h( x#连续型变量描述性统计8 E( _2 t t" T
windows()( F' Q* L& K" m0 z0 D
pairs(a[,c(6:10)]) #所有连续型变量间的散点图; i$ f: u) y$ ]7 @2 X: @
par(mfrow=c(2,2)) + w$ T' B: u/ M. t) g* @7 q
plot(rong,price) #每个连续型因变量与响应变量间的散点图2 v# Q8 j: p! `+ f1 W
plot(lv,price)
$ |& T% K; {/ v; Q* k% \# M4 H( rplot(area,price)
1 l8 @0 ?# b0 ? P( vplot(ratio,price)
, c; D8 K2 x8 s7 Jsummary(a[,c(6:10)]) #查看连续型变量的描述统计量
- J- V5 J4 T ]cor(a[,c(6:10)]) #查看连续型变量的相关系数
* y9 [/ q( {( k#属性变量描述性统计) ?/ ^/ J* r) q* {6 b6 u
windows()- Y( o$ u; k& L3 ^6 R& S
par(mfrow=c(2,3)) & C+ S" H" c: e; g
boxplot(price~dis) #每个属性变量关于响应变量的箱型图. k. w) C6 i7 J$ q8 o5 w; X
boxplot(price~wuye) ) b1 J# T- U( c
boxplot(price~fitment) 2 T _3 L B- ]* J [4 T& |* q
boxplot(price~ring) 9 ?0 p# |- J# a5 B" M3 B
boxplot(price~contype)) R$ I0 H1 c* H3 K- {
; v' u" S4 Q2 `0 O
7 J: Y5 ]7 O# c/ Y; c% D4 ~
- m6 O2 b2 ^) r. u9 [. \5 _$ X8 `9 w& s2 I d3 k1 p3 P
##模型建立$ m, Y3 ^- j: E4 @% \! Y( s1 y7 s
# ?/ }% l3 I( E% ?& ^3 k) G+ `8 l8 k3 {" i4 U$ q. s2 H$ z _" q
#在方差分析模型基础上加入连续型变量' m+ r j3 @: j8 d% F
lm1=lm(price~as.factor(dis)*as.factor(ring)+as.factor(wuye)+as.factor(fitment)+as.factor(contype)+rong+lv+area+ratio)
/ d/ \* X& q: n) U2 s& }anova(lm1) #方差分析% ?% |% x6 Q: ?3 G4 P! [$ u
summary(lm1) #模型参数估计等详细结果
/ |+ A' Z6 A8 v9 V# Q. {! I: A' _2 N7 Fwindows()9 u; M) m: x) ^: T+ u
par(mfrow=c(2,2))
9 | U% S3 r3 ^plot(lm1,which=c(1:4)) #回归诊断做残差图' V0 E# [2 D0 F% ]4 Z: H6 r
9 F$ ~7 J0 l8 M2 j6 ~
, q8 g2 H* J l5 S# q( v# z+ o1 K P8 r {7 r2 x+ y
4 n: ^; b7 v1 ]! |7 P) D8 ?1 a
##变量处理 v) m9 j z6 t: Y
/ f* h, Y" H8 l* D+ o$ Z9 k0 r% n! \8 p' E, h+ s
###对不显著的变量采用分组的方式希望能达到显著的效果
" ]& e5 s; w5 x, X. U/ s##对容积率的处理# l. M) l& ~! b* X% r
windows()
: N' Q+ l0 V( t' d# y4 Fn = 4
3 k, G6 V. a H: U$ qboxplot(price~ceiling(rong/n)) #容积率多分组下的箱型图
9 p8 m# R0 d; \table(ceiling(rong/n)) #容积率各分组下的样本数
: Z' K' W& W0 z, i3 tronggrp=1*(rong>n) #进行二分类
/ U ]/ Z" a6 P% m- s#ronggrp=ceiling(rong/n) ' U" I" W3 J$ v
table(ceiling(ronggrp)) #容积率二分类下的样本数
6 r. N7 y% V( [( q$ T2 M" C% }windows()
; l3 |! e9 Q$ P: n8 R$ ~boxplot(price~ceiling(ronggrp)) #容积率二分类下的房价箱型图0 ]) b% o% M" M$ R0 J
windows()4 U6 y3 j) t. d
par(mfrow=c(1,2))
, ~# Z3 Z' I4 Q1 G$ L, M- Qboxplot(rong~ring) #容积率与环线箱型图
8 v& K' c) x: I0 p0 J' m: Wboxplot(price~ring) #房价与环线箱型图 8 i. N8 K. r* F7 [3 w/ L
#加入容积率分组和容积率分组*所在环线交互因子的模型
* m8 X r& D. b9 c+ llm2=lm(price~as.factor(dis)*as.factor(ring)+as.factor(wuye)+as.factor(fitment)+as.factor(contype)+as.factor(ring)*as.factor(ronggrp)+lv+area+ratio)( |2 @" u* P9 b* l
anova(lm2) #方差分析. y; t: {2 p; `
summary(lm2) #模型参数估计等详细结果
1 H) Y, y/ ?$ bwindows()
# ]( Z# Y$ z) b" Vpar(mfrow=c(2,2))% ^9 n( b* }, _( E
plot(lm1,which=c(1:4)) #回归诊断
( d4 a5 D l; P, q$ k& W
& w9 G3 r/ _. w# A: }) ^3 N: ~7 Y# M* R& T3 w+ w
##对小区面积的处理
$ j) q7 o; d1 c* Isummary(area)7 N$ J3 M% b% z1 o- J/ K, h, x
plot(area,price)
4 N9 A) a; c$ swindows(). ? a$ }! X( [- w8 W7 w( z
n = 150000; }& N% d3 K, x: b: c
boxplot(price~ceiling(area/n)) . m1 ^9 P5 i+ i' b1 j5 \
table(ceiling(area/n))
" s, \) k1 g0 jareagrp=1*(area>n)
0 `1 b) j2 Z! p9 s3 X' P) Otable(ceiling(areagrp))8 a8 D' y5 w" K* V. L. ]
boxplot(price~ceiling(areagrp))
- _9 P# O2 N0 @# E9 c7 X#加入小区面积分组的模型
* i3 n: |% n% i3 Z# Xlm3=lm(price~as.factor(dis)*as.factor(ring)+as.factor(wuye)+as.factor(fitment)+as.factor(contype)+as.factor(ring)*as.factor(ronggrp)+lv+as.factor(areagrp)+ratio)# o+ ^5 L" j8 Q4 ]0 T# q$ h( n. h
anova(lm3) #方差分析; u. P9 `1 A8 d/ E
summary(lm3) #模型参数估计等详细结果7 L0 T+ f% h0 Y: Z9 l
windows()# N. u0 g( K! P8 [5 O; f4 R; g& u
par(mfrow=c(2,2))
1 Q2 B% H3 q; l/ A; g7 i# U! Rplot(lm3,which=c(1:4)) #回归诊断
7 F) n. H3 f' _/ G. T
X. G. K% { e, O6 u! y& c' A- B
##变量选择6 `# C" X" p$ D/ W0 a
4 A! V1 E7 i- R7 C- `' H9 q% U7 o9 S# `4 [" c+ B* t0 a" T# L
##AIC准则下的变量选择5 k: k5 G p8 R: v7 W6 G" b: [
lm4.aic=step(lm3,trace=F) #根据AIC准则选出最优模型,并赋值给lm.aic
2 A& t: F" | ]" E, r2 @summary(lm4.aic) #给出模型lm.aic中系数估计值、P值等细节
; Q- X0 K" q1 N8 x% C1 N' ?, [* S##BIC准则下的变量选择+ g) H2 M8 ^( a6 }1 u0 e2 w
lm5.bic=step(lm3,k=log(length(a[,1])),trace=F) #根据BIC准则选出最优模型,并赋值给lm.bic
. A+ G* i% `) o! m8 X; b# asummary(lm5.bic) #给出模型lm.bic中系数估计值、P值等细节
9 Z0 y* Y* o" ]2 X5 ^: W! Y t
/ `' G# H9 t8 C. Q, A/ H2 v W6 o
. _& D' \% ~8 u3 |#选用AIC准则下的模型进行回归诊断' s, h. D$ S9 k# m
windows()
7 {/ `' M1 C. q5 ~par(mfrow=c(2,2))
4 X, H" W( d( s- l# Nplot(lm4.aic,which=c(1:4)) + r5 s+ d2 M, N. C0 C2 l9 k
# w5 ^1 g3 j! `5 n; ?# x( b0 a: p3 O5 N9 t7 R" ]- r: L; d
+ p* ~) B& [; [( l2 K4 W. N; \" g1 D1 c9 s7 o9 O' M
##数据变换% d8 X7 P( t) Q/ u/ n
: A' Y s1 |8 x( g9 d2 `! n! e/ \: Z0 f u' c
#box-cox变换, z G `% Z2 \7 M4 Y
library(MASS)
|! Z% N. f& Qb=boxcox(price~as.factor(dis)*as.factor(ring)+as.factor(wuye)+as.factor(fitment)+as.factor(ring)*as.factor(ronggrp)+lv+ratio, data=a,lambda=seq(-3, 3, by=0.1))+ w! ~/ q/ P& | I6 }# j
I=which(b$y==max(b$y)) #定位似然函数最大的位置
3 Z! b6 s, M, A+ ]9 hlambda = b$x[I] #精确的λ值0 D; Z- F" x6 y% ]) o; Z
#λ接近于0,为模型简洁性,可以直接进行对数变换/ [) t2 k O% j7 Z$ H
logprice <- log(price)
" A2 r# r. C$ b, t4 O2 Vhist(logprice)
- E) M/ Q1 T0 [9 t( u! m" m8 `2 Y f' q
+ A) h( y5 d( m* M& Q" e! E##最终模型与诊断6 d2 ]: j+ l" d! @# ~/ d
# @5 O$ A1 E% O2 r$ Z; \' t7 b6 |7 U d! z
lm6=lm(logprice ~as.factor(dis)*as.factor(ring)+as.factor(wuye)+as.factor(fitment)+as.factor(ring)*as.factor(ronggrp)+lv+ratio)
1 R; f" _% x2 p" V/ W3 ^* X& Ewindows()
4 p- A$ z- r. n# k& i' s( [ }2 S7 T% Vpar(mfrow=c(2,2))
/ ^- s$ g* u' r O) Y- C gplot(lm6,which=c(1:4))! [: W7 \5 T/ n4 C! m, i
anova(lm6)9 p3 J9 z$ t4 k/ a N8 \9 l W
summary(lm6)
. u* u% R, v% {! x6 `: K. z* Q( u
( u/ x, b: W5 Y! L& S a: |7 T. `6 [* Z G
请关注数学中国网微博和数学中国公众号,联系QQ 3243710560( u3 _! R2 M" {
h8 `4 v$ t9 t
) N4 n, y8 P% S9 E5 V. ^5 G- u2 C4 ^
) Q( F: {9 n x' D( P% p
|
zan
|