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

杨利霞

2020 全国大学生数学建模竞赛C题思路+代码
题目链接：https://cloud.189.cn/t/ri2uUb7BRVJr


前言




又是一年数据挖掘题型，第一次接触这种题型还是在去年的mathorcup上，这种题的难度就在于指标的建立和数据的处理上。后面会出一份关于数据挖掘题型，我的相关经验，常用的工具和代码。

6 F  d# y: @$ y* Q* U7 Q0 q
下面的一，二问实际都在解决

- ]2 T" @% H4 v$ F: h+ O
贷不贷款？
贷款金额多少？
数据清洗
这道题的附件数据没有出现缺省或者异常数据，因此对于数据的预处理，更多的是根据问题的需求来做的。
3 V$ [1 Z; o' ^, r

, I3 _/ g9 O* o! C$ k将是否违约，违约设置为1，不违约设置为0
& k* j  v9 u% Y& @0 x0 _

信誉等级ABCD分别对应4,3,2,1
& f1 R2 \3 X, Z! @" P$ D
5 ~/ a  w9 {4 u* ^8 h
发票状态，有效发票为a，作废发票为b
0 q% D5 H, F2 s5 \5 ^5 z

我将销项和进项所有数据，以公司代码为区别，提取到了不同的sheet当中，对于该公司有效发票数，作废发票数，负数发票数，方便对数据观察。
! L$ d" m( R# s9 |
. M$ g: ^& A, F7 X. J
/ J# X; C0 C! t# B# 遍历所有sheet数据
for xsn in sn.sheet_names[1:]:
    # 读取文件
" M) {- o; l+ j    datas = pd.read_excel(file_pos, sheet_name=xsn)
    datas['date']=pd.to_datetime(datas['date'],format='%Y/%m/%d')
    datas.set_index('date', drop=True)
5 g+ l; f' u6 N0 _' _* C. e9 r    # 找到全部公司名称代号
+ P$ a4 ?! j2 t% s. c0 C  r    code_list = list(set(list((datas['code']))))
( I, H# M1 W0 E$ s( Z' v) f    for name in code_list:
9 l$ s0 y& D; l$ s& g& b: d        tmp_datas = datas[datas['code'] == name]
        tmp_datas.index = range(len(tmp_datas))
& E: E8 r& t9 M: k8 `8 \        # 转换日期未object类型
/ V/ r0 F# {! N5 M6 l        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):
& E2 |8 R' w, y+ M" G& k8 L% y7 c            tmp_datas['b_count'] = list(count1)[1]
        tmp2 = tmp_datas[tmp_datas['cost'] < 0]
        tmp_datas['neg_value_tax'] = len(tmp2)
        if xsn == sn.sheet_names[1]:
            tmp_datas.to_excel(writer1,sheet_name=name,index=False)
        else:
            tmp_datas.to_excel(writer2,sheet_name=name,index=False)
负数发票：在之前购买的物品，并开具了相关正向发票，后来退货所以开具了值为负数的发票，抵消前面正数发票的值。
# g' p8 H7 a$ k# Y5 H, j) d2 ^
8 s, a" G' P9 M  g3 e2 N; v
提取到信息：
/ {' V! g0 B/ m) E3 N" I" u' ^9 O

8 Q! t0 N7 L8 V6 S! M1 y8 A' b部分公司数据记录很少，或者时间跨度大，需要综合数据指标，抵消数据数量和跨度大的影响
有些负数发票，在之前找不到对应的正数发票，可能是因为在数据记录日期之前购买的，在之后退款，因此在附件中找不到记录。
' e/ F" q- p. N! d1 R! B问题一
% [0 |$ X* _3 a5 v! |6 w: r建立指标
' _: O' A6 j* E$ y4 _$ t& K进项发票作废率，进项负数发票率，进项每月平均交易额，进项每月交易次数，
5 v$ M* I/ g3 z8 K

销项发票作废率，销项负数发票率，销项每月平均交易额，销项每月交易次数，销售收入增长率

! I2 l& q" V6 K& q' D) x
提取出相关指标到附件

$ M* R/ o! J- V, [) h
for xsn in sn.sheet_names[1:]:
& m7 @$ w5 s8 H5 N7 n8 I    # 读取文件
8 m' n4 e( e' V3 j3 y9 m1 `2 H    datas = pd.read_excel(file_pos, sheet_name=xsn)
5 Z/ p" J! m7 l    code_list = list(set(list((datas['code']))))
! ?4 \3 M) `4 l5 \1 [    for name in code_list:
        tmp_datas = datas[datas['code'] == name]
        tmp_datas.index = range(len(tmp_datas))
6 c( m# ]2 A9 a2 V- R4 T        insert_datas.append(name)
        # 作废数
        cacel_count = len(tmp_datas[tmp_datas['tax_status'] == 'b'])
' a( s/ q2 {' {* O/ l        # 有效数
        valid_count = len(tmp_datas[tmp_datas['tax_status'] == 'a'])
        # 发票作废率
' x4 e) Z& D- Z7 K        count1 = (cacel_count / (cacel_count + valid_count))*100
        # 负数发票数
        neg_count = len(tmp_datas[tmp_datas['cost'] < 0])
        # 负数发票率
8 o9 Y! ?) @# O- e2 R9 x$ N        count2 = (neg_count / valid_count) * 100
        # 转换时间
        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
( w1 e& T0 }2 Q' l        # 时间最大值
        max_time = tmp_datas.iloc[0:,1].max()
$ _# q) E4 h! o5 x2 T2 B        # 时间最小值
7 h, m/ ~3 ^  T        min_time = tmp_datas.iloc[0:,1].min()
        # 时间差
        diff_time = months(max_time, min_time) + 1
2 I8 S' ]9 J2 i# }; Z# ~1 |5 t        # 有效票
& Y' |. E$ v8 T$ o) L" g6 t        valid_tax = tmp_datas[tmp_datas['tax_status'] == 'a']
8 b' Z" W. p8 W% Y5 W3 G        # 平均月交易额
1 g7 u! ?/ F7 X; X' P        avg_money = valid_tax['totle_cost'].sum() / diff_time
# s6 T; \9 S/ ~8 z( p8 b        # 平均每月交易次数
$ R2 _+ m& O# @1 a5 D6 m1 D0 m        trans_count = len(tmp_datas) / diff_time
" W5 q  `6 _2 ]        insert_datas += [count1, count2, avg_money, trans_count,]
% y  c( c! J& q; Q" Z  R2 ]+ K        if flag:
' N; ?- G/ y/ n- _& Q3 H            df1.loc[len(df1)] = insert_datas
            df1.to_excel(writer1,sheet_name='进项信息',index=False)
" _/ Y/ A7 c/ W( R        else:
            merge_time = tmp_datas.groupby(tmp_datas['date']).sum()
# s! P7 A% D; x- o            # 销售收入增长率
1 x: [: T  m: i& F/ \1 B            income_info= list((merge_time['cost'] - merge_time['cost'].shift(1)).fillna(1))
+ p. H! l3 ^2 j& l6 i+ |            diff_time_day = days(max_time,min_time)
, Q% k# r6 h# ~$ a8 K2 Q. x            income_tax = (sum(income_info) / diff_time_day)*100
            insert_datas.append(income_tax)
1 @+ k& d( i& L- c            df2.loc[len(df2)] = insert_datas
            df2.to_excel(writer1,sheet_name='销项信息',index=False)
$ }( D6 s" B/ x0 Q. k; ^. a        insert_datas = []
    flag = False




  L7 C0 L, c+ Q" v' f并将是否违约插入到最后一列


# 提取是否违约的列表
m = []
6 l" H# w% g2 [& j4 qfor name in code_list:
3 t5 P1 Z  ?  v1 Z  c    m.append(datas[datas['code']==name]['break_contract'].tolist()[0])
    df1.loc[:,len(df1)] = m
. b1 s2 H6 f, w7 b1 a- l; Y    df1.to_excel(writer3,sheet_name='sheet1',index=False)
0 S( j% `/ N+ C- A" {6 L建立模型
0 T% m: ^9 l/ n, y' R1 @Logistics违约率预测模型
使用Logistics违约预测模型，代入所有的指标数据为自变量，是否违约为因变量，预测出违约率。


X=datas[['进项发票作废率','进项负数发票率','进项每月平均交易额','进项每月交易次数','销项发票作废率','销项负数发票率','销项每月平均交易额','销项每月交易次数','销售收入增长率']]
) r$ s! _# a/ x4 Q* m1 {y=datas['是否违约']
; K& `/ e6 A( QX_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)
! I+ |, {) N1 ~+ o- E4 fa=model.predict_proba(X_validation)
) h; u0 M1 k* W' Z: \result=[]
for i in range(len(a)):
if a[1]>0.5:
0 K& D, c5 b  w  }: W; q7 \    result.append(1)
else:
% @. |  w9 H- Z4 S( z" A* |    result.append(0)
8 o( H; o" q9 A$ Zfrom sklearn import metrics
print('误差: %.4f' % (1-metrics.recall_score(y_validation,result,average='weighted')))
最终得到一张我们的分析表格






通过预测是否违约，我们就能解决贷不贷款的问题。


* t7 H6 m  S/ d6 F& M贷款金额
, A; j  x4 D& m9 Y/ A# P贷款金额的确认，根据该公司不违约率在所有公司中的权重，乘以总贷款金额确认：
1 H* L3 K" `; V8 U! b
9 a$ V6 b( H8 Y6 G/ t
2 u/ ^3 b/ E% @8 \& s3 ?

; G2 h7 i$ |* f2 Q
) T" J$ A1 {# R4 O& }1 ~% cri=1−Zi∑123j=1(1−Zi)×M
; Q9 e. a3 ?& Fri=1−Zi∑j=1123(1−Zi)×M
  F: y; E5 x! u5 x% b因此，我们得到的贷款金额是违约率和贷款总金额组成的关系式，这在第二问中能起到重要作用。

9 l; k3 s) X. m, M, i
贷款年利率
绘制出年利率与客户流失率图，可以分析出两者应该是有关系的。利用SPSS拟合出不同信誉等级，年利率与客户流失率的关系式。


$ T, x0 E% H! V信誉等级        R平方        关系式
A        0.9977        y = 37.97x^3-258.57x^2+640.944*x -1.121
B        0.9982        y = 33.995x^3-225.051x^2+552.829*x-1.017
3 _/ `( r" h7 F( l% [2 M( L2 iC        0.9982        y = 32.157x^3-207.386x^2+504.717*x-0.973
  I2 n" [* t& V. F# x( M" o! i银行获利=贷款金额x贷款年利率x(1-利率对于信誉评级客户流失率)


在贷款金额确认，贷款年利率范围在0.4~1.5的情况下，利用上面拟合的关系式，我们能够暴力跑出最优年利率。
  ]/ Z9 }$ {  ~
3 }  X+ S1 X8 \2 Q0 I7 ?0 ^8 u5 H
, {1 u7 J7 c, ?0 Rdouble turnover_rate(double x, char ch) {
        double y = 0, result = 0;
9 [; U  B4 A! r/ Z) U7 ?( ]2 ]        switch (ch) {
        case 'A':
, J$ n" D. m4 R                y = 37.969520 * pow(x, 3) - 258.570452 * pow(x, 2) + 640.944427 * x - 1.121484;
" l: j" d3 f! `% k/ L) B                result = x * (1 - y / 100.0);
: U. M; w% j7 v& r1 P6 p3 E% ~                break;
        case 'B':
6 j6 M$ C' M# J/ [, c! D. i                y = 33.994698 * pow(x, 3) - 225.050538 * pow(x, 2) + 552.829151 * x - 1.016503;
2 h) s& _0 d' }* h2 Y                result = x * (1 - y / 100.0);
$ L. Q, K! `0 [2 u9 |                break;
        case 'C':
4 d2 g( u+ J9 G                y = 32.156864 * pow(x, 3) - 207.385880 * pow(x, 2) + 504.716993 * x - 0.973497;
' }5 f4 I6 E6 i$ D0 x! x. a                result = x * (1 - y / 100.0);
                break;
0 Z" @8 F! U" |; V        default:
                cout << "输出有误！"  << ch << endl;
# x! a; d6 L* J0 F$ d  @, j        }
0 Z* C! m* Q0 A  k* V% S
4 `/ L, {% b% a. x        return result;
}
  D& Q4 i+ C" S4 {

; {, y0 ^4 D, k

/ m. c% i% U8 N% ]* r问题二
" {' Y- _. n# b, Y6 t利用代码，重新计算出各指标数据
9 G9 b$ j( b* `/ I$ ?代入Logistics违约率预测模型，预测出各公司的违约率
根据标准普尔评级建立，主标尺，对不同违约率进行A~D等级划分，信誉等级D不予贷款
将违约率代入，之前得到的公式，得到具体贷款金额
最优年利率沿用上一问
# g; C5 c) f& Z+ H1 x# 信用等级
9 d5 G: a6 x* ^+ _cs = []
0 N! f: Q" h; I. I. `5 P" T7 q# 最优年利率,客户流失率,利率值
tax = []
" [/ R2 H3 C# J# O$ |for i in m:
3 D0 O, L% ]$ ~! k! t    if i <= 0.0069264:
- R3 L  ]1 l% z9 k) Q( D        cs.append('A')
        tax.append([0.083,0.503173,0.0412366])
! |; v. a. o8 b9 v  J    elif i > 0.0069264 and i <= 0.22619:
        cs.append('B')
1 r7 f* t8 N/ ^9 J! _        tax.append([0.097,0.505215,0.0479942])
    elif i > 0.22619 and i <= 0.509915:
/ |2 `8 @+ w+ |9 ?8 }( E        cs.append('C')
        tax.append([0.1069,0.506501,0.052755])
5 G( Z, ^. T. Y7 S    elif i > 0.509915:
! t. ]/ X- Q  r        cs.append('D')
$ Y1 H) }$ x) _# w        tax.append([0.15,0,0])
    else:
' D6 n) l7 u2 `' ~8 g1 B        print('违规')
        
parr = []
" W4 q5 `. d! x- t) Zfor arr in list(a):
$ ^! W" u; e' y* e, o    parr.append(list(arr)[0])
4 |+ T1 U1 V9 x# X0 A7 d* Nsum_val = sum(parr)
" F- ?, ?2 C1 R6 j' s) x9 |amount = []
% X* {3 n( \  T3 D& yfor ival in parr:
8 v7 M7 M5 S- P( f) j8 G3 i    tmp = ival / sum_val * 100000000
0 ~: k( t9 l! T2 v2 ^1 ^, _    if ival < 1 - 0.509915:
' p% K5 B. a( T: _# L            amount.append(0)
    else:
            amount.append(tmp)


+ r) z& T, y, u, g
/ f7 i0 N* a6 `; m+ o( A) ~: w
. Z7 j. h2 k; O可以看到，贷款金额也都在10w~100w之内。
6 f+ Y8 m* J5 `8 D

5 W. I3 P" b' v+ |/ o* K问题三
这一问，我们做得有些匆忙了，有其他想法的可以按照自己的想法做做，这里只拿我们的做参照。

- T6 W& A5 R' d3 C# B6 j4 q3 |% a
7 B7 [3 x5 }9 A' k7 h3 N! C疫情对公司影响最大的就是每月平均销售额和每月平均销售数量，因此，
7 m: c/ n# N) f+ v5 G! q0 `
3 S; @( C: S3 R# R* |
5 P; b0 D% P. B# @* k对每月平均销售额和每月平均销售数量，分别取随机数，数量取10w组，其他指标数据值不变，是否违约数据根据第二问结果，设为初始值
; K6 n6 I' ]% ?$ `: o5 n* o4 M. A; ~代入Logistics模型中，预测出每一组的违约率
1 N$ W* w  W1 h, F) b- V. E1 ^9 `判断每个公司违约率的变化情况，根据变化情况来增/减贷款金额和年利率。
4 J1 z2 X- t+ E; [- L) |( r' Krand_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':[]}
df = DataFrame(data_form)

) N: |0 a; Z. p$ C. Y+ Sfor val in rand_num:
# j0 L$ `$ c6 z2 b9 L+ Z        a = rate_func(val)
        rate_list = get_rate(a)
        df.loc[len(df)]  = rate_list
6 e5 R) B" q) r) @* C: u        df.to_excel(writer,sheet_name='违约率变化',index=False)
销售数量和销售金额的随机数范围是：0~MAX

9 F8 m& S6 T9 o& k( v) h8 `7 R5 L


————————————————
% S5 H$ }/ ]$ c$ Y版权声明：本文为CSDN博主「Hk_Mayfly」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/qq_39542714/article/details/108592087


& `" V' z0 {/ j% `3 h* O$ a
* u6 n1 G- M: I: H/ W