2020 全国大学生数学建模竞赛C题思路+代码

杨利霞

2020 全国大学生数学建模竞赛C题思路+代码
9 T$ z( d" A. W$ c8 L. W题目链接：https://cloud.189.cn/t/ri2uUb7BRVJr
; q; q6 Y. H7 u" u6 ~8 g
4 p  N+ T+ _; ^% a
+ v+ y! B$ z+ [2 m" |) E前言

2 C3 }& ^/ n0 e' y5 V* W
6 b! t# [7 y2 `$ f6 Q# A
5 V' [  u3 `, w, J
又是一年数据挖掘题型，第一次接触这种题型还是在去年的mathorcup上，这种题的难度就在于指标的建立和数据的处理上。后面会出一份关于数据挖掘题型，我的相关经验，常用的工具和代码。


; F! d) `; x! `5 s  T下面的一，二问实际都在解决
) p1 K9 W3 t+ }" d4 a
0 n1 ]! C- K7 A. G& E. d+ K; v
贷不贷款？
贷款金额多少？
数据清洗
/ S& `& Y4 S. f, Q1 @- s6 T- A( i' K这道题的附件数据没有出现缺省或者异常数据，因此对于数据的预处理，更多的是根据问题的需求来做的。
% A1 ]' a- l8 A$ \8 r0 x( U, E( K. v
- r* r9 O3 W/ ~, u. Z8 s
将是否违约，违约设置为1，不违约设置为0


! K. k1 U" m2 O2 o信誉等级ABCD分别对应4,3,2,1
* V: M* z1 |( c  k
9 \1 R0 n. t3 l, y0 B; ]
发票状态，有效发票为a，作废发票为b
! {+ J7 d0 l$ [8 \) Z4 A, P
  P8 m- j$ N9 L9 s4 o9 w
我将销项和进项所有数据，以公司代码为区别，提取到了不同的sheet当中，对于该公司有效发票数，作废发票数，负数发票数，方便对数据观察。


# 遍历所有sheet数据
) J0 C. A; c; l1 m0 @, Cfor xsn in sn.sheet_names[1:]:
- W6 b# N* U* T" X- v7 s    # 读取文件
1 \6 r( w; S- |( w3 k    datas = pd.read_excel(file_pos, sheet_name=xsn)
- X6 w) e! j- L2 H! @9 w    datas['date']=pd.to_datetime(datas['date'],format='%Y/%m/%d')
    datas.set_index('date', drop=True)
    # 找到全部公司名称代号
1 M9 O, T( p1 z1 w* \' p    code_list = list(set(list((datas['code']))))
    for name in code_list:
: ~" O$ F4 G9 F8 @4 Z6 Y, |' G        tmp_datas = datas[datas['code'] == name]
0 [. G5 L, w2 n+ }+ J0 z+ T( z2 `        tmp_datas.index = range(len(tmp_datas))
( e' b/ B7 ?0 M# Y/ M0 z1 Z        # 转换日期未object类型
  G3 z6 [6 n3 ?8 I" i# b* k4 s        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
        count1 = tmp_datas['tax_status'].value_counts()
        tmp_datas['a_count'] = list(count1)[0]
, U6 b) f! O8 i7 ?! z        if(len(count1) > 1):
            tmp_datas['b_count'] = list(count1)[1]
$ S  o$ ^' f# i& q9 O" r        tmp2 = tmp_datas[tmp_datas['cost'] < 0]
        tmp_datas['neg_value_tax'] = len(tmp2)
/ V4 k% e5 z2 G5 t        if xsn == sn.sheet_names[1]:
. H, u% h0 d8 E6 ^6 Y) h  h            tmp_datas.to_excel(writer1,sheet_name=name,index=False)
; ?  f. N0 g5 V- b% t4 j        else:
  T8 a  F- P# `            tmp_datas.to_excel(writer2,sheet_name=name,index=False)
负数发票：在之前购买的物品，并开具了相关正向发票，后来退货所以开具了值为负数的发票，抵消前面正数发票的值。
+ f2 g9 R7 A; `

提取到信息：


! w0 A. Q. J2 H- v部分公司数据记录很少，或者时间跨度大，需要综合数据指标，抵消数据数量和跨度大的影响
8 {' J4 U1 u/ d1 {& R, _有些负数发票，在之前找不到对应的正数发票，可能是因为在数据记录日期之前购买的，在之后退款，因此在附件中找不到记录。
问题一
7 v9 t+ H8 K8 E9 E% Y1 J% W9 p+ A建立指标
  D' _; c4 [" `3 {进项发票作废率，进项负数发票率，进项每月平均交易额，进项每月交易次数，
8 u2 u. @/ v4 ?# Z7 S0 O
4 ~% S; g2 M/ S" A, ?: X
* Y9 u  y6 x, R7 E* \% l' _销项发票作废率，销项负数发票率，销项每月平均交易额，销项每月交易次数，销售收入增长率
% U0 m9 b8 p+ A  L, P+ a

  l2 ~: X4 Q% ~提取出相关指标到附件

& L5 w& a/ o4 |' Z. V
for xsn in sn.sheet_names[1:]:
4 }% `# ?; h9 l" N+ s3 w4 ~. ~    # 读取文件
    datas = pd.read_excel(file_pos, sheet_name=xsn)
    code_list = list(set(list((datas['code']))))
8 o4 ~* N5 o  E+ h8 M5 c; b# _    for name in code_list:
        tmp_datas = datas[datas['code'] == name]
  d! l( w- p1 }4 ?2 o) x        tmp_datas.index = range(len(tmp_datas))
        insert_datas.append(name)
: Q) P. z0 P# p+ }/ t2 F( ]7 L        # 作废数
        cacel_count = len(tmp_datas[tmp_datas['tax_status'] == 'b'])
4 z. v( F- t( M" X3 I        # 有效数
% {0 E0 p9 f! s8 D4 Y        valid_count = len(tmp_datas[tmp_datas['tax_status'] == 'a'])
- c& a6 P) Z! @% i6 j; l) ]$ g        # 发票作废率
        count1 = (cacel_count / (cacel_count + valid_count))*100
        # 负数发票数
/ @2 f* ?$ A9 S- W+ ?        neg_count = len(tmp_datas[tmp_datas['cost'] < 0])
        # 负数发票率
        count2 = (neg_count / valid_count) * 100
5 {# E( }9 c7 |& }& t/ @/ n' q        # 转换时间
1 G, P: e) s2 C3 I9 @        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
6 `0 ^0 w- ?" R' N7 I+ E7 u4 n3 u$ @        # 时间最大值
        max_time = tmp_datas.iloc[0:,1].max()
        # 时间最小值
        min_time = tmp_datas.iloc[0:,1].min()
        # 时间差
        diff_time = months(max_time, min_time) + 1
        # 有效票
        valid_tax = tmp_datas[tmp_datas['tax_status'] == 'a']
' a8 g) E; b4 n# E8 p7 M        # 平均月交易额
- K9 K/ K* r- b        avg_money = valid_tax['totle_cost'].sum() / diff_time
+ v3 r+ E; b- O  y        # 平均每月交易次数
        trans_count = len(tmp_datas) / diff_time
        insert_datas += [count1, count2, avg_money, trans_count,]
        if flag:
            df1.loc[len(df1)] = insert_datas
# I% M2 \, k2 O% [; k            df1.to_excel(writer1,sheet_name='进项信息',index=False)
+ E8 V; N9 C* i+ N* i% k1 ?        else:
7 b1 C- W" c+ O1 J/ J- E" u# J$ d            merge_time = tmp_datas.groupby(tmp_datas['date']).sum()
5 Q. k8 t# m8 O: ]& M1 P            # 销售收入增长率
& {* K& ?5 j7 C* ]            income_info= list((merge_time['cost'] - merge_time['cost'].shift(1)).fillna(1))
            diff_time_day = days(max_time,min_time)
4 O. C9 {3 C+ f: u) x; J( n            income_tax = (sum(income_info) / diff_time_day)*100
/ ?$ z5 ?8 f" W( J9 |6 u            insert_datas.append(income_tax)
            df2.loc[len(df2)] = insert_datas
            df2.to_excel(writer1,sheet_name='销项信息',index=False)
) D0 Z6 m+ r0 q1 k        insert_datas = []
: m5 B5 q, w& A: `9 r$ h0 e    flag = False
! f1 H9 H" }% a9 e9 W& H" A
) n3 s0 y9 a& w7 [5 |8 `, q
8 i+ w( ~# b1 ?2 j1 ^5 n& N

* V( i! Z4 ?  _1 v3 D并将是否违约插入到最后一列


4 e2 q9 V/ i1 T% D# 提取是否违约的列表
m = []
  c- X/ F2 W' i! d* m5 c3 h7 |for name in code_list:
& |$ x! m& M. j' N    m.append(datas[datas['code']==name]['break_contract'].tolist()[0])
    df1.loc[:,len(df1)] = m
) |. X. A7 z* \: X/ k% X    df1.to_excel(writer3,sheet_name='sheet1',index=False)
# c' b5 y2 k$ c建立模型
Logistics违约率预测模型
* _5 T9 j- @. ?) O使用Logistics违约预测模型，代入所有的指标数据为自变量，是否违约为因变量，预测出违约率。

$ ]8 K# m& ^3 L
X=datas[['进项发票作废率','进项负数发票率','进项每月平均交易额','进项每月交易次数','销项发票作废率','销项负数发票率','销项每月平均交易额','销项每月交易次数','销售收入增长率']]
y=datas['是否违约']
X_train, X_test, y_train, y_test = train_test_split(X,y, test_size=0.2,random_state=2020)
X_validation, X_test, y_validation, y_test = train_test_split(X_test,y_test, test_size=0.1,random_state=2020)
model = LogisticRegression()
model.fit(X_train,y_train)
  k) R3 f9 b/ B) T! pa=model.predict_proba(X_validation)
result=[]
- P3 F2 p' l2 }for i in range(len(a)):
) w9 A3 o8 \/ c" E9 h& sif a[1]>0.5:
! C1 }# W- J8 ?& e    result.append(1)
" @. c/ S2 D: melse:
    result.append(0)
' N) n0 L) ~: G! Z; [from sklearn import metrics
print('误差: %.4f' % (1-metrics.recall_score(y_validation,result,average='weighted')))
最终得到一张我们的分析表格
& w( M& ~: i3 I' E! V
/ t, h1 Z. N( u+ O# u+ E




& w; x7 {9 ^* p- c! ~通过预测是否违约，我们就能解决贷不贷款的问题。


) k* d6 K( M2 l- o; q8 K贷款金额
贷款金额的确认，根据该公司不违约率在所有公司中的权重，乘以总贷款金额确认：


0 q0 L0 y! N+ L
7 ^, Z9 y) J. U$ E  [: [

0 U) j6 {/ O9 Dri=1−Zi∑123j=1(1−Zi)×M
% v3 L0 u5 n4 P! X8 _# m" `1 I" vri=1−Zi∑j=1123(1−Zi)×M
% h& P: w* i, s: s) j因此，我们得到的贷款金额是违约率和贷款总金额组成的关系式，这在第二问中能起到重要作用。


1 v8 U: D- }; W2 i/ Q贷款年利率
绘制出年利率与客户流失率图，可以分析出两者应该是有关系的。利用SPSS拟合出不同信誉等级，年利率与客户流失率的关系式。
- H4 G( H" b; y) M7 g% x
1 l/ T2 P' D$ c# Q1 a$ V
. E  t6 ~9 D8 g- f5 ~4 s5 W8 _. z信誉等级        R平方        关系式
A        0.9977        y = 37.97x^3-258.57x^2+640.944*x -1.121
) C+ c2 a/ O& ?) qB        0.9982        y = 33.995x^3-225.051x^2+552.829*x-1.017
6 k8 S  v9 @' H1 @  v% ~C        0.9982        y = 32.157x^3-207.386x^2+504.717*x-0.973
银行获利=贷款金额x贷款年利率x(1-利率对于信誉评级客户流失率)

+ `) t5 D; Q7 h5 R4 g1 ~$ G
' x/ ]4 |6 k+ _% a6 |在贷款金额确认，贷款年利率范围在0.4~1.5的情况下，利用上面拟合的关系式，我们能够暴力跑出最优年利率。

' D# A, u* u7 M5 v( W
4 v8 d' R5 r: A4 J3 Jdouble turnover_rate(double x, char ch) {
        double y = 0, result = 0;
% p: D- @5 ]/ c: L0 v* B        switch (ch) {
        case 'A':
  S# v4 g' a5 L( F2 `4 R; g0 \                y = 37.969520 * pow(x, 3) - 258.570452 * pow(x, 2) + 640.944427 * x - 1.121484;
9 j9 w# f! p; Y( i                result = x * (1 - y / 100.0);
+ Q# s4 `' b  U                break;
) r# `' s( e+ Y: s        case 'B':
) P( B# a5 q( a0 Z7 f7 z                y = 33.994698 * pow(x, 3) - 225.050538 * pow(x, 2) + 552.829151 * x - 1.016503;
                result = x * (1 - y / 100.0);
                break;
# q6 y$ F% e/ T, \2 R        case 'C':
7 e/ F# W" w  H$ D8 B0 d2 O                y = 32.156864 * pow(x, 3) - 207.385880 * pow(x, 2) + 504.716993 * x - 0.973497;
                result = x * (1 - y / 100.0);
                break;
, ?: m/ x2 T5 }1 p5 q        default:
                cout << "输出有误！"  << ch << endl;
        }
+ ^- v2 g+ G8 d, N& q3 t+ t
( {$ A0 N! q( t+ [0 C        return result;
' q+ G" p2 N6 H}
& |, C/ y+ f; \1 m0 s3 [9 y

# P; q$ V4 Z" h8 c5 ~

问题二
1 Y( [( R9 T5 h/ c利用代码，重新计算出各指标数据
代入Logistics违约率预测模型，预测出各公司的违约率
, e' s) N( t: t4 E1 z6 b9 j根据标准普尔评级建立，主标尺，对不同违约率进行A~D等级划分，信誉等级D不予贷款
# b: R1 N7 f+ Y6 n- g9 ^1 T+ f将违约率代入，之前得到的公式，得到具体贷款金额
最优年利率沿用上一问
# 信用等级
: v6 Z  e8 |$ g9 O0 ]- Kcs = []
8 ]' T% o7 b& |% _4 q# 最优年利率,客户流失率,利率值
. z+ o; A" Y0 [, q. a+ btax = []
for i in m:
    if i <= 0.0069264:
; C+ ?7 ]6 j1 r' ^# @        cs.append('A')
) r7 K2 x4 [% i+ {. e8 e. ^        tax.append([0.083,0.503173,0.0412366])
0 A- K% ^% y5 a+ Q    elif i > 0.0069264 and i <= 0.22619:
        cs.append('B')
        tax.append([0.097,0.505215,0.0479942])
; ~7 A- E8 N; |! d& b1 x' A    elif i > 0.22619 and i <= 0.509915:
9 @9 f- n/ }0 Z* v        cs.append('C')
" g! X/ s5 L2 m0 o1 R. s        tax.append([0.1069,0.506501,0.052755])
    elif i > 0.509915:
        cs.append('D')
        tax.append([0.15,0,0])
7 A3 D4 [: B; b    else:
        print('违规')
        
parr = []
for arr in list(a):
    parr.append(list(arr)[0])
% ?$ V* R  P0 |# t) O' @( Psum_val = sum(parr)
2 F  k$ W! z! a) b+ mamount = []
0 U; s  y4 T! v" Z4 G% X; \7 [for ival in parr:
8 Y8 R* _: O8 r    tmp = ival / sum_val * 100000000
/ g6 a2 X) |' W8 P8 ?9 t! H7 y    if ival < 1 - 0.509915:
            amount.append(0)
! `; z+ h: w. e    else:
8 E7 m% ~. G7 T- s5 p$ F            amount.append(tmp)
! q  f6 U& Q) O1 P, X8 ]6 F

* [& r" ]0 ?/ I1 \' |( E
! S, ]. b! k3 w. Z( V! ?9 C; V- W6 T
3 I9 _1 V& ?) I. R: P! l可以看到，贷款金额也都在10w~100w之内。


8 s, i- q( r6 R6 u+ E; |问题三
这一问，我们做得有些匆忙了，有其他想法的可以按照自己的想法做做，这里只拿我们的做参照。


疫情对公司影响最大的就是每月平均销售额和每月平均销售数量，因此，
  h" h) n( X+ S/ g* t6 s# Q) Q

. Y3 ~8 G5 G' f) E$ c6 F9 A  |! \对每月平均销售额和每月平均销售数量，分别取随机数，数量取10w组，其他指标数据值不变，是否违约数据根据第二问结果，设为初始值
代入Logistics模型中，预测出每一组的违约率
判断每个公司违约率的变化情况，根据变化情况来增/减贷款金额和年利率。
8 J; v( {, [, L8 e+ grand_num = pro_rand()

data_form = {'E377':[], 'E311':[], 'E297':[], 'E386':[], 'E233':[], 'E367':[], 'E194':[], 'E196':[], 'E249':[], 'E205':[], 'E159':[], 'E300':[], 'E346':[], 'E395':[], 'E360':[], 'E200':[], 'E195':[], 'E247':[], 'E317':[], 'E303':[], 'E135':[], 'E217':[], 'E298':[], 'E325':[], 'E335':[], 'E157':[], 'E287':[], 'E270':[], 'E286':[], 'E212':[], 'E261':[], 'E423':[], 'E390':[], 'E189':[], 'E129':[], 'E366':[], 'E191':[], 'E246':[], 'E406':[], 'E357':[], 'E349':[], 'E387':[], 'E372':[], 'E385':[], 'E209':[], 'E267':[], 'E414':[], 'E316':[], 'E341':[], 'E281':[], 'E206':[], 'E363':[], 'E166':[], 'E190':[], 'E140':[], 'E130':[], 'E225':[], 'E347':[], 'E356':[], 'E185':[], 'E192':[], 'E207':[], 'E234':[], 'E136':[], 'E379':[], 'E274':[], 'E383':[], 'E242':[], 'E361':[], 'E408':[], 'E198':[], 'E204':[], 'E413':[], 'E389':[], 'E253':[], 'E226':[], 'E231':[], 'E182':[], 'E318':[], 'E392':[], 'E275':[], 'E425':[], 'E388':[], 'E305':[], 'E155':[], 'E348':[], 'E400':[], 'E256':[], 'E351':[], 'E201':[], 'E345':[], 'E278':[], 'E306':[], 'E308':[], 'E215':[], 'E382':[], 'E407':[], 'E376':[], 'E291':[], 'E369':[], 'E139':[], 'E296':[], 'E260':[], 'E145':[], 'E257':[], 'E216':[], 'E125':[], 'E312':[], 'E354':[], 'E380':[], 'E326':[], 'E268':[], 'E186':[], 'E179':[], 'E223':[], 'E319':[], 'E422':[], 'E224':[], 'E151':[], 'E293':[], 'E302':[], 'E309':[], 'E378':[], 'E373':[], 'E364':[], 'E144':[], 'E162':[], 'E208':[], 'E399':[], 'E355':[], 'E197':[], 'E375':[], 'E368':[], 'E334':[], 'E172':[], 'E254':[], 'E352':[], 'E230':[], 'E313':[], 'E255':[], 'E273':[], 'E220':[], 'E169':[], 'E126':[], 'E237':[], 'E153':[], 'E164':[], 'E344':[], 'E284':[], 'E320':[], 'E328':[], 'E143':[], 'E338':[], 'E290':[], 'E251':[], 'E314':[], 'E330':[], 'E424':[], 'E183':[], 'E248':[], 'E181':[], 'E271':[], 'E403':[], 'E359':[], 'E324':[], 'E304':[], 'E213':[], 'E412':[], 'E337':[], 'E235':[], 'E371':[], 'E283':[], 'E391':[], 'E370':[], 'E173':[], 'E158':[], 'E171':[], 'E299':[], 'E310':[], 'E174':[], 'E295':[], 'E402':[], 'E419':[], 'E327':[], 'E294':[], 'E339':[], 'E203':[], 'E285':[], 'E410':[], 'E241':[], 'E152':[], 'E228':[], 'E245':[], 'E263':[], 'E404':[], 'E160':[], 'E292':[], 'E397':[], 'E329':[], 'E163':[], 'E301':[], 'E365':[], 'E178':[], 'E288':[], 'E175':[], 'E211':[], 'E252':[], 'E276':[], 'E156':[], 'E232':[], 'E321':[], 'E128':[], 'E150':[], 'E148':[], 'E161':[], 'E222':[], 'E142':[], 'E134':[], 'E238':[], 'E396':[], 'E239':[], 'E227':[], 'E147':[], 'E188':[], 'E243':[], 'E401':[], 'E362':[], 'E418':[], 'E272':[], 'E333':[], 'E405':[], 'E416':[], 'E358':[], 'E202':[], 'E259':[], 'E180':[], 'E409':[], 'E244':[], 'E394':[], 'E282':[], 'E210':[], 'E124':[], 'E398':[], 'E374':[], 'E277':[], 'E265':[], 'E421':[], 'E342':[], 'E322':[], 'E165':[], 'E340':[], 'E381':[], 'E331':[], 'E280':[], 'E393':[], 'E353':[], 'E221':[], 'E219':[], 'E323':[], 'E336':[], 'E199':[], 'E170':[], 'E137':[], 'E168':[], 'E350':[], 'E141':[], 'E214':[], 'E177':[], 'E250':[], 'E269':[], 'E167':[], 'E307':[], 'E420':[], 'E411':[], 'E262':[], 'E176':[], 'E218':[], 'E229':[], 'E315':[], 'E264':[], 'E127':[], 'E132':[], 'E154':[], 'E133':[], 'E332':[], 'E138':[], 'E258':[], 'E266':[], 'E236':[], 'E415':[], 'E146':[], 'E193':[], 'E187':[], 'E384':[], 'E417':[], 'E289':[], 'E279':[], 'E131':[], 'E184':[], 'E240':[], 'E149':[], 'E343':[]}
; Q' [0 N( v/ |" Fdf = DataFrame(data_form)
* @, Y. `* \  L5 d+ G/ X: _- B
( `- x5 F! h/ a) rfor val in rand_num:
        a = rate_func(val)
        rate_list = get_rate(a)
        df.loc[len(df)]  = rate_list
        df.to_excel(writer,sheet_name='违约率变化',index=False)
' A8 W& L; L8 _8 S1 I2 b销售数量和销售金额的随机数范围是：0~MAX
4 v& N' a7 U. \, L

1 ~' k" w, ?9 x! |: d
$ {9 v- ~1 g. V2 H* A
( g) I" P2 T/ F8 B& e# R( N0 c————————————————
" |& B; U" v3 b0 K* {/ C6 d# d$ z5 o版权声明：本文为CSDN博主「Hk_Mayfly」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/qq_39542714/article/details/108592087