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

杨利霞

2020 全国大学生数学建模竞赛C题思路+代码
3 q+ S% U2 [+ u, v题目链接：https://cloud.189.cn/t/ri2uUb7BRVJr
5 w& c) O; R! P+ e4 b' X! Q/ k

前言

' }9 z  g& s3 d. q0 j2 G


! y. v* n! A% Z' V又是一年数据挖掘题型，第一次接触这种题型还是在去年的mathorcup上，这种题的难度就在于指标的建立和数据的处理上。后面会出一份关于数据挖掘题型，我的相关经验，常用的工具和代码。
7 d9 w# H$ ~6 \, t2 x) h
" o* }, ?7 [( }+ Z4 G% h
) m8 E+ B4 t: |8 n3 q下面的一，二问实际都在解决

) C5 B4 f/ H$ a! Z/ @# H. |& y
' O2 g8 s5 ^. S: O3 s  q2 x贷不贷款？
贷款金额多少？
数据清洗
, h5 K0 t6 Z; E1 B  R5 ?  u7 k. O这道题的附件数据没有出现缺省或者异常数据，因此对于数据的预处理，更多的是根据问题的需求来做的。
& m' w+ F0 `! S& e# Y

将是否违约，违约设置为1，不违约设置为0
* r, }1 t* z  ~
7 Y! \( k, N1 y7 C6 o
+ Y' J/ W0 W0 M) G8 R* U% ^2 x信誉等级ABCD分别对应4,3,2,1
; O7 }, H/ C# q% A' a8 x5 G

9 t( k+ X" v( `! Y0 @! }: u发票状态，有效发票为a，作废发票为b

  p! k/ s" ^# M* |, a! z# M& S5 m) o
我将销项和进项所有数据，以公司代码为区别，提取到了不同的sheet当中，对于该公司有效发票数，作废发票数，负数发票数，方便对数据观察。
2 G1 S" n# H. n

# 遍历所有sheet数据
9 J; N- h# v, U6 ^! zfor xsn in sn.sheet_names[1:]:
    # 读取文件
    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)
0 l& b9 o# p1 |7 y9 W) y  q    # 找到全部公司名称代号
3 B* O" W3 U6 V: z& V) Y  b    code_list = list(set(list((datas['code']))))
: `7 u: G3 E% Q    for name in code_list:
        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()
0 {2 P6 o  H; V: _, [7 l        tmp_datas['a_count'] = list(count1)[0]
* B+ l, M$ _+ ]  ?! T! a, |- _8 E        if(len(count1) > 1):
            tmp_datas['b_count'] = list(count1)[1]
! `7 }* F/ F3 z1 S        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)
2 k6 r2 C3 k* P1 w) N        else:
6 v% W5 v% h' X& Z  |6 a            tmp_datas.to_excel(writer2,sheet_name=name,index=False)
' Q8 e3 t7 G6 D( U4 V1 k负数发票：在之前购买的物品，并开具了相关正向发票，后来退货所以开具了值为负数的发票，抵消前面正数发票的值。

+ ^4 N7 `: k; S" ]3 b8 \) K
& ?$ Q% v1 Y: x! C) p提取到信息：
* z: \. q! O  R$ E9 g
5 N( x& o8 V3 I) u6 j: t
3 i4 i6 J0 ?/ @部分公司数据记录很少，或者时间跨度大，需要综合数据指标，抵消数据数量和跨度大的影响
有些负数发票，在之前找不到对应的正数发票，可能是因为在数据记录日期之前购买的，在之后退款，因此在附件中找不到记录。
问题一
建立指标
1 N0 _/ a. u8 ]3 ~) ~) z3 r进项发票作废率，进项负数发票率，进项每月平均交易额，进项每月交易次数，
3 p( r* Y; |9 E3 _, W  t

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

, c9 G4 X- {( b
+ Z% o% O( ~) ^0 z/ N提取出相关指标到附件


for xsn in sn.sheet_names[1:]:
    # 读取文件
    datas = pd.read_excel(file_pos, sheet_name=xsn)
. {, J# s! U5 y' u    code_list = list(set(list((datas['code']))))
4 D$ |2 ]; @+ [9 z, M9 T6 W7 q    for name in code_list:
0 E2 c" q. {8 o        tmp_datas = datas[datas['code'] == name]
        tmp_datas.index = range(len(tmp_datas))
# ^5 V( W) C$ u. _  ?        insert_datas.append(name)
        # 作废数
, x& _# a1 c( W7 R1 O        cacel_count = len(tmp_datas[tmp_datas['tax_status'] == 'b'])
        # 有效数
, n' Y" O- y. `, h5 I2 K        valid_count = len(tmp_datas[tmp_datas['tax_status'] == 'a'])
2 X# D, m! w& Y% F7 _/ s        # 发票作废率
        count1 = (cacel_count / (cacel_count + valid_count))*100
) `: X2 y5 g1 N6 i; N, n        # 负数发票数
1 @+ R' _6 j: e/ [" h3 N- Y6 r        neg_count = len(tmp_datas[tmp_datas['cost'] < 0])
7 J2 V3 D0 K+ E% C/ q# k        # 负数发票率
        count2 = (neg_count / valid_count) * 100
        # 转换时间
9 D; h( t9 D0 n$ D        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
        # 时间最大值
        max_time = tmp_datas.iloc[0:,1].max()
        # 时间最小值
/ e3 c3 Q$ T0 s3 V6 F        min_time = tmp_datas.iloc[0:,1].min()
        # 时间差
        diff_time = months(max_time, min_time) + 1
        # 有效票
, {5 i+ A2 j' ~& B        valid_tax = tmp_datas[tmp_datas['tax_status'] == 'a']
- _' z1 n" d% q$ F        # 平均月交易额
        avg_money = valid_tax['totle_cost'].sum() / diff_time
        # 平均每月交易次数
. a: o$ |+ B. {! N1 \7 `& x, J  J        trans_count = len(tmp_datas) / diff_time
8 q# T% b4 A% f/ A5 `        insert_datas += [count1, count2, avg_money, trans_count,]
        if flag:
* U$ S6 C- B4 s            df1.loc[len(df1)] = insert_datas
/ j% d+ h5 }% H9 u1 U" |& ^4 g            df1.to_excel(writer1,sheet_name='进项信息',index=False)
        else:
: C1 z  P! b* K8 D/ p3 u6 V            merge_time = tmp_datas.groupby(tmp_datas['date']).sum()
            # 销售收入增长率
5 d: P  a8 \4 N/ P! H            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
) s" A; G2 ~1 H4 o; Q8 R            insert_datas.append(income_tax)
3 M8 ^1 Y0 D3 l! c" q- [: ]            df2.loc[len(df2)] = insert_datas
            df2.to_excel(writer1,sheet_name='销项信息',index=False)
) ?5 r9 w$ O" x7 c* W& b; C        insert_datas = []
5 x1 u+ X+ n. V# i+ {    flag = False
+ q8 v6 a" ^8 h

+ ]6 Q0 ^/ ?+ x3 ~$ W# W
) R2 D" F7 n. o* C% Y
并将是否违约插入到最后一列

' ^- O7 Z6 ~; C* W
# 提取是否违约的列表
" B; i7 ]! ^2 `9 O6 um = []
) ^3 A6 T) ~0 R+ V/ u$ jfor name in code_list:
( Q* ?7 P7 }! P! A8 q9 V    m.append(datas[datas['code']==name]['break_contract'].tolist()[0])
    df1.loc[:,len(df1)] = m
    df1.to_excel(writer3,sheet_name='sheet1',index=False)
! h: b' s  s% j$ G/ D& M6 m$ |5 @建立模型
Logistics违约率预测模型
- t- c  D" e3 K6 M/ p/ W使用Logistics违约预测模型，代入所有的指标数据为自变量，是否违约为因变量，预测出违约率。


X=datas[['进项发票作废率','进项负数发票率','进项每月平均交易额','进项每月交易次数','销项发票作废率','销项负数发票率','销项每月平均交易额','销项每月交易次数','销售收入增长率']]
y=datas['是否违约']
- ]; b1 Q! {6 C1 Y% NX_train, X_test, y_train, y_test = train_test_split(X,y, test_size=0.2,random_state=2020)
( J, @6 C0 c7 I# p! x/ K# HX_validation, X_test, y_validation, y_test = train_test_split(X_test,y_test, test_size=0.1,random_state=2020)
6 n6 Q  \, r; [9 E$ D8 pmodel = LogisticRegression()
1 T: |& w! v/ d* L$ T8 [, ymodel.fit(X_train,y_train)
a=model.predict_proba(X_validation)
8 ^7 w; `; `& J' _% y8 K; u% iresult=[]
% j- Y) m/ D7 J" K, ?for i in range(len(a)):
if a[1]>0.5:
    result.append(1)
# b& ?* |# M/ helse:
    result.append(0)
from sklearn import metrics
print('误差: %.4f' % (1-metrics.recall_score(y_validation,result,average='weighted')))
1 ^5 K8 y0 D1 P8 @% u- y最终得到一张我们的分析表格
# |7 a" X9 Z* N2 k# d6 ~4 Z

- P; h( E) }6 r. t
8 x" g% {- d% c8 r) G' A


# v0 o1 T5 q5 U$ D通过预测是否违约，我们就能解决贷不贷款的问题。


贷款金额
贷款金额的确认，根据该公司不违约率在所有公司中的权重，乘以总贷款金额确认：

' H  S8 H- F& w! o# B! W7 {, Z: Z
8 A0 r0 E+ r. J( _5 s( k
5 g* I! p4 }3 f( w5 _7 p2 M

ri=1−Zi∑123j=1(1−Zi)×M
) C. W8 \. F9 a' u- ?ri=1−Zi∑j=1123(1−Zi)×M
- U/ [' U- w; M6 o9 O因此，我们得到的贷款金额是违约率和贷款总金额组成的关系式，这在第二问中能起到重要作用。

8 Y2 M6 m. q" n) T- @
) E+ X/ d: n1 I! H) {贷款年利率
绘制出年利率与客户流失率图，可以分析出两者应该是有关系的。利用SPSS拟合出不同信誉等级，年利率与客户流失率的关系式。


% E5 @: p+ k! Q4 L信誉等级        R平方        关系式
( k' e/ e6 W% k4 B4 c5 {6 Y1 WA        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
  ]/ R/ V4 [+ X" _, ~8 `# L  j( hC        0.9982        y = 32.157x^3-207.386x^2+504.717*x-0.973
' u! L9 _3 X! D/ K0 v" v银行获利=贷款金额x贷款年利率x(1-利率对于信誉评级客户流失率)
) `! W2 G& \' p$ `( l: I
" K* H. W) W. z. L8 Y
6 G( N- |0 M2 a8 H0 O# T$ V在贷款金额确认，贷款年利率范围在0.4~1.5的情况下，利用上面拟合的关系式，我们能够暴力跑出最优年利率。
# w5 {+ r: b9 \, ]- x

double turnover_rate(double x, char ch) {
        double y = 0, result = 0;
* J1 y: [* Y. \, ^) d( }% ^        switch (ch) {
        case 'A':
                y = 37.969520 * pow(x, 3) - 258.570452 * pow(x, 2) + 640.944427 * x - 1.121484;
3 p$ w2 a9 ^& G0 [                result = x * (1 - y / 100.0);
                break;
        case 'B':
                y = 33.994698 * pow(x, 3) - 225.050538 * pow(x, 2) + 552.829151 * x - 1.016503;
                result = x * (1 - y / 100.0);
4 u& w' A$ y# M( W7 b+ _# ?& w/ [  g                break;
        case 'C':
                y = 32.156864 * pow(x, 3) - 207.385880 * pow(x, 2) + 504.716993 * x - 0.973497;
4 T0 M0 J* w$ ~7 b# o, V                result = x * (1 - y / 100.0);
/ I4 k; [+ o+ I0 d; T4 ?" j                break;
        default:
; {3 z3 l) L. e0 ]5 V$ Y2 F. {                cout << "输出有误！"  << ch << endl;
& k9 z) o3 \% ^  e6 U; j4 r; U' E* l        }
& _. c" v& ^& F' @( o
  j$ V# L/ t- q! @        return result;
}


/ q4 q% {2 ^# @+ a; P) _
- N; [# ]# ?6 c) Q4 d
& Q7 P+ c. q! x9 [/ P4 I1 ^3 l问题二
+ O5 `  T1 o: P* }利用代码，重新计算出各指标数据
. |0 M7 T/ P) ^3 f) h2 \# k, X代入Logistics违约率预测模型，预测出各公司的违约率
9 s; ^' x9 h, o$ G' B根据标准普尔评级建立，主标尺，对不同违约率进行A~D等级划分，信誉等级D不予贷款
将违约率代入，之前得到的公式，得到具体贷款金额
最优年利率沿用上一问
( d' f% b' R0 k8 h; u( x# 信用等级
0 Q0 a0 G0 R! l* n4 O& |$ Wcs = []
# 最优年利率,客户流失率,利率值
* R8 ]! A- N# rtax = []
& A; N, [/ O7 k# vfor i in m:
    if i <= 0.0069264:
, D- |6 q- w, B7 K- |! ?' c        cs.append('A')
        tax.append([0.083,0.503173,0.0412366])
    elif i > 0.0069264 and i <= 0.22619:
" O( p& D" a  P        cs.append('B')
        tax.append([0.097,0.505215,0.0479942])
    elif i > 0.22619 and i <= 0.509915:
        cs.append('C')
        tax.append([0.1069,0.506501,0.052755])
! @- j% K7 y7 Y. r( ~1 W. H! [    elif i > 0.509915:
        cs.append('D')
        tax.append([0.15,0,0])
    else:
9 @( C# w6 d' ~4 D2 k+ V; v4 q        print('违规')
        
4 @2 }; O* p2 y6 \+ T: P7 cparr = []
for arr in list(a):
    parr.append(list(arr)[0])
sum_val = sum(parr)
amount = []
for ival in parr:
2 E& F- D" ]4 M% I9 Z    tmp = ival / sum_val * 100000000
/ \* h* @4 R/ e" E  c& ~+ n    if ival < 1 - 0.509915:
            amount.append(0)
    else:
            amount.append(tmp)
! d9 e+ R- B' T: {7 n  M0 @
: b( H# n* E! y
1 h7 l8 k& Z& o8 s: `' Q( J% v6 F  q, H
: P  V! j" [) _& p0 I$ f' y2 {: b+ Y
7 a5 V$ }( A& ?4 m3 v2 l可以看到，贷款金额也都在10w~100w之内。
" u$ o# i! I' }8 v
' }# |! s6 [" M+ N3 L
问题三
这一问，我们做得有些匆忙了，有其他想法的可以按照自己的想法做做，这里只拿我们的做参照。


! \( I( v6 O- C6 J0 ]9 m2 F疫情对公司影响最大的就是每月平均销售额和每月平均销售数量，因此，

3 @; K' N2 Q8 v* ^
对每月平均销售额和每月平均销售数量，分别取随机数，数量取10w组，其他指标数据值不变，是否违约数据根据第二问结果，设为初始值
代入Logistics模型中，预测出每一组的违约率
; |5 e8 N( l. k+ f! M. Y判断每个公司违约率的变化情况，根据变化情况来增/减贷款金额和年利率。
) ?4 g' F! q+ a$ y" Qrand_num = pro_rand()
, X. V! e2 @$ O$ }. O7 @1 }
, j1 o" {. W2 pdata_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':[]}
7 `1 f4 i+ C- a2 x6 s4 Odf = DataFrame(data_form)
+ Z  a2 b: H. Y& c, y$ b
for val in rand_num:
# y& `  y! F+ |        a = rate_func(val)
        rate_list = get_rate(a)
        df.loc[len(df)]  = rate_list
' b. S, ?% }! V        df.to_excel(writer,sheet_name='违约率变化',index=False)
. L6 B/ s1 w, w. n4 M0 L销售数量和销售金额的随机数范围是：0~MAX
9 ~, h& a. l9 J, b( G( B6 k7 ?* ]3 W
2 F9 R, a9 Y- s* d( b4 ~


1 P2 F! [" H% n————————————————
版权声明：本文为CSDN博主「Hk_Mayfly」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
( P( N. K( |# v) R原文链接：https://blog.csdn.net/qq_39542714/article/details/108592087


0 p- a- [4 S  c+ C& f" o