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

杨利霞

2020 全国大学生数学建模竞赛C题思路+代码
题目链接：https://cloud.189.cn/t/ri2uUb7BRVJr
/ _& B; V- f1 h+ z6 T5 `
3 A4 p* k0 H( P& w: V( O" R/ ^5 b
前言


; L9 G% ?% n; X6 A1 C

又是一年数据挖掘题型，第一次接触这种题型还是在去年的mathorcup上，这种题的难度就在于指标的建立和数据的处理上。后面会出一份关于数据挖掘题型，我的相关经验，常用的工具和代码。
  b/ ]2 O' ]2 c! x+ ~: a6 K1 C" m

! c, x+ m+ {$ l; q- N下面的一，二问实际都在解决
9 M" f' k9 l7 ~, x' G% }4 n5 N/ u
% Q+ U, ]0 z& F% O! F4 `
贷不贷款？
( H  \4 x. d& C7 W: L3 S贷款金额多少？
$ c9 B6 P( Q. `, c4 m3 w数据清洗
) ~" {3 @! t- E. I4 C这道题的附件数据没有出现缺省或者异常数据，因此对于数据的预处理，更多的是根据问题的需求来做的。
) N* G& e. T# r" C) ?" N1 e
+ M3 r# Q8 \1 O5 @: G% b
将是否违约，违约设置为1，不违约设置为0
, U5 T1 @+ W3 {! F0 I
5 B# A) J" c( z+ H4 A4 e( r
$ c3 C( Z7 p/ G* b6 B: {( _信誉等级ABCD分别对应4,3,2,1
" {; c; N3 |' x! d1 g  O
' J( ]# l7 Z' a/ N7 [
发票状态，有效发票为a，作废发票为b
# y: T5 [8 a; R: n' d5 A  [
( W# ~4 B1 p* t3 x; A8 ?! Q4 K
% T; ~7 T# \3 I; ^& u我将销项和进项所有数据，以公司代码为区别，提取到了不同的sheet当中，对于该公司有效发票数，作废发票数，负数发票数，方便对数据观察。


! [: t- |& y- H4 N8 a  z# 遍历所有sheet数据
for xsn in sn.sheet_names[1:]:
    # 读取文件
0 _/ N0 p. z  y- @7 A    datas = pd.read_excel(file_pos, sheet_name=xsn)
1 V8 t) U" n; i    datas['date']=pd.to_datetime(datas['date'],format='%Y/%m/%d')
  L& l6 M- o! [    datas.set_index('date', drop=True)
, {" n& r3 `4 \) m( f% x. J! O" I    # 找到全部公司名称代号
    code_list = list(set(list((datas['code']))))
    for name in code_list:
+ K$ ?7 U; s& Q! Q' V" N3 O        tmp_datas = datas[datas['code'] == name]
        tmp_datas.index = range(len(tmp_datas))
        # 转换日期未object类型
        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]
        if(len(count1) > 1):
# q/ {6 B2 O8 M. y( ~) X& G) S            tmp_datas['b_count'] = list(count1)[1]
9 C2 y4 U, T( H' ~        tmp2 = tmp_datas[tmp_datas['cost'] < 0]
( G: L, ~. Q! \& s        tmp_datas['neg_value_tax'] = len(tmp2)
        if xsn == sn.sheet_names[1]:
            tmp_datas.to_excel(writer1,sheet_name=name,index=False)
        else:
+ F. Y( g8 g" B5 w, V, n9 E$ t* Z            tmp_datas.to_excel(writer2,sheet_name=name,index=False)
负数发票：在之前购买的物品，并开具了相关正向发票，后来退货所以开具了值为负数的发票，抵消前面正数发票的值。
4 }! a, I* k( F0 p
3 e7 j1 I& L" a" M: @, n$ f
提取到信息：

) y) \4 F7 f! O: D7 T" |" a
部分公司数据记录很少，或者时间跨度大，需要综合数据指标，抵消数据数量和跨度大的影响
8 m3 R; m6 _8 R1 M% \& E. C( e有些负数发票，在之前找不到对应的正数发票，可能是因为在数据记录日期之前购买的，在之后退款，因此在附件中找不到记录。
, s7 R1 J! Q1 X* G5 ^/ v2 H问题一
& L: M* B" p# |  d1 J% E/ w8 v建立指标
3 A; h' w0 ~+ n5 K  y进项发票作废率，进项负数发票率，进项每月平均交易额，进项每月交易次数，


  X: d. Q$ F/ t& t- T6 y销项发票作废率，销项负数发票率，销项每月平均交易额，销项每月交易次数，销售收入增长率


提取出相关指标到附件


for xsn in sn.sheet_names[1:]:
    # 读取文件
    datas = pd.read_excel(file_pos, sheet_name=xsn)
! h8 o' x) u) J  K( [    code_list = list(set(list((datas['code']))))
    for name in code_list:
        tmp_datas = datas[datas['code'] == name]
. i- z! V! u8 Y+ K3 C' E        tmp_datas.index = range(len(tmp_datas))
        insert_datas.append(name)
        # 作废数
        cacel_count = len(tmp_datas[tmp_datas['tax_status'] == 'b'])
9 l9 @0 B0 q; r0 Z8 s( b        # 有效数
( m6 r! z% q# Q& V* h        valid_count = len(tmp_datas[tmp_datas['tax_status'] == 'a'])
1 Q, K8 @$ Q! [. k        # 发票作废率
        count1 = (cacel_count / (cacel_count + valid_count))*100
* T8 j, v' }2 Y2 F+ B& L        # 负数发票数
        neg_count = len(tmp_datas[tmp_datas['cost'] < 0])
        # 负数发票率
- d' |: C' Y' C6 y( \* H        count2 = (neg_count / valid_count) * 100
0 J1 m* `% N9 _! G        # 转换时间
        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
        # 时间最大值
" a) L7 f9 [6 l3 s! d        max_time = tmp_datas.iloc[0:,1].max()
0 m$ J1 c* G) g3 E, ^        # 时间最小值
        min_time = tmp_datas.iloc[0:,1].min()
        # 时间差
        diff_time = months(max_time, min_time) + 1
" h3 P* P" o* m# ^% S        # 有效票
, e0 I, I0 t' U( C; W# Q; t4 d        valid_tax = tmp_datas[tmp_datas['tax_status'] == 'a']
        # 平均月交易额
        avg_money = valid_tax['totle_cost'].sum() / diff_time
7 G  _, c( P9 W1 _        # 平均每月交易次数
        trans_count = len(tmp_datas) / diff_time
        insert_datas += [count1, count2, avg_money, trans_count,]
9 r+ X' J2 B" w        if flag:
3 K: u5 h/ F6 q' }7 F  A$ }            df1.loc[len(df1)] = insert_datas
" q4 l" m! ]0 V, ?9 H- C2 Y            df1.to_excel(writer1,sheet_name='进项信息',index=False)
8 E- {- M# c# Q' @- A, b1 T        else:
            merge_time = tmp_datas.groupby(tmp_datas['date']).sum()
            # 销售收入增长率
% W# M! [- g; G0 a! B3 e, ^            income_info= list((merge_time['cost'] - merge_time['cost'].shift(1)).fillna(1))
            diff_time_day = days(max_time,min_time)
            income_tax = (sum(income_info) / diff_time_day)*100
; U+ P- z6 P; w& u, R" V            insert_datas.append(income_tax)
; i1 [0 q* Y+ ?/ X' x            df2.loc[len(df2)] = insert_datas
9 X) l; W1 K1 V            df2.to_excel(writer1,sheet_name='销项信息',index=False)
        insert_datas = []
+ {+ C. }4 T8 t( t; k1 C    flag = False
- u2 c3 ]( M, g8 a5 y4 B3 N
( T2 i% k: o+ F; p) K2 O; j4 F* C& g


并将是否违约插入到最后一列
$ _8 R& c' Q% M4 X7 k5 s* M

' v1 Q7 `7 ?/ i4 O, f# 提取是否违约的列表
" \+ v( J; T0 `, n$ U5 Lm = []
for name in code_list:
$ b* \4 s3 ]) k    m.append(datas[datas['code']==name]['break_contract'].tolist()[0])
  N3 a9 j1 a% p    df1.loc[:,len(df1)] = m
5 E/ H7 q4 m; d: Q4 y6 y1 K    df1.to_excel(writer3,sheet_name='sheet1',index=False)
建立模型
$ F' A7 u3 s$ y- ?; i0 m6 C& ^Logistics违约率预测模型
使用Logistics违约预测模型，代入所有的指标数据为自变量，是否违约为因变量，预测出违约率。
  w6 n. J7 r& k4 u/ I
, z& G+ E" G/ m8 P% R
3 B! P8 o' m7 K. J7 n* n" zX=datas[['进项发票作废率','进项负数发票率','进项每月平均交易额','进项每月交易次数','销项发票作废率','销项负数发票率','销项每月平均交易额','销项每月交易次数','销售收入增长率']]
: ?9 @8 J: _+ E. a  j7 ^, q9 ry=datas['是否违约']
: p7 U7 K& ?9 U& k" n/ DX_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)
: d/ H3 F% M. k7 N( _2 E6 x' Xmodel = LogisticRegression()
model.fit(X_train,y_train)
; O( C: ]% r( m1 Z1 Xa=model.predict_proba(X_validation)
result=[]
; q6 q& h- C( b$ N" \for i in range(len(a)):
if a[1]>0.5:
+ e3 q% @* z0 M6 O6 N7 |    result.append(1)
8 j7 z/ q) ^/ Y# V, L1 lelse:
    result.append(0)
* c$ _( M* X- Z- K9 T8 W" Xfrom sklearn import metrics
8 w( R) m0 ]1 z7 j- I1 \! _print('误差: %.4f' % (1-metrics.recall_score(y_validation,result,average='weighted')))
3 {& M; ]1 C% {$ j2 p8 J, }最终得到一张我们的分析表格


" P; Q0 j6 X! g( n% V1 u5 {

9 u; E5 E: p5 c
3 W# f  q1 S$ F& P/ c5 Y* y
通过预测是否违约，我们就能解决贷不贷款的问题。

9 J- o8 C# F- \; z
( g& d. _3 B& p. _7 Q9 E; V! Y* o贷款金额
贷款金额的确认，根据该公司不违约率在所有公司中的权重，乘以总贷款金额确认：
7 f! P8 L, E3 \. t1 r

- w- E' W3 N8 p$ h


" F6 r: p0 c6 V4 M: Xri=1−Zi∑123j=1(1−Zi)×M
ri=1−Zi∑j=1123(1−Zi)×M
因此，我们得到的贷款金额是违约率和贷款总金额组成的关系式，这在第二问中能起到重要作用。


贷款年利率
绘制出年利率与客户流失率图，可以分析出两者应该是有关系的。利用SPSS拟合出不同信誉等级，年利率与客户流失率的关系式。


信誉等级        R平方        关系式
, k) ~  v* \7 a: Y  |9 oA        0.9977        y = 37.97x^3-258.57x^2+640.944*x -1.121
0 [- e+ m3 s5 G" Q  IB        0.9982        y = 33.995x^3-225.051x^2+552.829*x-1.017
! U1 U. g1 ~6 h- o" Z1 z' QC        0.9982        y = 32.157x^3-207.386x^2+504.717*x-0.973
7 c$ W( D) ]* W银行获利=贷款金额x贷款年利率x(1-利率对于信誉评级客户流失率)

$ p" K- Q' l7 M: C
4 r/ q9 A0 k' k; a$ H8 M$ l在贷款金额确认，贷款年利率范围在0.4~1.5的情况下，利用上面拟合的关系式，我们能够暴力跑出最优年利率。
; M. t4 I0 r/ B3 {

& ^) u* `0 n! Q  ?6 tdouble turnover_rate(double x, char ch) {
: F1 Y2 {; t' ]; @5 l9 P/ B; m        double y = 0, result = 0;
        switch (ch) {
1 K# ?0 L/ {' i5 E        case 'A':
8 t, Q$ `1 ~9 T4 x8 |  _$ n                y = 37.969520 * pow(x, 3) - 258.570452 * pow(x, 2) + 640.944427 * x - 1.121484;
( |# X; }4 \; T, k" |                result = x * (1 - y / 100.0);
* ^6 Q+ F! d4 z4 L                break;
: D) a& P* j7 p6 u$ L, f% j7 d- Z        case 'B':
# g6 y* a1 P; d                y = 33.994698 * pow(x, 3) - 225.050538 * pow(x, 2) + 552.829151 * x - 1.016503;
- o) }$ @& {3 g. b0 Z" b                result = x * (1 - y / 100.0);
7 d. `+ W+ n- U0 y                break;
        case 'C':
4 ~9 |0 r4 `8 T  e                y = 32.156864 * pow(x, 3) - 207.385880 * pow(x, 2) + 504.716993 * x - 0.973497;
! D0 Z2 ?: W) n% M( P( `7 q, i                result = x * (1 - y / 100.0);
                break;
        default:
( i; H" P( v. q                cout << "输出有误！"  << ch << endl;
        }

        return result;
}

* ]4 M& {: h. M  W3 y+ i

  F! k( P( r3 n( ?4 C/ A# n2 \
问题二
利用代码，重新计算出各指标数据
代入Logistics违约率预测模型，预测出各公司的违约率
根据标准普尔评级建立，主标尺，对不同违约率进行A~D等级划分，信誉等级D不予贷款
将违约率代入，之前得到的公式，得到具体贷款金额
最优年利率沿用上一问
# k( `2 g: N! s% @6 h+ x3 O* X, H# 信用等级
9 }$ s- X' ]" \cs = []
# 最优年利率,客户流失率,利率值
tax = []
for i in m:
    if i <= 0.0069264:
3 d1 i) E& U8 c( `, g6 v3 D        cs.append('A')
        tax.append([0.083,0.503173,0.0412366])
' ^& d5 N4 }3 O/ F; t    elif i > 0.0069264 and i <= 0.22619:
        cs.append('B')
5 m4 d0 N# p5 M) r8 y        tax.append([0.097,0.505215,0.0479942])
    elif i > 0.22619 and i <= 0.509915:
$ v3 `) E6 p6 K  |+ R  o! L        cs.append('C')
! Z# d. a0 D: C6 k& }8 D) R. K) k        tax.append([0.1069,0.506501,0.052755])
    elif i > 0.509915:
- {9 f) u. @! D& o2 Q7 z; q8 A/ Z( i- d        cs.append('D')
        tax.append([0.15,0,0])
    else:
        print('违规')
        
parr = []
for arr in list(a):
9 E) _& ?; F0 F4 T    parr.append(list(arr)[0])
/ d8 z) H' G" I% x, p7 {sum_val = sum(parr)
4 H! g; c+ v3 K& d7 E. Q! B+ Mamount = []
for ival in parr:
    tmp = ival / sum_val * 100000000
6 y% R5 z. Z! t& n    if ival < 1 - 0.509915:
& V) s) K( K0 L" c* z5 m( [            amount.append(0)
    else:
            amount.append(tmp)

( E- s! }' \0 Z9 E/ H, P
- C; w; A' _2 i

* ^3 j0 z( W: q: X  s6 P, m' h可以看到，贷款金额也都在10w~100w之内。
, J& P' a. d8 @/ R& O
$ E4 b5 g1 M1 J8 t3 A1 h
1 i; G$ ^2 s4 o/ ?问题三
7 q: W0 p  @. D$ Q这一问，我们做得有些匆忙了，有其他想法的可以按照自己的想法做做，这里只拿我们的做参照。
$ t& N# F& n; j1 N

+ B6 C0 o2 o6 V! `( o; G; d; p疫情对公司影响最大的就是每月平均销售额和每月平均销售数量，因此，


对每月平均销售额和每月平均销售数量，分别取随机数，数量取10w组，其他指标数据值不变，是否违约数据根据第二问结果，设为初始值
/ W5 b% Z) a+ X; J, H$ p' H代入Logistics模型中，预测出每一组的违约率
判断每个公司违约率的变化情况，根据变化情况来增/减贷款金额和年利率。
rand_num = pro_rand()
. Z- E4 w+ ^$ c8 L* \! Y
9 x& v: A0 l8 ^. i4 J: tdata_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':[]}
df = DataFrame(data_form)
) N$ a( M% a5 J: h( v8 }; V  C
for val in rand_num:
  t( R% {3 e. H' \4 L# h        a = rate_func(val)
        rate_list = get_rate(a)
        df.loc[len(df)]  = rate_list
+ H" F8 I. C2 \; }9 \/ J        df.to_excel(writer,sheet_name='违约率变化',index=False)
: `2 `2 Z$ C- l4 I( B, @销售数量和销售金额的随机数范围是：0~MAX


6 q7 d$ h' K+ f/ |# n

————————————————
版权声明：本文为CSDN博主「Hk_Mayfly」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
; }9 J& X( ^5 j原文链接：https://blog.csdn.net/qq_39542714/article/details/108592087

' a- C4 p* [* C+ l6 u2 E
# T6 ?( m" m0 ?, d! u0 R
; F( c, b6 p4 f7 u) r" f