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

杨利霞

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


前言



2 t1 T: t- d* S1 w$ Z
) P7 C: n9 [# h2 j又是一年数据挖掘题型，第一次接触这种题型还是在去年的mathorcup上，这种题的难度就在于指标的建立和数据的处理上。后面会出一份关于数据挖掘题型，我的相关经验，常用的工具和代码。
2 W- ^4 Q% M3 l4 c

1 C" w+ N4 W8 b9 X  M6 y下面的一，二问实际都在解决
) \2 i0 Q0 F9 h$ g8 c
& N2 L# t4 f. ]
贷不贷款？
贷款金额多少？
, x0 H9 n( ^! D2 r% W数据清洗
这道题的附件数据没有出现缺省或者异常数据，因此对于数据的预处理，更多的是根据问题的需求来做的。

# J0 Z3 C+ o* b. ]& H3 z
将是否违约，违约设置为1，不违约设置为0
4 b# Z! G7 r; x# r

信誉等级ABCD分别对应4,3,2,1

6 Z* {8 T6 v$ A1 }
* m4 s* `- x- `  G& F发票状态，有效发票为a，作废发票为b

4 o# l( u: a2 K  Q
我将销项和进项所有数据，以公司代码为区别，提取到了不同的sheet当中，对于该公司有效发票数，作废发票数，负数发票数，方便对数据观察。
# V9 N$ r  d3 J6 w& ^- O9 q1 i
1 ]) K& `) J6 ]+ c
# 遍历所有sheet数据
for xsn in sn.sheet_names[1:]:
0 t, Z( ^! j# l% w5 Y    # 读取文件
& }7 @/ K$ k4 |    datas = pd.read_excel(file_pos, sheet_name=xsn)
9 u* {5 y% h- t/ q% V& N    datas['date']=pd.to_datetime(datas['date'],format='%Y/%m/%d')
" E# F5 s/ m3 I6 h/ Y5 W7 O    datas.set_index('date', drop=True)
    # 找到全部公司名称代号
    code_list = list(set(list((datas['code']))))
    for name in code_list:
        tmp_datas = datas[datas['code'] == name]
: m$ v' v) D' S6 x8 c        tmp_datas.index = range(len(tmp_datas))
        # 转换日期未object类型
* s6 W5 }" u, b+ a        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
# B" o; v; Z7 \# R) [/ w( U4 f  k        count1 = tmp_datas['tax_status'].value_counts()
7 [8 ?" U$ [3 J+ b" @. _        tmp_datas['a_count'] = list(count1)[0]
        if(len(count1) > 1):
, N  c" C* {$ [7 |% U" r9 C            tmp_datas['b_count'] = list(count1)[1]
        tmp2 = tmp_datas[tmp_datas['cost'] < 0]
8 X9 y6 p0 S9 b2 N0 P        tmp_datas['neg_value_tax'] = len(tmp2)
, W" _" {3 n# Q2 w5 @/ ~        if xsn == sn.sheet_names[1]:
7 C% ^* v  }4 l6 z/ \3 O1 [& g" H6 ]8 x            tmp_datas.to_excel(writer1,sheet_name=name,index=False)
        else:
            tmp_datas.to_excel(writer2,sheet_name=name,index=False)
* J; a  i& o0 I( g负数发票：在之前购买的物品，并开具了相关正向发票，后来退货所以开具了值为负数的发票，抵消前面正数发票的值。
) Z& ~+ K" B2 C. ^

! \3 r+ V* R( a4 |提取到信息：
$ n; T+ D7 x! q% O6 W/ b$ x
& j5 t* \9 k) j6 [+ s
9 ~/ n$ x$ X0 x! L' T& D: q部分公司数据记录很少，或者时间跨度大，需要综合数据指标，抵消数据数量和跨度大的影响
. |' a5 p/ W  {9 p9 R# w% h) p有些负数发票，在之前找不到对应的正数发票，可能是因为在数据记录日期之前购买的，在之后退款，因此在附件中找不到记录。
* W. k; ]9 ]- n* F  s$ _0 Y2 \) x问题一
  k, a3 q" I' m" r$ O建立指标
; R, z: A3 e( Y/ L% K+ q9 ~进项发票作废率，进项负数发票率，进项每月平均交易额，进项每月交易次数，
, J5 l7 B$ r' [' H9 Y- l

, S) c9 _2 D8 s; J销项发票作废率，销项负数发票率，销项每月平均交易额，销项每月交易次数，销售收入增长率


提取出相关指标到附件
# {4 c& E. U. v1 r, W
. w; L  z* O5 v; e+ l; r
for xsn in sn.sheet_names[1:]:
    # 读取文件
6 _1 s' Z( e. N- d4 [! q7 t    datas = pd.read_excel(file_pos, sheet_name=xsn)
    code_list = list(set(list((datas['code']))))
    for name in code_list:
& f( y2 V$ {) q        tmp_datas = datas[datas['code'] == name]
        tmp_datas.index = range(len(tmp_datas))
  m, L% Z8 c7 H4 j        insert_datas.append(name)
/ W( w; r" X, h( u" ~" K; A: N; J        # 作废数
4 z* ~% d3 c4 t/ w6 N4 P        cacel_count = len(tmp_datas[tmp_datas['tax_status'] == 'b'])
. q, V7 S2 @6 O3 S% T" W0 V8 a/ ?        # 有效数
        valid_count = len(tmp_datas[tmp_datas['tax_status'] == 'a'])
        # 发票作废率
0 |; R7 U4 }# s! q" g        count1 = (cacel_count / (cacel_count + valid_count))*100
        # 负数发票数
- _" z$ t. \* U7 _. Z/ v( T) A        neg_count = len(tmp_datas[tmp_datas['cost'] < 0])
        # 负数发票率
        count2 = (neg_count / valid_count) * 100
        # 转换时间
6 K6 g6 q& _" s3 d2 Z; [& a        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
        # 时间最大值
        max_time = tmp_datas.iloc[0:,1].max()
        # 时间最小值
) d6 W9 C) c/ R- j  P        min_time = tmp_datas.iloc[0:,1].min()
6 N! U& ?3 x* o/ m4 P        # 时间差
        diff_time = months(max_time, min_time) + 1
        # 有效票
        valid_tax = tmp_datas[tmp_datas['tax_status'] == 'a']
1 K5 C1 d+ S; t3 l) [! u# Q7 r; u1 }        # 平均月交易额
$ a) G8 E4 X9 i7 \; w  X        avg_money = valid_tax['totle_cost'].sum() / diff_time
( X3 w8 Z% U% ^* G) U        # 平均每月交易次数
& N, y" b% a  c2 W* i8 Q! g1 a        trans_count = len(tmp_datas) / diff_time
        insert_datas += [count1, count2, avg_money, trans_count,]
# x; q5 b* A$ D; F  p; E  [* e- F' s) S4 O        if flag:
" B$ p6 b5 P1 e6 H: m, m            df1.loc[len(df1)] = insert_datas
  j3 I/ Z3 s+ Q7 `& ]. n            df1.to_excel(writer1,sheet_name='进项信息',index=False)
$ |& k8 X. P: p4 B* U( Z2 ^3 Q        else:
4 `& \* n; \" }0 A; F0 O' a% j            merge_time = tmp_datas.groupby(tmp_datas['date']).sum()
3 U1 R  r6 ^9 `% J1 f            # 销售收入增长率
- T* S0 E/ e0 C% d            income_info= list((merge_time['cost'] - merge_time['cost'].shift(1)).fillna(1))
            diff_time_day = days(max_time,min_time)
/ r! _+ \8 m* E: \) r; 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 = []
& T6 B: O. O2 M, F1 R! S. |    flag = False
7 ?6 ]: U8 h0 }8 S: M

( e) a8 r: ]8 m9 ^0 B5 R5 h
( K+ Y" |! W' [3 f
并将是否违约插入到最后一列

7 R+ W+ q8 f* A! h+ L
# 提取是否违约的列表
m = []
. v* c7 D3 z' e, f2 U# Mfor name in code_list:
5 j4 L) ~6 K! b, H4 {+ `8 H4 e: J- @    m.append(datas[datas['code']==name]['break_contract'].tolist()[0])
! P8 D6 g6 m! Y# T* ~: h    df1.loc[:,len(df1)] = m
    df1.to_excel(writer3,sheet_name='sheet1',index=False)
建立模型
  E7 L8 }% x; `: H& V8 R3 K: b6 iLogistics违约率预测模型
使用Logistics违约预测模型，代入所有的指标数据为自变量，是否违约为因变量，预测出违约率。

0 }, o) X# @2 k5 }
0 P8 N7 m% h  T6 m' a- n! yX=datas[['进项发票作废率','进项负数发票率','进项每月平均交易额','进项每月交易次数','销项发票作废率','销项负数发票率','销项每月平均交易额','销项每月交易次数','销售收入增长率']]
& Z" u/ x8 l& b3 o' {y=datas['是否违约']
X_train, X_test, y_train, y_test = train_test_split(X,y, test_size=0.2,random_state=2020)
! U& c9 _# `. Z% C* rX_validation, X_test, y_validation, y_test = train_test_split(X_test,y_test, test_size=0.1,random_state=2020)
model = LogisticRegression()
* \, I6 X7 i1 G/ @model.fit(X_train,y_train)
% Z8 K+ F* Y. F2 j% `: Va=model.predict_proba(X_validation)
result=[]
0 t! f& J  |/ c2 Z. g, Jfor i in range(len(a)):
if a[1]>0.5:
    result.append(1)
else:
    result.append(0)
/ M/ F! U2 j# d) vfrom sklearn import metrics
. @! T6 {6 K7 ]  Jprint('误差: %.4f' % (1-metrics.recall_score(y_validation,result,average='weighted')))
8 Z; j+ H* K* H, p  U4 I2 R最终得到一张我们的分析表格
6 q2 M* a2 R8 ~4 p( w: M. Z

3 ?9 L/ G4 C4 w  a# m


; ?# T4 x$ S4 I. \
, E$ H! m3 J# r- r. ~) r( N通过预测是否违约，我们就能解决贷不贷款的问题。

/ {+ `6 H% F  h/ ?
贷款金额
贷款金额的确认，根据该公司不违约率在所有公司中的权重，乘以总贷款金额确认：




, v5 q9 c& V2 Y
( m1 Z0 h3 {. n. q! C$ \ri=1−Zi∑123j=1(1−Zi)×M
ri=1−Zi∑j=1123(1−Zi)×M
因此，我们得到的贷款金额是违约率和贷款总金额组成的关系式，这在第二问中能起到重要作用。

9 d# P: T* {3 @+ E# N
贷款年利率
6 l' I  r# k+ @8 u, `4 }绘制出年利率与客户流失率图，可以分析出两者应该是有关系的。利用SPSS拟合出不同信誉等级，年利率与客户流失率的关系式。

7 N' P3 Z: U2 ]% E
信誉等级        R平方        关系式
2 I& {" ]7 V1 n/ sA        0.9977        y = 37.97x^3-258.57x^2+640.944*x -1.121
2 Y& V; t; f/ F7 N9 l8 f4 |B        0.9982        y = 33.995x^3-225.051x^2+552.829*x-1.017
/ ^, A# E2 g' O& H$ C* P+ U0 }C        0.9982        y = 32.157x^3-207.386x^2+504.717*x-0.973
银行获利=贷款金额x贷款年利率x(1-利率对于信誉评级客户流失率)
& R9 \' W7 w6 n5 E- o3 _
9 S7 M' h6 ]. T" {8 p3 x, ?) O( V
: Y* O" N2 A8 I: p在贷款金额确认，贷款年利率范围在0.4~1.5的情况下，利用上面拟合的关系式，我们能够暴力跑出最优年利率。
( `7 F: E- y. B5 N1 M% [

5 K4 o, h3 }4 m& |! K/ ]0 Qdouble turnover_rate(double x, char ch) {
        double y = 0, result = 0;
        switch (ch) {
6 ?  i% h: I  D  j" W        case 'A':
- U, N0 e7 u( _9 c: I                y = 37.969520 * pow(x, 3) - 258.570452 * pow(x, 2) + 640.944427 * x - 1.121484;
                result = x * (1 - y / 100.0);
                break;
' X4 h* Y' ^8 A, Y        case 'B':
                y = 33.994698 * pow(x, 3) - 225.050538 * pow(x, 2) + 552.829151 * x - 1.016503;
                result = x * (1 - y / 100.0);
                break;
        case 'C':
% d- ]+ r* X/ b+ K7 K                y = 32.156864 * pow(x, 3) - 207.385880 * pow(x, 2) + 504.716993 * x - 0.973497;
                result = x * (1 - y / 100.0);
                break;
0 ?. O) N  s1 T9 ^+ g4 ^' D        default:
                cout << "输出有误！"  << ch << endl;
        }

        return result;
4 m: |' L9 b! I}

0 a8 \1 M" [. [4 c+ l" v


问题二
利用代码，重新计算出各指标数据
代入Logistics违约率预测模型，预测出各公司的违约率
4 i7 W6 l/ D/ f; l( y; W' L$ h0 |根据标准普尔评级建立，主标尺，对不同违约率进行A~D等级划分，信誉等级D不予贷款
8 P, y- ]# {( O6 e" e将违约率代入，之前得到的公式，得到具体贷款金额
最优年利率沿用上一问
# 信用等级
cs = []
# 最优年利率,客户流失率,利率值
4 o2 ?  Q/ _. Q& w8 Qtax = []
3 ], O2 H1 o/ T+ S$ Dfor i in m:
0 H& i9 ]/ t/ `/ z5 ~    if i <= 0.0069264:
6 s3 L" `: I! H" \5 Y! ]: H  I5 Y        cs.append('A')
0 y6 l  p( q$ o' L        tax.append([0.083,0.503173,0.0412366])
    elif i > 0.0069264 and i <= 0.22619:
        cs.append('B')
        tax.append([0.097,0.505215,0.0479942])
    elif i > 0.22619 and i <= 0.509915:
        cs.append('C')
2 z! k+ Z" c* H5 N* W0 h( v        tax.append([0.1069,0.506501,0.052755])
2 K; t1 M5 j; ?    elif i > 0.509915:
        cs.append('D')
        tax.append([0.15,0,0])
2 e! s# e2 v: \8 M, q    else:
1 L; J( X, y4 ~. v% r        print('违规')
  N9 f" E0 p0 W1 r/ ]        
parr = []
for arr in list(a):
* Y7 @. h  _- I/ h+ r0 d    parr.append(list(arr)[0])
sum_val = sum(parr)
: s& ^7 m4 N1 V  N+ I3 Wamount = []
for ival in parr:
    tmp = ival / sum_val * 100000000
* X9 D" ?! N2 v7 x8 B    if ival < 1 - 0.509915:
            amount.append(0)
    else:
' V8 Q3 S/ ]9 q: [3 |$ r9 {8 F            amount.append(tmp)

0 R, Z0 v7 i" m; n* {
- J( @+ i  V, A* C# J9 X, |* L
! C. s! @( \- p4 ?  m/ B$ {
可以看到，贷款金额也都在10w~100w之内。
! P) A, p- H% L8 r- I
$ B! U  q$ Y7 ]2 W" g$ b
问题三
这一问，我们做得有些匆忙了，有其他想法的可以按照自己的想法做做，这里只拿我们的做参照。
; y2 t/ s4 `4 g& s& A4 L  k5 K8 L& A

1 E+ B. M4 r8 K2 t6 P6 P疫情对公司影响最大的就是每月平均销售额和每月平均销售数量，因此，
. g3 ?! z, Z: g! ?, h$ ~

对每月平均销售额和每月平均销售数量，分别取随机数，数量取10w组，其他指标数据值不变，是否违约数据根据第二问结果，设为初始值
代入Logistics模型中，预测出每一组的违约率
6 b% P3 d7 u  w7 ~- V6 w# J# m判断每个公司违约率的变化情况，根据变化情况来增/减贷款金额和年利率。
% j0 D1 W5 N6 z8 t, c5 n4 Yrand_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)

for val in rand_num:
/ X. {: ~" P5 T        a = rate_func(val)
        rate_list = get_rate(a)
7 J( C; I' ~5 }( |# o. j/ s        df.loc[len(df)]  = rate_list
        df.to_excel(writer,sheet_name='违约率变化',index=False)
& g8 |$ \( O+ c销售数量和销售金额的随机数范围是：0~MAX

7 Z+ |: ?+ A" x# ^/ x& Y. ]

* p, k+ ?5 g# t3 G# e2 T$ t7 o
————————————————
版权声明：本文为CSDN博主「Hk_Mayfly」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/qq_39542714/article/details/108592087


3 i9 x& `" t. Z
* E9 @9 }! E* j