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

杨利霞

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

4 S5 P- U+ t# e  I
前言
" B+ J+ [3 t; y( ]

) M0 b% G" a% K8 D$ j9 x

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


6 X; T+ @$ ~* u) `/ n/ i下面的一，二问实际都在解决


: f1 y/ O/ `) ^7 B8 \5 p6 N贷不贷款？
贷款金额多少？
8 N* K3 g$ m) g6 E# v数据清洗
这道题的附件数据没有出现缺省或者异常数据，因此对于数据的预处理，更多的是根据问题的需求来做的。
2 d# V) @+ X$ X" @2 i& J2 `( O
, U) ]3 Q4 N4 N' a# f/ O1 P
将是否违约，违约设置为1，不违约设置为0

1 N4 _( M" c' N: T
7 S- d. F: G! Y信誉等级ABCD分别对应4,3,2,1
8 E0 h/ O7 f' g. Q1 O
3 }' R+ R. S( F$ A3 }0 c
9 r( s3 z* Q# ]6 Y7 r  l发票状态，有效发票为a，作废发票为b
0 Y/ N# D- Z" @7 V: Y( Z/ O; n

我将销项和进项所有数据，以公司代码为区别，提取到了不同的sheet当中，对于该公司有效发票数，作废发票数，负数发票数，方便对数据观察。


) ?. [; c2 @' n; {4 O( \# E5 \2 D# u# 遍历所有sheet数据
for xsn in sn.sheet_names[1:]:
' {# d# K$ ~# I, a    # 读取文件
4 P( r+ R' o2 x8 p    datas = pd.read_excel(file_pos, sheet_name=xsn)
    datas['date']=pd.to_datetime(datas['date'],format='%Y/%m/%d')
7 N/ v" H; M# N    datas.set_index('date', drop=True)
/ [% j  M0 ~( z$ F) J/ C% H" B5 k$ X    # 找到全部公司名称代号
' r$ ?7 ~, ~1 [/ ~9 `! V* \    code_list = list(set(list((datas['code']))))
% I6 @& t* j0 P+ g/ I4 R4 g) W5 |    for name in code_list:
        tmp_datas = datas[datas['code'] == name]
: K0 o' b* l! y1 \$ u2 o# T        tmp_datas.index = range(len(tmp_datas))
        # 转换日期未object类型
, t$ Z0 s9 Y. a9 ?4 d0 t        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
        count1 = tmp_datas['tax_status'].value_counts()
8 c. b# Y8 ~0 T8 u: h/ W4 c- v4 D" P        tmp_datas['a_count'] = list(count1)[0]
; M, l% C9 j- G  |8 a; E% @        if(len(count1) > 1):
            tmp_datas['b_count'] = list(count1)[1]
( M) m9 M: }! F, L9 g        tmp2 = tmp_datas[tmp_datas['cost'] < 0]
        tmp_datas['neg_value_tax'] = len(tmp2)
1 {, J, M+ j* m' s        if xsn == sn.sheet_names[1]:
$ e) U# u) o3 x' P            tmp_datas.to_excel(writer1,sheet_name=name,index=False)
        else:
            tmp_datas.to_excel(writer2,sheet_name=name,index=False)
负数发票：在之前购买的物品，并开具了相关正向发票，后来退货所以开具了值为负数的发票，抵消前面正数发票的值。

& v0 c, o6 c! q
提取到信息：


部分公司数据记录很少，或者时间跨度大，需要综合数据指标，抵消数据数量和跨度大的影响
有些负数发票，在之前找不到对应的正数发票，可能是因为在数据记录日期之前购买的，在之后退款，因此在附件中找不到记录。
5 f  k& h% i" p0 M问题一
$ _9 o8 f# H( j$ N# h建立指标
4 R' u- @9 `* M) d8 {进项发票作废率，进项负数发票率，进项每月平均交易额，进项每月交易次数，
: P5 E# S+ ~7 m: C: m( d* H$ l
& l1 F& o) \2 B  b  @% A
( O# B5 [: @9 }3 ]销项发票作废率，销项负数发票率，销项每月平均交易额，销项每月交易次数，销售收入增长率
  K7 ], H% E' }/ k/ K$ ?
% F6 F5 _' M' r! o. J
提取出相关指标到附件

* [, b/ [$ h! f. @
for xsn in sn.sheet_names[1:]:
    # 读取文件
    datas = pd.read_excel(file_pos, sheet_name=xsn)
$ S! a; ]9 g: j: A    code_list = list(set(list((datas['code']))))
; b% C, m0 A1 X8 X    for name in code_list:
        tmp_datas = datas[datas['code'] == name]
        tmp_datas.index = range(len(tmp_datas))
        insert_datas.append(name)
        # 作废数
        cacel_count = len(tmp_datas[tmp_datas['tax_status'] == 'b'])
3 ]& M+ F4 Z7 \( P7 S6 N        # 有效数
        valid_count = len(tmp_datas[tmp_datas['tax_status'] == 'a'])
        # 发票作废率
        count1 = (cacel_count / (cacel_count + valid_count))*100
7 P1 m: h7 t. Z" u( v        # 负数发票数
* M9 `: L8 C# w3 C; C: j" z# j        neg_count = len(tmp_datas[tmp_datas['cost'] < 0])
        # 负数发票率
' }+ q' b& w; i  x! M4 v% H        count2 = (neg_count / valid_count) * 100
4 i7 f" d, v4 u        # 转换时间
/ y. V. {) w9 ^        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
; S+ h7 A  F9 B, w' i        # 时间最大值
        max_time = tmp_datas.iloc[0:,1].max()
4 A+ N/ N+ x' L9 r6 P$ n        # 时间最小值
        min_time = tmp_datas.iloc[0:,1].min()
        # 时间差
0 a& b7 Q0 z: I% m  j6 W, N        diff_time = months(max_time, min_time) + 1
        # 有效票
        valid_tax = tmp_datas[tmp_datas['tax_status'] == 'a']
$ @4 V0 ^) ]% @+ Z        # 平均月交易额
( H0 V8 x! {) ?7 g        avg_money = valid_tax['totle_cost'].sum() / diff_time
        # 平均每月交易次数
        trans_count = len(tmp_datas) / diff_time
        insert_datas += [count1, count2, avg_money, trans_count,]
' _  p0 F# Y  h        if flag:
+ w4 C9 ?4 b. s; ~7 Q            df1.loc[len(df1)] = insert_datas
6 o+ w/ \$ r( n4 h  [2 F! n2 P9 V            df1.to_excel(writer1,sheet_name='进项信息',index=False)
5 p, v: n6 n- j, Z$ t+ j        else:
8 J: V: L# A# f! J            merge_time = tmp_datas.groupby(tmp_datas['date']).sum()
            # 销售收入增长率
            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
            insert_datas.append(income_tax)
" D% L6 h/ N3 u! U            df2.loc[len(df2)] = insert_datas
            df2.to_excel(writer1,sheet_name='销项信息',index=False)
        insert_datas = []
    flag = False
  A- y) n4 m) I* e3 v( _
- E6 i: E$ P* N3 `3 g" ^8 @" n

  j1 i  Q; b4 N0 v" O* O
并将是否违约插入到最后一列
) u, T3 |6 y! E

4 Q5 j) [- j  z, z1 M# 提取是否违约的列表
m = []
for name in code_list:
    m.append(datas[datas['code']==name]['break_contract'].tolist()[0])
    df1.loc[:,len(df1)] = m
& t; J7 V0 W. d3 o$ L! h    df1.to_excel(writer3,sheet_name='sheet1',index=False)
/ _: r% i6 W/ I' H建立模型
Logistics违约率预测模型
使用Logistics违约预测模型，代入所有的指标数据为自变量，是否违约为因变量，预测出违约率。
3 T4 k3 W8 f6 O% A
5 k" ~/ D; m4 m
/ I, k$ V5 {  V) K" H7 r1 jX=datas[['进项发票作废率','进项负数发票率','进项每月平均交易额','进项每月交易次数','销项发票作废率','销项负数发票率','销项每月平均交易额','销项每月交易次数','销售收入增长率']]
y=datas['是否违约']
; a  S  k! y5 O0 W. K+ pX_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()
& T. i0 y7 t1 A% a8 a  ~" k6 gmodel.fit(X_train,y_train)
+ X& l0 ]/ q+ w, C+ d: O3 b0 I; m6 sa=model.predict_proba(X_validation)
result=[]
+ Y, G: A6 u% ffor i in range(len(a)):
+ E9 ^# l% [3 Qif a[1]>0.5:
    result.append(1)
" Q. E# ?. e; i3 S9 ~( ?else:
9 j% |8 r/ N! T  N3 z1 _; Z: H    result.append(0)
( t" Q/ _$ l" dfrom sklearn import metrics
print('误差: %.4f' % (1-metrics.recall_score(y_validation,result,average='weighted')))
- R- J; A1 n. `; K最终得到一张我们的分析表格
' i1 z9 S3 Y" d


5 v" m' _! o/ S3 V
/ g0 G& j# c0 y6 J4 m% E. F

. n; H1 |) W2 [* d) R通过预测是否违约，我们就能解决贷不贷款的问题。


1 J( |# W. F. }7 Z) i8 N! w' E贷款金额
贷款金额的确认，根据该公司不违约率在所有公司中的权重，乘以总贷款金额确认：


; J. P9 g8 x+ G- R3 W. B- N  O) w6 }
+ V9 W# B5 x0 b  }* f2 l0 E
* I" _9 r& Q  G0 J
ri=1−Zi∑123j=1(1−Zi)×M
9 n; p+ ?/ |- l- Gri=1−Zi∑j=1123(1−Zi)×M
/ I9 _. g: ~( |' g( _1 b" ^因此，我们得到的贷款金额是违约率和贷款总金额组成的关系式，这在第二问中能起到重要作用。
$ z" O& [4 J  e; C

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

7 p# c7 y$ D! X$ R$ |! S5 ~
% b/ `" H9 I5 \/ d* q: N/ z0 ^3 h信誉等级        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
7 l* P3 j5 m1 t1 u, @5 UC        0.9982        y = 32.157x^3-207.386x^2+504.717*x-0.973
! R. {+ j9 L" a/ ^/ v+ y银行获利=贷款金额x贷款年利率x(1-利率对于信誉评级客户流失率)

) A( k% Y% J0 i: P, s
; g0 m0 o: @+ ~% r- i  b) f在贷款金额确认，贷款年利率范围在0.4~1.5的情况下，利用上面拟合的关系式，我们能够暴力跑出最优年利率。


4 Z/ {! G6 N/ C* |# O6 pdouble turnover_rate(double x, char ch) {
        double y = 0, result = 0;
1 C% P) n9 S* j3 J& a        switch (ch) {
% g* z5 ?+ W2 n        case 'A':
                y = 37.969520 * pow(x, 3) - 258.570452 * pow(x, 2) + 640.944427 * x - 1.121484;
8 l5 [) w$ |2 |# \) k! h$ a* k  ]/ \                result = x * (1 - y / 100.0);
                break;
$ k3 _; _2 c8 e  n        case 'B':
+ ^/ b! F- F6 u0 Z! q& P                y = 33.994698 * pow(x, 3) - 225.050538 * pow(x, 2) + 552.829151 * x - 1.016503;
                result = x * (1 - y / 100.0);
! J: [$ n3 S5 h% t                break;
        case 'C':
                y = 32.156864 * pow(x, 3) - 207.385880 * pow(x, 2) + 504.716993 * x - 0.973497;
6 [5 k. m! h1 P5 ^' ]- F2 k) R) g                result = x * (1 - y / 100.0);
                break;
        default:
                cout << "输出有误！"  << ch << endl;
+ W( f- N9 N( A  O/ Z+ n        }

        return result;
}
' d# g+ E+ E2 ~  `- z, E& Z

& l1 W! G3 ^' [2 p# J& T6 Z" v- D# X

问题二
利用代码，重新计算出各指标数据
代入Logistics违约率预测模型，预测出各公司的违约率
根据标准普尔评级建立，主标尺，对不同违约率进行A~D等级划分，信誉等级D不予贷款
, X3 q" `5 L/ V  o* D7 i" s! C将违约率代入，之前得到的公式，得到具体贷款金额
最优年利率沿用上一问
# 信用等级
cs = []
) I$ m8 A. q! s7 p! V# 最优年利率,客户流失率,利率值
7 T0 J% [6 U0 A4 Z, l8 H6 \tax = []
for i in m:
* o1 b" }$ K0 y' e) N    if i <= 0.0069264:
$ U7 ?- n% V: v' x        cs.append('A')
        tax.append([0.083,0.503173,0.0412366])
    elif i > 0.0069264 and i <= 0.22619:
        cs.append('B')
9 E, T. {1 s4 M; @        tax.append([0.097,0.505215,0.0479942])
8 {/ V7 S0 d4 }3 m2 B    elif i > 0.22619 and i <= 0.509915:
        cs.append('C')
& f  g/ A7 \6 \5 Q* P        tax.append([0.1069,0.506501,0.052755])
7 Q- D: g6 a' R1 C( H' U- O    elif i > 0.509915:
        cs.append('D')
        tax.append([0.15,0,0])
    else:
        print('违规')
        
: L  A4 W. r/ p; ~, C) kparr = []
4 ]8 }5 a- F/ rfor arr in list(a):
    parr.append(list(arr)[0])
sum_val = sum(parr)
' x3 ?2 b6 R) A8 Y2 _" Jamount = []
% b) r8 F! X& s! \for ival in parr:
    tmp = ival / sum_val * 100000000
    if ival < 1 - 0.509915:
9 P. w6 f+ g% U/ I9 V            amount.append(0)
    else:
            amount.append(tmp)

. K( p. i5 j/ f7 k9 I4 M
0 F. P9 m* S' t0 I" P$ z! Y
. S; X/ j* }+ v1 _! s. p
可以看到，贷款金额也都在10w~100w之内。

. s6 ~. w3 a/ U8 f
  _8 U0 D, z8 \问题三
/ U: N6 v4 o9 H4 Q5 X( Z- T. E这一问，我们做得有些匆忙了，有其他想法的可以按照自己的想法做做，这里只拿我们的做参照。


/ S. ]) Z  q4 z" V疫情对公司影响最大的就是每月平均销售额和每月平均销售数量，因此，
1 r5 L9 _# _( O# {0 I5 b9 n
  S" B5 [3 S# s
! K5 n2 q3 O. B# V; P$ {  u# d- v对每月平均销售额和每月平均销售数量，分别取随机数，数量取10w组，其他指标数据值不变，是否违约数据根据第二问结果，设为初始值
代入Logistics模型中，预测出每一组的违约率
6 P. C# d' H* R1 J& C判断每个公司违约率的变化情况，根据变化情况来增/减贷款金额和年利率。
% i- ]! M# T* }$ P7 A5 erand_num = pro_rand()
  L# A% Z# L9 _
$ Z4 L2 [  F, o/ Zdata_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':[]}
3 l4 l* H4 v( V# T" F! Fdf = DataFrame(data_form)

for val in rand_num:
        a = rate_func(val)
6 N: P' q& Q5 M5 S2 M        rate_list = get_rate(a)
: D6 Z/ Z' V- L% c9 A$ v; G        df.loc[len(df)]  = rate_list
! x2 N! h. O" C/ l3 ?        df.to_excel(writer,sheet_name='违约率变化',index=False)
销售数量和销售金额的随机数范围是：0~MAX
2 z  o* h: \) n4 C9 P7 I# v% ~


# C( X+ d4 a- L* m% e" }
————————————————
版权声明：本文为CSDN博主「Hk_Mayfly」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/qq_39542714/article/details/108592087
% A) d; A% _5 _6 a
0 X' V% P# G9 ^; D% q