2020 全国大学生数学建模竞赛C题思路+代码 5 |, w' |1 v( M8 O" c) m7 O1 |题目链接:https://cloud.189.cn/t/ri2uUb7BRVJr0 ~8 ]+ K: z# E' o- l# ~% m
8 u4 c4 y' |5 K; ]! r- L, L% n
9 u" }2 w$ l1 L E
前言5 \# F9 F& F1 \3 r) L" H
, y7 L' ~$ Y, x# A & W9 c' [- S' g6 k5 d# {/ `" F7 y1 O
- O0 {6 y- z; L' i1 {! X" i
又是一年数据挖掘题型,第一次接触这种题型还是在去年的mathorcup上,这种题的难度就在于指标的建立和数据的处理上。后面会出一份关于数据挖掘题型,我的相关经验,常用的工具和代码。( Q) e2 {6 n4 J5 E9 a, s6 }
* v. `* ?) a4 r
8 Y3 u8 H v% _+ J# w下面的一,二问实际都在解决 6 E- Q, h2 z+ }8 W' L4 C4 c* X* D+ O& _8 e" j. ~
; M0 m: u' D3 Y A- \$ [贷不贷款? 3 y. [" T) P6 x \# m% S! k# U; z贷款金额多少? : o# c' @; L7 `" E/ |数据清洗. q. [/ _/ {) h: e8 F4 t6 z
这道题的附件数据没有出现缺省或者异常数据,因此对于数据的预处理,更多的是根据问题的需求来做的。 ! I5 O6 j1 i3 l1 L4 U: [# T3 _$ T1 T2 q( u& K& ?
0 A5 W5 U; c$ N3 j5 W
将是否违约,违约设置为1,不违约设置为0 8 j a. Q# G* U $ E' p4 }6 c. r/ O1 M' }0 ]% E( A
信誉等级ABCD分别对应4,3,2,17 v% a+ X5 B5 {4 ?* e% Q
( s- W2 G* @+ S# O0 C
7 n& j% f7 g- w8 ]0 U发票状态,有效发票为a,作废发票为b ' f1 Y3 C! W8 _- d/ x+ Z5 Y3 V- r# ` + e: K2 H1 y7 I8 @0 M7 A6 W5 q/ f- c1 Q: ]" o' V( R7 H) O- `8 h" U
我将销项和进项所有数据,以公司代码为区别,提取到了不同的sheet当中,对于该公司有效发票数,作废发票数,负数发票数,方便对数据观察。 # A% l. ]; v, b8 w$ }6 x5 l9 }) y
. f/ i9 G: q5 E; U1 q( t; e
# 遍历所有sheet数据3 X8 R, B! e Q5 h: m9 w
for xsn in sn.sheet_names[1:]: 3 W6 r% S, @! _5 N. k6 _* [. c( @ # 读取文件 1 J) b" |0 e I3 \. |+ I, k8 H# n$ @ datas = pd.read_excel(file_pos, sheet_name=xsn)5 @( J8 I; _- A( D% [
datas['date']=pd.to_datetime(datas['date'],format='%Y/%m/%d'), ], I5 P4 T+ S1 M
datas.set_index('date', drop=True) / p" h& _: @. m; y9 i# Q # 找到全部公司名称代号 9 H% r" w. `, F9 |% F3 z code_list = list(set(list((datas['code'])))) ! T( c4 R& m% Q6 G! Z for name in code_list: + Z0 D. \$ M, N7 c( c u' z+ | tmp_datas = datas[datas['code'] == name] {. A* V0 A, h( k7 T
tmp_datas.index = range(len(tmp_datas))$ j/ M" J% r& H4 c. Q! q/ l% N
# 转换日期未object类型& `1 g2 w- N3 W1 G
tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']] . Z. z+ b4 c3 t# j$ e count1 = tmp_datas['tax_status'].value_counts() $ p" C' C" m6 x0 X! T tmp_datas['a_count'] = list(count1)[0] % C% Q4 }2 O$ n0 i5 {& i if(len(count1) > 1): S& [1 L* V! S1 f" n5 _$ A tmp_datas['b_count'] = list(count1)[1] * a% W( I+ ~4 H. W tmp2 = tmp_datas[tmp_datas['cost'] < 0] 8 \0 E W" f) Y6 H tmp_datas['neg_value_tax'] = len(tmp2)- E8 d' H! v5 ]% n* ~
if xsn == sn.sheet_names[1]:1 x! `, C6 F% z4 h# i
tmp_datas.to_excel(writer1,sheet_name=name,index=False)8 }8 |- I, m+ w$ q; x9 F
else:0 G6 [4 U' p& P. ~% _+ i5 u& S
tmp_datas.to_excel(writer2,sheet_name=name,index=False) 1 V$ e6 J6 u+ n1 F$ m负数发票:在之前购买的物品,并开具了相关正向发票,后来退货所以开具了值为负数的发票,抵消前面正数发票的值。 9 I1 `. {! R6 @8 R% z 6 `5 c$ I; f% L1 @$ K5 E. z1 ~- v$ _* Z- L# D
提取到信息: ' B2 F8 `- w4 V8 Y1 }# G # U6 t4 }8 e2 C2 b) Q1 D ) w: y" Z& f$ T部分公司数据记录很少,或者时间跨度大,需要综合数据指标,抵消数据数量和跨度大的影响 / m6 D3 U! ]' }: J9 v% _2 g有些负数发票,在之前找不到对应的正数发票,可能是因为在数据记录日期之前购买的,在之后退款,因此在附件中找不到记录。( ?2 M, X& J' p" R& ]) K
问题一 $ Z# |4 n r; G: ^建立指标 ( [* A7 t! _% i% a% L进项发票作废率,进项负数发票率,进项每月平均交易额,进项每月交易次数, 2 n# f7 F3 g. G5 d- z% o2 ~2 E& s+ _
( {4 p( k1 l E% u销项发票作废率,销项负数发票率,销项每月平均交易额,销项每月交易次数,销售收入增长率7 n& t+ L1 R U. Q5 K- w. k$ l9 l0 e4 G
5 M; r$ [0 q I2 P3 r o+ W: d7 l4 A U, D, F% q
提取出相关指标到附件& C, S }8 L/ E! F+ \% C+ \# o8 N+ {% A$ y
' j( v7 Q( f% ` : S5 y, u2 T- Q9 c( t1 k1 |for xsn in sn.sheet_names[1:]:7 d5 W* K6 n) _* a9 m0 h5 s
# 读取文件' A9 D8 m7 s% ~; J
datas = pd.read_excel(file_pos, sheet_name=xsn) & G: ?% I9 \- P code_list = list(set(list((datas['code']))))8 @: W0 p" m- o1 e
for name in code_list: ( V8 A5 ? \' v" Z: N/ ^2 ]3 W tmp_datas = datas[datas['code'] == name] 8 C6 |6 D5 y* i8 b* v tmp_datas.index = range(len(tmp_datas)) * m9 s( d0 [0 T2 g insert_datas.append(name), k# g; Q) D6 X1 H- e
# 作废数 # h) u% [5 |6 l+ [4 G cacel_count = len(tmp_datas[tmp_datas['tax_status'] == 'b'])% g+ ^/ e% `/ O8 \
# 有效数 7 g5 f- P, J/ i, ]! L) L valid_count = len(tmp_datas[tmp_datas['tax_status'] == 'a']) ) c8 L8 ?" q3 L; ? # 发票作废率/ ?' r; u8 g3 i0 |. b. b' t
count1 = (cacel_count / (cacel_count + valid_count))*1003 @ N _$ ~3 p5 U0 _
# 负数发票数9 G0 M6 }/ J- |
neg_count = len(tmp_datas[tmp_datas['cost'] < 0])/ ^1 @' o4 S3 [$ F1 D9 r' j' Q* J; N
# 负数发票率9 E& ?/ O) h! z/ K
count2 = (neg_count / valid_count) * 100 ! R2 l7 N; v; y) [ # 转换时间 1 @% I1 `1 p0 ~! j! x/ ^ tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']] 8 c o; {: o# ` # 时间最大值 8 m1 c/ l$ `! x, H, r7 a+ i max_time = tmp_datas.iloc[0:,1].max()" @. R; D6 c2 N% [1 O( i
# 时间最小值 b( G3 w- P: @* K+ b
min_time = tmp_datas.iloc[0:,1].min()1 P7 N1 W% _2 @' `3 q& S; w, M4 ~: ~
# 时间差( D, a5 M2 z7 ^( }3 ~2 g( M. h8 H
diff_time = months(max_time, min_time) + 1 5 \2 ^. c2 h: k4 c# g # 有效票 ; P# I5 G; {0 ]" W2 l @8 a valid_tax = tmp_datas[tmp_datas['tax_status'] == 'a']; p, ]( h- n% w
# 平均月交易额 0 r; `. {/ R: L+ O0 B+ v- t# S! c avg_money = valid_tax['totle_cost'].sum() / diff_time ( Z7 K8 L( j8 D2 ]& z L, _ # 平均每月交易次数* O+ r; P. f9 ~: y3 p8 ]$ ]: S
trans_count = len(tmp_datas) / diff_time$ R- u& y; z+ H: f/ E
insert_datas += [count1, count2, avg_money, trans_count,] q6 s1 T( t$ f+ l2 ?0 R
if flag:9 w( x9 n1 W+ ~" p- s" C: P
df1.loc[len(df1)] = insert_datas3 q$ z6 _, e y {; m- B* y
df1.to_excel(writer1,sheet_name='进项信息',index=False) 5 x, M" E- y' \3 e else: ( p- U6 h% p/ O' @- t5 h merge_time = tmp_datas.groupby(tmp_datas['date']).sum() 2 G' W, q/ C2 Q+ _3 p. w9 y # 销售收入增长率 4 J( | q: G/ r income_info= list((merge_time['cost'] - merge_time['cost'].shift(1)).fillna(1))# i' N/ [" S" X# {
diff_time_day = days(max_time,min_time) 6 M/ \5 H: q+ c) o( p) r6 i income_tax = (sum(income_info) / diff_time_day)*100 ; d V4 Z6 G7 b! j9 v insert_datas.append(income_tax) / @- M, i2 S# v" q( E# W df2.loc[len(df2)] = insert_datas" P3 ^2 R+ t2 S! F) T/ ^
df2.to_excel(writer1,sheet_name='销项信息',index=False)% R5 {+ B. c1 [' K
insert_datas = []+ a' Y- k+ F! w3 ]( {
flag = False * C+ G! [9 G" X: N 6 W' V' A3 k; ]0 ]- k , [& Y$ P" S, S5 M+ R, J% t5 W) D e) y8 r1 x& I# p' o8 w d
) a8 }. n! ]. ]; T并将是否违约插入到最后一列 1 b* R% u: y4 N ( j9 {. w$ ^' u: j$ {* k) V" |. k% ~0 \/ W
# 提取是否违约的列表( ]0 e' L" M6 }: e G! t4 F5 F
m = [] 1 D, W Z% P1 Z H7 Efor name in code_list:/ f* r1 i6 k( K2 x8 `& `' n2 G& ]
m.append(datas[datas['code']==name]['break_contract'].tolist()[0])! m: C! h* \) J
df1.loc[:,len(df1)] = m. ~0 p; ~4 f' g
df1.to_excel(writer3,sheet_name='sheet1',index=False)0 I0 C# O% G8 ~+ N) F
建立模型! i. s$ i0 h9 r/ I
Logistics违约率预测模型* \' C7 B. v* u+ m! X
使用Logistics违约预测模型,代入所有的指标数据为自变量,是否违约为因变量,预测出违约率。+ x) u! @# N4 K! K
6 a5 a* ]) a! f. e m3 V) o) m+ ?4 i4 J/ L! R
X=datas[['进项发票作废率','进项负数发票率','进项每月平均交易额','进项每月交易次数','销项发票作废率','销项负数发票率','销项每月平均交易额','销项每月交易次数','销售收入增长率']]9 T3 R0 g- O5 `: x5 K+ v q
y=datas['是否违约']8 E0 P8 u- B7 g$ r4 K( t b
X_train, X_test, y_train, y_test = train_test_split(X,y, test_size=0.2,random_state=2020). M i! |' f, j d4 u) r- }
X_validation, X_test, y_validation, y_test = train_test_split(X_test,y_test, test_size=0.1,random_state=2020)6 x6 I% s6 V; o* K
model = LogisticRegression() # \$ x c% w0 x. n5 Zmodel.fit(X_train,y_train) ( D4 G6 B/ v1 |9 r- {6 T$ va=model.predict_proba(X_validation) 3 |$ _7 x1 S9 _8 F. ?' v. mresult=[] 3 |( p( p5 J, D2 [7 i2 Kfor i in range(len(a)):! K* p) i1 M. C& u
if a[1]>0.5:0 d0 D: x$ L- g0 f3 I
result.append(1)) g: [" _* N) i
else:8 _2 f. A& p( ?1 T+ h1 Y1 N
result.append(0)8 F. `3 ~+ f9 D; t4 j* u6 _
from sklearn import metrics % [4 T g. c/ N" S2 m! Gprint('误差: %.4f' % (1-metrics.recall_score(y_validation,result,average='weighted'))) : `, O' d" X2 Z. x% M) t最终得到一张我们的分析表格 & w4 p# I) u) N8 l; i7 f! u$ K1 ?5 R& k
6 ]5 G6 s0 L: H, k! }4 }
( I* H2 X3 }; u3 [4 n6 k4 M) F, g9 ~1 N" c* K
, Q0 q& D4 ^% B d- [/ x& T R6 L/ e: ?/ r
通过预测是否违约,我们就能解决贷不贷款的问题。 " A$ d/ [4 A; o5 s; K# u M4 X9 b, v
& }, }% S; j$ R0 o
贷款金额 8 y. H5 g/ g! C- a8 A6 `贷款金额的确认,根据该公司不违约率在所有公司中的权重,乘以总贷款金额确认:7 J2 U5 p; o+ j5 [1 u) d* s
! H; x/ j( U: D7 h" ?6 G" B2 A2 z4 g5 r8 `
8 U! |# e& w0 s
/ R7 I/ E" i& {
, } I8 k' u& y4 I0 ?9 E
ri=1−Zi∑123j=1(1−Zi)×M ! N! t$ e! ~# b- ?+ o2 Hri=1−Zi∑j=1123(1−Zi)×M ' m0 l6 g, n& c1 {5 C0 k6 f' t因此,我们得到的贷款金额是违约率和贷款总金额组成的关系式,这在第二问中能起到重要作用。 ( [6 F2 o y8 f. Z! T d9 ~% K" O
5 Q+ |" B' c, K2 P* ?! T$ r
贷款年利率# C0 U" ]1 b3 i9 C7 T" F
绘制出年利率与客户流失率图,可以分析出两者应该是有关系的。利用SPSS拟合出不同信誉等级,年利率与客户流失率的关系式。 - a7 ?% p& k, U5 c, K4 T+ U" s' U% e- c3 x9 E
, g- P& J) b! y5 v$ E
信誉等级 R平方 关系式 3 Q5 G/ D# }% ]# D: bA 0.9977 y = 37.97x^3-258.57x^2+640.944*x -1.121/ f j: Y) u# U9 _( g
B 0.9982 y = 33.995x^3-225.051x^2+552.829*x-1.017! p! K# f# u( b( Z- l" _6 e
C 0.9982 y = 32.157x^3-207.386x^2+504.717*x-0.973 2 G* k1 n$ v1 `% a0 d2 e+ {4 E Q银行获利=贷款金额x贷款年利率x(1-利率对于信誉评级客户流失率)" l$ w5 r$ O/ u/ y" X; R: X
! ?9 n' ~# Q# ?2 Y/ k
: E- m7 c% J1 l8 k
在贷款金额确认,贷款年利率范围在0.4~1.5的情况下,利用上面拟合的关系式,我们能够暴力跑出最优年利率。, h7 N! m. ~- {. G B* @
2 ~# U, E5 f% c3 i ! O+ o. Q2 L) d9 W6 ^5 ?* S1 bdouble turnover_rate(double x, char ch) { , g) y8 r5 X! r" {$ @ double y = 0, result = 0; & L/ \6 i, L/ j! j5 q' \1 Y# i" v switch (ch) { 9 z4 _) p" f0 L2 `% x case 'A': . a' Y: @2 a% R3 R1 t# [ y = 37.969520 * pow(x, 3) - 258.570452 * pow(x, 2) + 640.944427 * x - 1.121484; $ f8 N* t' C% V0 \5 J$ m result = x * (1 - y / 100.0); 8 p) A% P& a6 G+ | break;0 ]! r: M# S ] j# v
case 'B':$ {& v! X* Q; O: C
y = 33.994698 * pow(x, 3) - 225.050538 * pow(x, 2) + 552.829151 * x - 1.016503;& C; S& I, ]; e. w+ s
result = x * (1 - y / 100.0);. O+ G' L* t0 f& W# k9 _! z! v( z
break;/ i7 ~8 W- c* Z9 w+ T
case 'C': & a1 X& x& X/ O: l3 ~: E( ~3 r y = 32.156864 * pow(x, 3) - 207.385880 * pow(x, 2) + 504.716993 * x - 0.973497;. k' ~* I& }5 `4 W- Y. R3 `
result = x * (1 - y / 100.0); ; P! o# m- d; x. v" G( n break;0 ?. W2 N1 N5 g: G6 x$ c" y3 d
default: - \* Y4 c$ @* r: F1 S cout << "输出有误!" << ch << endl;" P2 c$ j, r2 @ C8 A, y% a0 S
} ( q5 s7 B! p# k& e , O _7 e. N6 |# Z( F9 c
return result; + u, l1 k6 g* J! k, g}0 a% X1 n% A. a6 i+ F
% S A( m9 R) D1 s# k
. F8 d- B0 s! t# Z5 _6 X D% p7 ]% {+ J. l" r- [2 e8 {$ U- L4 y2 A0 A
问题二0 P& `! a& p* V# L
利用代码,重新计算出各指标数据 4 T, G! ?$ _5 u/ u6 E* ?% a代入Logistics违约率预测模型,预测出各公司的违约率% h# T( ~0 q" V$ I7 a& x4 {- P
根据标准普尔评级建立,主标尺,对不同违约率进行A~D等级划分,信誉等级D不予贷款3 m- ~7 n5 G; C& }5 w
将违约率代入,之前得到的公式,得到具体贷款金额$ I3 T7 Q' d& U0 `0 o
最优年利率沿用上一问 4 Z2 o, k$ f6 Q/ M8 {# 信用等级+ ]( ]5 b. t: r2 K) R& I# ?2 M) e
cs = [] # d& m/ n1 N5 h# 最优年利率,客户流失率,利率值; p. a- p( H; J' M9 ]* U0 ?
tax = []/ e% R- M. l0 z6 \
for i in m: : K' X" r) Z; Y if i <= 0.0069264: - z h0 Q) E9 d0 } cs.append('A'); M( O; h8 ]: v0 P- b
tax.append([0.083,0.503173,0.0412366]) 9 j( @2 q; h" E" ~9 D elif i > 0.0069264 and i <= 0.22619:8 ^# B5 X; u' f; `4 B/ w
cs.append('B') " f% M+ I; k5 Y7 f/ a d; W tax.append([0.097,0.505215,0.0479942])" ^) V' O* J. |1 l. i
elif i > 0.22619 and i <= 0.509915: ; `! p; i @- x8 U4 a3 p4 E& X# K cs.append('C') " } J) h- L7 ?: ?" ]3 A0 M tax.append([0.1069,0.506501,0.052755]) 9 r" N5 x) J4 s x7 ]9 U9 L elif i > 0.509915: ! g/ q; r! W/ t$ P% N cs.append('D')5 b/ ^7 f8 o' B+ f: ~
tax.append([0.15,0,0]) 1 p% d7 u2 {2 Q- m) Q6 V- a else: ! n, E& d5 U' Q* b print('违规') / L5 k4 O2 k6 a* S; X3 N. ?$ l + l: y* l) O8 L! J- T5 e" K+ U
parr = []# G6 O9 [8 y7 m8 _- ^# F" R
for arr in list(a):+ t; j$ x- E9 N! d3 |3 e+ M$ I
parr.append(list(arr)[0])9 _+ d# V( J$ A) } c
sum_val = sum(parr); p2 ]6 r7 N- t# o! Z8 r' w
amount = [] s" E4 m) j- xfor ival in parr: * Q9 [0 T5 E% L tmp = ival / sum_val * 1000000000 B: t- g8 f7 r9 z
if ival < 1 - 0.509915:) s. o9 U: n5 J
amount.append(0); F {# `9 q) \/ R3 H& Z
else: 8 F1 G; @2 t# W4 I amount.append(tmp)6 [2 g+ i; j2 d; b
0 }2 j8 k' Q" Q