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

杨利霞

2020 全国大学生数学建模竞赛C题思路+代码
/ K% S* T# x, b5 p0 x8 @题目链接：https://cloud.189.cn/t/ri2uUb7BRVJr
3 }% A4 Z: B3 s- x/ g$ @' x

前言
. d" l( K. e, {9 q4 ]) ^$ D4 c( o

2 E. a: S# K1 d, u8 T* C
8 ~' ?) k1 _2 k3 D' Q
9 k3 u: g( o6 d7 V6 g0 O- [& R: t7 p又是一年数据挖掘题型，第一次接触这种题型还是在去年的mathorcup上，这种题的难度就在于指标的建立和数据的处理上。后面会出一份关于数据挖掘题型，我的相关经验，常用的工具和代码。
# G2 o7 V2 w2 }6 d; V# i& |

下面的一，二问实际都在解决

/ E9 v4 g7 Q8 m+ o
; b4 Z" r( L) G# J7 @贷不贷款？
贷款金额多少？
6 z7 s( ^2 @  d: v  U) O/ l. X数据清洗
这道题的附件数据没有出现缺省或者异常数据，因此对于数据的预处理，更多的是根据问题的需求来做的。
( o$ ^0 @( W3 A% o6 x# y; c

0 [8 I1 D5 _7 C, d' H( f' _8 h将是否违约，违约设置为1，不违约设置为0

1 v# j2 m3 B% X0 G! _
信誉等级ABCD分别对应4,3,2,1


9 `/ Y4 y' Z' Z+ b4 i; e发票状态，有效发票为a，作废发票为b

) R7 E( n6 z" M" H% N& s" U
我将销项和进项所有数据，以公司代码为区别，提取到了不同的sheet当中，对于该公司有效发票数，作废发票数，负数发票数，方便对数据观察。

- ~" s1 a+ f; j3 f; p
# 遍历所有sheet数据
for xsn in sn.sheet_names[1:]:
4 u" U/ p' v5 v, |2 [    # 读取文件
    datas = pd.read_excel(file_pos, sheet_name=xsn)
* g9 J, B, g' u/ ~/ g    datas['date']=pd.to_datetime(datas['date'],format='%Y/%m/%d')
    datas.set_index('date', drop=True)
    # 找到全部公司名称代号
    code_list = list(set(list((datas['code']))))
7 U7 c, F+ ?. O& _    for name in code_list:
        tmp_datas = datas[datas['code'] == name]
        tmp_datas.index = range(len(tmp_datas))
( u) c6 u$ |1 N0 l$ ^2 s1 W* y6 Z        # 转换日期未object类型
        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
. V% Y" O6 Q2 V! O! b        count1 = tmp_datas['tax_status'].value_counts()
        tmp_datas['a_count'] = list(count1)[0]
        if(len(count1) > 1):
) q) B' O' ?; L2 [            tmp_datas['b_count'] = list(count1)[1]
0 l$ f& o7 P% T$ [! t/ m        tmp2 = tmp_datas[tmp_datas['cost'] < 0]
/ }% p0 ?* ^7 `# [& N4 p/ d* E* e        tmp_datas['neg_value_tax'] = len(tmp2)
7 u- u9 {8 b1 X( H4 @0 s* D        if xsn == sn.sheet_names[1]:
; q- m2 X& D/ }            tmp_datas.to_excel(writer1,sheet_name=name,index=False)
. r- M3 i) S% {% g) I, u% V" o        else:
* u* `2 ~  z4 W) q# e            tmp_datas.to_excel(writer2,sheet_name=name,index=False)
0 w- E2 G- }+ l负数发票：在之前购买的物品，并开具了相关正向发票，后来退货所以开具了值为负数的发票，抵消前面正数发票的值。
, }( i5 |$ Y7 y

提取到信息：

: y' ~+ S* P. b. X9 E4 w5 _
部分公司数据记录很少，或者时间跨度大，需要综合数据指标，抵消数据数量和跨度大的影响
有些负数发票，在之前找不到对应的正数发票，可能是因为在数据记录日期之前购买的，在之后退款，因此在附件中找不到记录。
问题一
建立指标
* b/ f* s& K- t* e进项发票作废率，进项负数发票率，进项每月平均交易额，进项每月交易次数，
8 u5 c) \8 R' ?8 B5 [' G5 h7 f' Z
4 `$ v- O; Q1 S' S9 l: j& Q
销项发票作废率，销项负数发票率，销项每月平均交易额，销项每月交易次数，销售收入增长率

/ J& H1 P+ J  e, D
提取出相关指标到附件
6 k9 D2 [5 j4 W; J  f

for xsn in sn.sheet_names[1:]:
    # 读取文件
; K/ F( X* I$ ^    datas = pd.read_excel(file_pos, sheet_name=xsn)
    code_list = list(set(list((datas['code']))))
    for name in code_list:
: e3 h! l3 |/ R3 W4 H        tmp_datas = datas[datas['code'] == name]
        tmp_datas.index = range(len(tmp_datas))
        insert_datas.append(name)
+ P) E6 i( Q3 u4 T, V1 C3 @1 T/ u        # 作废数
        cacel_count = len(tmp_datas[tmp_datas['tax_status'] == 'b'])
* p% R5 p4 v1 |  Y  }: S2 J4 h/ }) S        # 有效数
        valid_count = len(tmp_datas[tmp_datas['tax_status'] == 'a'])
4 @/ M2 i1 l; z9 c        # 发票作废率
) K1 I1 `4 M) n" j$ d5 W        count1 = (cacel_count / (cacel_count + valid_count))*100
: h% X- F/ g/ Z. W- A+ e! O        # 负数发票数
* g" x- [* w7 H) T* r7 n9 [        neg_count = len(tmp_datas[tmp_datas['cost'] < 0])
        # 负数发票率
        count2 = (neg_count / valid_count) * 100
        # 转换时间
5 |9 J- O  u# L- c        tmp_datas['date'] = [x.strftime('%Y/%m/%d') for x in tmp_datas['date']]
! X, q% g1 p( g4 D        # 时间最大值
        max_time = tmp_datas.iloc[0:,1].max()
/ ?& @+ a+ C3 M4 I5 n; u! z        # 时间最小值
        min_time = tmp_datas.iloc[0:,1].min()
4 M( P1 D* c. j, b& s! `        # 时间差
! l7 d. S% K4 D/ y! N/ ]        diff_time = months(max_time, min_time) + 1
        # 有效票
- M" n! ?1 ]. B3 C: y        valid_tax = tmp_datas[tmp_datas['tax_status'] == 'a']
        # 平均月交易额
8 s; q5 i# g8 K! S9 w! u        avg_money = valid_tax['totle_cost'].sum() / diff_time
        # 平均每月交易次数
1 g* V, ^+ [2 `$ g& o5 D* u9 D        trans_count = len(tmp_datas) / diff_time
        insert_datas += [count1, count2, avg_money, trans_count,]
        if flag:
            df1.loc[len(df1)] = insert_datas
            df1.to_excel(writer1,sheet_name='进项信息',index=False)
7 n4 b! B( i2 z  f* E5 r        else:
            merge_time = tmp_datas.groupby(tmp_datas['date']).sum()
            # 销售收入增长率
6 F# t5 X* r; c            income_info= list((merge_time['cost'] - merge_time['cost'].shift(1)).fillna(1))
3 q' ^. {* X7 J9 ]' q            diff_time_day = days(max_time,min_time)
            income_tax = (sum(income_info) / diff_time_day)*100
            insert_datas.append(income_tax)
) K- b* x" _+ O' H& Y            df2.loc[len(df2)] = insert_datas
' k; e. {7 U  u- X; m, w            df2.to_excel(writer1,sheet_name='销项信息',index=False)
        insert_datas = []
    flag = False


  R* K; |$ I! Z6 E

$ r9 \* d4 {# a" M- c并将是否违约插入到最后一列


# 提取是否违约的列表
5 }9 }+ w) B9 `, Rm = []
for name in code_list:
3 D2 _) e- [3 a+ I    m.append(datas[datas['code']==name]['break_contract'].tolist()[0])
    df1.loc[:,len(df1)] = m
0 M- k+ G# \4 u& T8 h    df1.to_excel(writer3,sheet_name='sheet1',index=False)
建立模型
Logistics违约率预测模型
使用Logistics违约预测模型，代入所有的指标数据为自变量，是否违约为因变量，预测出违约率。
- a/ T4 z0 R; X4 ~( E

% U- i' P( r5 Y& NX=datas[['进项发票作废率','进项负数发票率','进项每月平均交易额','进项每月交易次数','销项发票作废率','销项负数发票率','销项每月平均交易额','销项每月交易次数','销售收入增长率']]
8 |' t/ y1 {, W5 p9 M$ n9 L  Ny=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)
model = LogisticRegression()
model.fit(X_train,y_train)
) v: Z! V8 t2 ]0 b  ~2 s$ B9 w6 U2 Aa=model.predict_proba(X_validation)
4 V8 n+ p5 S" d! h: P" nresult=[]
3 t! n8 t2 K. j2 A/ d# T! Cfor i in range(len(a)):
- Y' X2 a4 A. o* N# L- `8 q5 k7 i6 xif a[1]>0.5:
3 b5 g5 M: t8 k2 D6 ~& L    result.append(1)
else:
    result.append(0)
( ?( R3 _" a2 k$ T. Lfrom sklearn import metrics
print('误差: %.4f' % (1-metrics.recall_score(y_validation,result,average='weighted')))
( P; \* O+ w2 c3 O$ E最终得到一张我们的分析表格


; J  h! U5 B6 ]+ h7 x
0 [% X. ?, b" v
0 h2 ?9 W" h" Y' m- x
  d# \& r- @* M$ a: q
通过预测是否违约，我们就能解决贷不贷款的问题。

4 F' l; Q4 B9 Q0 P5 w$ m
$ _( z; c  q9 A( J* {* P贷款金额
贷款金额的确认，根据该公司不违约率在所有公司中的权重，乘以总贷款金额确认：
- i- n2 E" _- P0 Z4 L' ]& D
+ l# h, X; W) ]# r* X: h& D/ _% b
  A% r7 A$ Z0 z' H% @


' n' O8 [9 o: m7 |% ~ri=1−Zi∑123j=1(1−Zi)×M
ri=1−Zi∑j=1123(1−Zi)×M
因此，我们得到的贷款金额是违约率和贷款总金额组成的关系式，这在第二问中能起到重要作用。

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

* }$ K1 L6 W" k9 F+ j& j5 T
信誉等级        R平方        关系式
A        0.9977        y = 37.97x^3-258.57x^2+640.944*x -1.121
1 _! w* U' }1 Z4 T. eB        0.9982        y = 33.995x^3-225.051x^2+552.829*x-1.017
C        0.9982        y = 32.157x^3-207.386x^2+504.717*x-0.973
银行获利=贷款金额x贷款年利率x(1-利率对于信誉评级客户流失率)


在贷款金额确认，贷款年利率范围在0.4~1.5的情况下，利用上面拟合的关系式，我们能够暴力跑出最优年利率。


double turnover_rate(double x, char ch) {
0 i- ^7 S7 s% B        double y = 0, result = 0;
        switch (ch) {
# @6 @2 f) ~& R3 n        case 'A':
/ g8 h( s# W% {6 n- ?0 g" E                y = 37.969520 * pow(x, 3) - 258.570452 * pow(x, 2) + 640.944427 * x - 1.121484;
                result = x * (1 - y / 100.0);
* d# J6 g- T+ ^                break;
+ J3 n8 {  E4 f# t; ^7 n        case 'B':
                y = 33.994698 * pow(x, 3) - 225.050538 * pow(x, 2) + 552.829151 * x - 1.016503;
* u9 e; f4 Q' i* o2 i8 g& O                result = x * (1 - y / 100.0);
, _( B1 I9 f5 F; [+ ~1 k+ }  h                break;
; \! E1 B* n" v& P/ T        case 'C':
                y = 32.156864 * pow(x, 3) - 207.385880 * pow(x, 2) + 504.716993 * x - 0.973497;
" t  `* O: U; i6 t6 d# o                result = x * (1 - y / 100.0);
                break;
, @7 S, z5 @$ d3 C  t8 a7 }        default:
# L  r2 ^7 e2 i& g                cout << "输出有误！"  << ch << endl;
/ i) V0 D7 f! I+ F        }
4 Z6 I& F: ~) A$ W; U
        return result;
}
' X) h  d+ n& U% N( J( W4 U
; l6 ~* O: ?( H2 r6 |; G. Q
' f5 O* y* O" ?% ^! V

" q, v6 {1 V  d7 {0 @问题二
利用代码，重新计算出各指标数据
8 D0 I. S# b2 p4 T1 I5 q6 _代入Logistics违约率预测模型，预测出各公司的违约率
1 E; \4 w9 }1 S, @! @/ K( n3 F根据标准普尔评级建立，主标尺，对不同违约率进行A~D等级划分，信誉等级D不予贷款
3 u$ O+ |* w6 L- c# G将违约率代入，之前得到的公式，得到具体贷款金额
7 N: T! t' ?3 A9 E. L最优年利率沿用上一问
# 信用等级
cs = []
# 最优年利率,客户流失率,利率值
tax = []
! k1 c# y/ C" }/ Lfor i in m:
    if i <= 0.0069264:
        cs.append('A')
- {# I% P/ B3 |0 T; `' G        tax.append([0.083,0.503173,0.0412366])
/ \, L; V& q3 V. s9 o    elif i > 0.0069264 and i <= 0.22619:
        cs.append('B')
        tax.append([0.097,0.505215,0.0479942])
. K' f3 M! z. {  o3 D    elif i > 0.22619 and i <= 0.509915:
        cs.append('C')
# j7 P  S% W0 e1 H" z" O6 n( {8 M        tax.append([0.1069,0.506501,0.052755])
    elif i > 0.509915:
        cs.append('D')
5 N* V, I5 o6 j( @9 m! u4 w. H        tax.append([0.15,0,0])
& u) \3 U) s$ `5 g' [% N    else:
2 A0 M" a+ Z9 j9 Q  {( v' _! j        print('违规')
+ a5 `8 |' v9 A9 m* _/ _        
parr = []
) I( R, B* N" p& Pfor arr in list(a):
    parr.append(list(arr)[0])
sum_val = sum(parr)
5 ?; C' R& i) Y9 B% ?0 yamount = []
for ival in parr:
& D1 H$ z2 |5 H' k    tmp = ival / sum_val * 100000000
2 _  O: c6 {1 P* b$ l8 f2 C) v! Q    if ival < 1 - 0.509915:
9 Z3 w1 ^' t: Y! u            amount.append(0)
% B( @' `1 ]( |( D2 r  O    else:
            amount.append(tmp)

$ X+ e9 e/ ]( ~0 m- r8 |  x

% _8 h6 p, g6 ^( R% U' v$ N; @
) t! s: O) X  e) @1 U可以看到，贷款金额也都在10w~100w之内。


问题三
, u* G% I& D! N% i$ T这一问，我们做得有些匆忙了，有其他想法的可以按照自己的想法做做，这里只拿我们的做参照。


疫情对公司影响最大的就是每月平均销售额和每月平均销售数量，因此，


$ q3 p* v, w& t& k2 O; r对每月平均销售额和每月平均销售数量，分别取随机数，数量取10w组，其他指标数据值不变，是否违约数据根据第二问结果，设为初始值
8 h8 y$ R: `! N' \2 z1 C" e* N' y1 E4 k" A代入Logistics模型中，预测出每一组的违约率
' E( e( t: U1 J% v7 s判断每个公司违约率的变化情况，根据变化情况来增/减贷款金额和年利率。
& O, M- h7 z) n9 ~% `4 Nrand_num = pro_rand()
4 K5 H& H5 U2 U- H9 d9 A
& ~: y- w" L% |7 kdata_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)

8 E/ i5 Z$ L2 E( D# b0 k) F1 ffor val in rand_num:
; Q4 _$ h2 n1 ^1 n; [        a = rate_func(val)
        rate_list = get_rate(a)
        df.loc[len(df)]  = rate_list
        df.to_excel(writer,sheet_name='违约率变化',index=False)
销售数量和销售金额的随机数范围是：0~MAX
* e* p  ~# G1 ?# M' O


9 W3 H& F1 @! L- d7 A
————————————————
! P! D. F% u$ {/ ]版权声明：本文为CSDN博主「Hk_Mayfly」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
* C4 h1 i0 W7 I* K  r% Y- i原文链接：https://blog.csdn.net/qq_39542714/article/details/108592087



2 B" E2 T! M% N! k) c4 c* |; Q3 h