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

杨利霞

2020 全国大学生数学建模竞赛C题思路+代码
# [1 f1 ~) S6 h' X0 `) Z- ]" E1 |题目链接：https://cloud.189.cn/t/ri2uUb7BRVJr


. \! E* K) W5 |$ t6 W& k前言



3 D: g, z! c* v  w* d* P
, p( z4 D# ^" Q, u又是一年数据挖掘题型，第一次接触这种题型还是在去年的mathorcup上，这种题的难度就在于指标的建立和数据的处理上。后面会出一份关于数据挖掘题型，我的相关经验，常用的工具和代码。

4 C4 q3 }1 k7 R$ I( Y6 a1 B  i
下面的一，二问实际都在解决
5 N, ^2 J' j. d1 p8 H

6 y1 F5 a: N0 `6 d. f  ^贷不贷款？
贷款金额多少？
数据清洗
这道题的附件数据没有出现缺省或者异常数据，因此对于数据的预处理，更多的是根据问题的需求来做的。

2 s+ ?$ I, p1 ~2 |+ q: W
将是否违约，违约设置为1，不违约设置为0
! X  S7 ?7 s8 D, C- P( K) G0 C

' ^) ^; M- A0 F& Q: H9 I# i1 n0 \信誉等级ABCD分别对应4,3,2,1
4 Q6 [4 n# ~& g7 q4 F. u8 Z4 a

发票状态，有效发票为a，作废发票为b
2 \/ g9 P) M1 g- }& w; k
: f3 J3 b* Z6 d7 z" c! U  w2 X& n
我将销项和进项所有数据，以公司代码为区别，提取到了不同的sheet当中，对于该公司有效发票数，作废发票数，负数发票数，方便对数据观察。


# 遍历所有sheet数据
3 _: U7 ^1 k3 ]* @7 d6 E6 ^  }. m" Ufor xsn in sn.sheet_names[1:]:
# v" k6 f) ]! G    # 读取文件
    datas = pd.read_excel(file_pos, sheet_name=xsn)
: M! `+ C+ Q# w; k6 z" w; v    datas['date']=pd.to_datetime(datas['date'],format='%Y/%m/%d')
( n. R2 l# t6 G# x% U! J* t    datas.set_index('date', drop=True)
    # 找到全部公司名称代号
$ u2 W5 E; E" s    code_list = list(set(list((datas['code']))))
    for name in code_list:
        tmp_datas = datas[datas['code'] == name]
        tmp_datas.index = range(len(tmp_datas))
  E. e, f% r* s6 a2 d: r        # 转换日期未object类型
        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
: Y. o7 D! Q6 E6 X) J8 H4 r        count1 = tmp_datas['tax_status'].value_counts()
        tmp_datas['a_count'] = list(count1)[0]
        if(len(count1) > 1):
' n: ^! f! q  D5 a2 g' l            tmp_datas['b_count'] = list(count1)[1]
8 a) t  ~0 k3 D) V        tmp2 = tmp_datas[tmp_datas['cost'] < 0]
        tmp_datas['neg_value_tax'] = len(tmp2)
        if xsn == sn.sheet_names[1]:
3 W  G+ ^% `; {2 K( q1 o$ L7 P            tmp_datas.to_excel(writer1,sheet_name=name,index=False)
& {& V& [. m* L$ A+ g# v+ j0 S1 b        else:
            tmp_datas.to_excel(writer2,sheet_name=name,index=False)
0 w; n* _& {) k+ ]" S; R2 O负数发票：在之前购买的物品，并开具了相关正向发票，后来退货所以开具了值为负数的发票，抵消前面正数发票的值。


提取到信息：

3 c# p* b" M3 P8 z0 w/ L$ i+ b
: l+ M6 y. y: Q( _部分公司数据记录很少，或者时间跨度大，需要综合数据指标，抵消数据数量和跨度大的影响
有些负数发票，在之前找不到对应的正数发票，可能是因为在数据记录日期之前购买的，在之后退款，因此在附件中找不到记录。
问题一
建立指标
' ^$ z, z1 x9 b5 Y! R; F7 z6 E进项发票作废率，进项负数发票率，进项每月平均交易额，进项每月交易次数，
: i3 H% U' @; @, L" @& ^

# |4 ]# n. S* \" u2 i+ _销项发票作废率，销项负数发票率，销项每月平均交易额，销项每月交易次数，销售收入增长率
9 n$ z, |6 ]% _. d- z  D) {3 s+ E
- M) p+ i& `7 S6 ^$ H
提取出相关指标到附件


for xsn in sn.sheet_names[1:]:
    # 读取文件
    datas = pd.read_excel(file_pos, sheet_name=xsn)
0 j# R2 V- }# @9 _    code_list = list(set(list((datas['code']))))
    for name in code_list:
1 Z$ U" ]  ^, K8 d' E) {' U        tmp_datas = datas[datas['code'] == name]
' C6 F5 @8 R9 r7 S) g+ s  Z        tmp_datas.index = range(len(tmp_datas))
        insert_datas.append(name)
        # 作废数
        cacel_count = len(tmp_datas[tmp_datas['tax_status'] == 'b'])
        # 有效数
* w& ^3 H8 a) y, A1 n' ^5 O! S        valid_count = len(tmp_datas[tmp_datas['tax_status'] == 'a'])
3 o7 t. `+ f" z- m        # 发票作废率
        count1 = (cacel_count / (cacel_count + valid_count))*100
7 e7 L/ R+ Z0 ~. K' z- c        # 负数发票数
7 P4 T% X8 e# C7 s/ X6 D" \6 o        neg_count = len(tmp_datas[tmp_datas['cost'] < 0])
        # 负数发票率
5 T. u- D1 o: I# a5 r        count2 = (neg_count / valid_count) * 100
        # 转换时间
6 X- d; ]" a4 i+ H. v/ a( T        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
; s0 A$ x- h# Y: M        # 时间最大值
        max_time = tmp_datas.iloc[0:,1].max()
7 x+ T$ s8 }) X! p0 \' w        # 时间最小值
        min_time = tmp_datas.iloc[0:,1].min()
, b% G9 V) e- ]9 `3 [! x  E        # 时间差
        diff_time = months(max_time, min_time) + 1
        # 有效票
6 I5 q1 w  w; y! q        valid_tax = tmp_datas[tmp_datas['tax_status'] == 'a']
        # 平均月交易额
4 U7 O- C" h5 j1 {$ g0 a- z        avg_money = valid_tax['totle_cost'].sum() / diff_time
        # 平均每月交易次数
0 M0 g0 b  Y/ B: j        trans_count = len(tmp_datas) / diff_time
        insert_datas += [count1, count2, avg_money, trans_count,]
+ l" g( X9 \! q8 _+ O7 F        if flag:
            df1.loc[len(df1)] = insert_datas
1 f) I: f/ Z, Z% E; L            df1.to_excel(writer1,sheet_name='进项信息',index=False)
* {7 K5 c  _4 ^        else:
            merge_time = tmp_datas.groupby(tmp_datas['date']).sum()
8 n$ l9 n( l7 p4 B  H! l9 R            # 销售收入增长率
$ u- u/ T8 e$ r" Q  l7 X8 \            income_info= list((merge_time['cost'] - merge_time['cost'].shift(1)).fillna(1))
            diff_time_day = days(max_time,min_time)
: p% ~5 K" [( n5 U4 g* U4 g            income_tax = (sum(income_info) / diff_time_day)*100
            insert_datas.append(income_tax)
            df2.loc[len(df2)] = insert_datas
            df2.to_excel(writer1,sheet_name='销项信息',index=False)
        insert_datas = []
    flag = False
* Z  b7 B% y3 W) P: W4 u8 S
9 q/ {9 Q" R6 l  Y* U9 p
: L! P9 ]  ?5 e  [
6 S! @" N, b  c, b
, ]6 h& l5 t3 V; U并将是否违约插入到最后一列

5 b6 ^9 H$ C* y
' _" U$ B7 X* H  O9 ~: ~: t# 提取是否违约的列表
m = []
: M/ d% }9 v9 v( d8 J/ zfor name in code_list:
) {; j* q  o! o- |, Y% `( {( F; D    m.append(datas[datas['code']==name]['break_contract'].tolist()[0])
1 J: K- @. G- p3 o1 ]" }    df1.loc[:,len(df1)] = m
    df1.to_excel(writer3,sheet_name='sheet1',index=False)
. e9 j# v* l8 z% S) d建立模型
( q+ k; q4 v0 E; l/ z9 v" ILogistics违约率预测模型
使用Logistics违约预测模型，代入所有的指标数据为自变量，是否违约为因变量，预测出违约率。

$ B1 E. S1 h3 T
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)
8 Y% i9 w% f* a! a" e( l; |4 Lmodel = LogisticRegression()
8 N$ E+ ?. ~+ G# C5 t" Pmodel.fit(X_train,y_train)
a=model.predict_proba(X_validation)
result=[]
- C" U' Z) ]0 P6 I2 |! Efor i in range(len(a)):
; K. d* k2 ~4 q+ B1 B1 k" mif a[1]>0.5:
, D$ V" e: `  Y6 H' F( `    result.append(1)
2 Y0 [5 J) x/ Nelse:
    result.append(0)
$ W, _1 v# x, O% E' ffrom sklearn import metrics
print('误差: %.4f' % (1-metrics.recall_score(y_validation,result,average='weighted')))
最终得到一张我们的分析表格




) ?. }* b1 X3 g* y! o# r, U, q

通过预测是否违约，我们就能解决贷不贷款的问题。
3 Z: _) e; p/ p) S

  ~8 p$ ^0 x  {# H5 h2 X9 a! m贷款金额
贷款金额的确认，根据该公司不违约率在所有公司中的权重，乘以总贷款金额确认：


' Y* D: U* t$ o0 _
6 A5 d1 {0 N& P9 w9 E
9 J8 a! D% ]% X2 @; ]' s6 q
6 I' i& j: h# l8 y) Lri=1−Zi∑123j=1(1−Zi)×M
ri=1−Zi∑j=1123(1−Zi)×M
1 H8 {4 J  K( C6 v因此，我们得到的贷款金额是违约率和贷款总金额组成的关系式，这在第二问中能起到重要作用。

' n: L" d* C; B: U6 B. o: A, S
贷款年利率
, x' s5 g# L& q; U1 ^4 L& N7 T绘制出年利率与客户流失率图，可以分析出两者应该是有关系的。利用SPSS拟合出不同信誉等级，年利率与客户流失率的关系式。
7 D  g& j6 R0 r
0 L; N" @* S4 P
信誉等级        R平方        关系式
A        0.9977        y = 37.97x^3-258.57x^2+640.944*x -1.121
# L$ Y# x5 \! w& o% _* ?; xB        0.9982        y = 33.995x^3-225.051x^2+552.829*x-1.017
" @: U% C  G( Z3 M  sC        0.9982        y = 32.157x^3-207.386x^2+504.717*x-0.973
  o& b" y* S1 |9 q# k银行获利=贷款金额x贷款年利率x(1-利率对于信誉评级客户流失率)
* L) }4 q: }$ f

/ l7 `. [# S' `+ b( c在贷款金额确认，贷款年利率范围在0.4~1.5的情况下，利用上面拟合的关系式，我们能够暴力跑出最优年利率。
9 S# W$ H" g, g: F' U

double turnover_rate(double x, char ch) {
        double y = 0, result = 0;
3 Z4 H: ~5 k% c9 q        switch (ch) {
- b% P0 ^& }6 u; Q        case 'A':
: w  u5 T, Z% M                y = 37.969520 * pow(x, 3) - 258.570452 * pow(x, 2) + 640.944427 * x - 1.121484;
! Z3 X7 i5 L' Z6 ^9 }! @                result = x * (1 - y / 100.0);
; n) Y- f' |' _2 d3 e: e                break;
; K6 \3 _0 O; I# k        case 'B':
" `* N: [& B* f8 m$ d7 T6 Q4 C8 r                y = 33.994698 * pow(x, 3) - 225.050538 * pow(x, 2) + 552.829151 * x - 1.016503;
                result = x * (1 - y / 100.0);
                break;
4 U- @/ q2 D( b( U* A% U0 c1 }+ j        case 'C':
                y = 32.156864 * pow(x, 3) - 207.385880 * pow(x, 2) + 504.716993 * x - 0.973497;
% E- C6 y; h. O1 J                result = x * (1 - y / 100.0);
8 T% h) I  \  U                break;
' u* ?8 K, K1 c5 L% x1 U" J% z        default:
2 s8 M/ K& {; S% z                cout << "输出有误！"  << ch << endl;
% |; g; j- Q( L* G* P        }
: g" I3 l' d, ?5 w6 L) j$ H0 ?* E
, U$ Z" |: p7 J7 n$ v4 Y% i+ t        return result;
}

  Z9 R. a' `1 _  e


; y) j3 r. K+ i6 p$ {- a问题二
利用代码，重新计算出各指标数据
代入Logistics违约率预测模型，预测出各公司的违约率
根据标准普尔评级建立，主标尺，对不同违约率进行A~D等级划分，信誉等级D不予贷款
8 a; x. E  [: V8 G将违约率代入，之前得到的公式，得到具体贷款金额
6 P( d9 ]' K* N  w3 x! Q最优年利率沿用上一问
9 f& c4 ]- s7 b# 信用等级
8 s7 Z( F, i7 ?/ S. `% Vcs = []
; ?) @% r, l: h# 最优年利率,客户流失率,利率值
tax = []
for i in m:
# I7 y' }. y6 U: b9 S& k% s    if i <= 0.0069264:
        cs.append('A')
4 j7 E& o& R( v2 V        tax.append([0.083,0.503173,0.0412366])
4 k3 Q- Y4 @3 B) {% S1 X  ~    elif i > 0.0069264 and i <= 0.22619:
        cs.append('B')
/ A/ P4 m4 ]! z3 ~. n        tax.append([0.097,0.505215,0.0479942])
0 M7 E% a$ i4 \% P0 K    elif i > 0.22619 and i <= 0.509915:
        cs.append('C')
! h$ x" U1 ?+ h; f7 R: O        tax.append([0.1069,0.506501,0.052755])
: v( N& t& Q2 I/ B' |( o. @! l6 @    elif i > 0.509915:
        cs.append('D')
7 ]. @/ K9 |+ h1 f  C) }5 N5 k4 T  \        tax.append([0.15,0,0])
# y) T* P% i" d6 h0 M7 k+ s    else:
- d8 q( b# u* p6 J; S/ F        print('违规')
, p# L6 ~2 q( T        
parr = []
for arr in list(a):
    parr.append(list(arr)[0])
sum_val = sum(parr)
; M' L% B  p" s" Ramount = []
for ival in parr:
    tmp = ival / sum_val * 100000000
    if ival < 1 - 0.509915:
! E! F+ Z( {$ ^$ c' i8 h; b            amount.append(0)
    else:
            amount.append(tmp)
: ]+ Y1 b% @+ f" O% H- Q& ~% ?
' l  r* H# Z3 x+ `

) [, J7 w* ?7 m6 m
: v, @, L- H# Y" Q可以看到，贷款金额也都在10w~100w之内。

7 E2 a9 I9 Y) f, t; Z( Q5 u
问题三
' G6 S, s$ l+ v6 _: Q8 c这一问，我们做得有些匆忙了，有其他想法的可以按照自己的想法做做，这里只拿我们的做参照。
. {7 H6 @' d  J. ^
9 n# O% k- }1 N
疫情对公司影响最大的就是每月平均销售额和每月平均销售数量，因此，

2 V: c/ }- g1 D& m& p4 d% d
对每月平均销售额和每月平均销售数量，分别取随机数，数量取10w组，其他指标数据值不变，是否违约数据根据第二问结果，设为初始值
& b! J- b" R) X0 }, Q( H/ ]9 Y代入Logistics模型中，预测出每一组的违约率
  _9 V4 `/ @8 J& f/ @判断每个公司违约率的变化情况，根据变化情况来增/减贷款金额和年利率。
rand_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)
9 y4 v& @* g( o* }3 u
for val in rand_num:
5 h+ i" c1 \7 K# i$ O& q        a = rate_func(val)
, b4 n7 n3 w- w4 X  o5 Z9 g. _9 _        rate_list = get_rate(a)
& N' k6 E3 O& P' |+ Z$ ]        df.loc[len(df)]  = rate_list
& ^" v; y8 P- H) m5 M        df.to_excel(writer,sheet_name='违约率变化',index=False)
% ?: w5 a" a/ W" x/ A! a销售数量和销售金额的随机数范围是：0~MAX


0 b- D# [6 p/ r1 [

9 |4 @2 G3 u; ]# o" \————————————————
版权声明：本文为CSDN博主「Hk_Mayfly」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
: U( l% h& R! g原文链接：https://blog.csdn.net/qq_39542714/article/details/108592087

1 m+ j2 B3 ~9 v, Y/ }% P, D. `  O! V
1 x) j% g1 F6 f% b