Python商业数据分析代码

1.商业数据分析及其工具

importpandasaspd

#读取本地数据

df=pd.readcsv(，/2glkx/data/al2-l.csv*)

#读取网上数据

importpandasaspd

data_url=

“https:〃raw.githubusercontent.com/alstat/Analysis-with-Programming/master/2014

/Python/Numerical-Descriptions-of-the-Data/data.csv〃

df=pd.readcsv(dataurl)

print(df.head())

print(df.tail())

print(df.columns)

Index([u,Abra，,u'Apayao',u'Benguet',Ifugao,,u'Kalinga'],dtype='object')

#Extractingrownamesortheindex

print(df.index)

Rangeindex(start=0,stop=79,step=l)

#Transposedata

print(df.T)

print(df.ix[:,0].headO)现在改版ix用loc

print(df.ix[10:20,0:3])

print(df.drop(df.columns[[2,3]],axis=1).headO)

print(df.describe())

fromscipyimportstatsasss

#Performonesamplet-testusing1500asthetruemean

print(ss.ttestlsamp(a=df.ix[:»'Abra'],popmean=15000))

Ttest_lsampResult(statistic=-l.1281738488299586,pvalue=0.26270472069109496)

print(ss.ttest_lsamp(a二df,popmean=15000))

importmatplotlib.pyplotaspit

pit.show(df.plot(kind='box'))

1

importseabornassns

#Dotheboxplot

pit.show(sns.boxplot(df))

defadd2int(x,y):

returnx+y

print(add2int(2,2))

2Python商业数据的存取

importpandasaspd

importnumpyasnp

a=[，apple，,'pear','watch','money']

b=[[l,2,3,4,5],[5,7,8,9,0],[1,3,5,7,9],[2,4,6,8,0]]

d二diet(zip(a,b))

d

p=pd.DataFrame(d)

P

p.to_csv(，F:\\2glkx\\data\\IBM.csv')

pd.readcsv('F:\\2glkx\\data\\al2-l.csv')

importpandasaspd

importnumpyasnp

df=pd.readexcel(*F:\\2glkx\\data\\al2-2.xls')

df.head()

importpandasaspd

importnumpyasnp

#读取数据并创建数据表，名称为data。

data=pd.DataFrame(pd.readexcel('G:\\2glkx\\data\\al2-2.xls'))

#查看数据表前5行的内容

data,head()

importtushareasts

df=ts.get_hist_data(,000875))#从网上取数据

#直接保存

#df.to_csv(5:/2glkx/data/000875.csv')

#选择数据保存

df.to_csv(*F:/2glkx/data/000875.csv',columns=['open','high','low','close'])

importpandasaspd

importnumpyasnp

df=pd.read_csv(，F:/2glkx/data/000875.csv')

2

df.head()

importtushareasts

importos

filename=*F:/2glkx/data/bigfile.csv

forcodein「000875','600848','000981']:

df=ts.gethistdata(code)

ifos.path,exists(filename):

#df.tocsv(fi1ename,mode='a,header=None)

df.to_csv(filename,mode='a')

else:

df.to_csv(filename)

importpandasaspd

importnumpyasnp

df=pd.readexcelCF:/2glkx/data/000875.xls')

df.head()

importtushareasts#需先安装tushare程序包

#此程序包的安装命令：pipinstalltushare

importpandasaspd

importnumpyasnp象中

data=pd.DataFrame()

datal=ts.get_hist_data(，600000,)

#需要修改上面的时间

datal=datal['close']

datal=datal[:：-l]#按日期从远到近结束

data]'600000']=datal

data2=ts.gethistdata('000980')

data2=data2[，close，]

data2=data2[:：-l]

data]'000980']=data2

data3=ts.gethistdata(500098T)

data3=data3[，close']

data3=data3[:：-l]

datal7000981']=data3

data.info()#查看数据情况

data=data.dropna()

data.info()

data,head()

data,tail()

3

data=datatf600000,,‘000981']]

data,head()

data.ix[1:4]#现在ix改为loc

data.iloc[:2,:3]

importtushareasts

importpandasaspd

pd.setoptionexpandframerepr，,False)#显示所有歹ll

ts.set_token(,yourtoken,)#获取token号，需要先注册

pro=ts.pro_api()

stock_data=pro.daily(ts_code=，000001.SZ',start_date=,2010010T,

enddate』20190101，)

stock_data.head()

importpandas_datareader.dataasweb

importdatetime

start=datetime,datetime(2017,1,1)#获取数据的时间段一起始时间

end=datetime.date.today()#获取数据的时间段-结束时间

stock=web.DataReader(，z600797.SSZ，,“yahoo",start,end)

#获取浙大网新2017年1月1日至今的股票数

stock,head()#打印DataFrame数据前5行

importpandasaspd

frompandasimportSeries,DataFrame

importnumpyasnp

importmatplotlib.pyplotaspit

frompandas_datareaderimportdata,wb

fromdatetimeimportdatetime

end=datetime.now()

start=datetime(end.year-1,end.month,end.day)

df=data.DataReader(J600797.SS?,'yahoo',start,end)

dffAdjClose'].plot(legend=True,figsize=(10,4))

pit.show()

importnumpyasnp

importpandasaspd

importpandas_datareader.dataasweb

importdatetime

#获取600797.SS浙大网新数据

4

df_csvsave=

web.DataReader(，z600018.SS","yahoo”,datetime.datetime(2019,1,1),datetime,date,to

day())

print(df_csvsave)

dfcsvsave,tocsv(r'F:\2glkx\data\600018.csv',columns=df_csvsave.columns,index=

True)

importtushareasts

importpandasaspd

pd.set_option(,expand_frame_repr,,False)#显示所有歹U

ts.settokenyourtoken，)#获取token号

pro=ts.pro_api()

code_list=[5000001.SZ\'600000.SH','000002.SZ']

stock_data=pd.DataFrame()

forcodeincodelist:

print(code)

df=pro.daily(tscode=code,startdate=，2018010T,enddate='20180104')

stock_data=stock_data.append(df,ignore_index=True)

print(stockdata)

grouped=stock_data.groupby('ts_code')

print(grouped)

3Python商业数据的图形绘制与可视化

importmatplotlib.pyplotaspit

importpandasaspd

importnumpyasnp

df=pd.read_excel(，ZG:/2glkx/data/al3-1.xls〃)

#或者df=pd.read_excel(*G:\\2glkx\\data\\al3-l.xls，)

df.head()

fig=plt.figure()

ax=fig.addsubplot(1,1,1)

ax.hist(df['Age'],bins=7)

pit.show()

fig=plt.figure()

ax=fig.add_subplot(1,1,1)

ax.hist(df['Age'],bins=7)

pit.title(JAgedistribution')

5

pit.show()

fig=plt.figure()

ax=fig.add_subplot(1,1,1)

ax.hist(df['Age'],bins=7)

pit.titleCAgedistribution,)

pit.xlabel('Age')

pit.ylabelC#Employee,)

pit.show()

fig=pit.figure()

ax=fig.addsubplot(1,1,1)

ax.scatter(df[,Age'],df['Sales'])

#Youcanalsoaddmorevariablesheretorepresentcolorandsize,

pit.title(5Age&SalesScatterofEmployee,)

#Variable

pit.xlabelCAgef)

pit.ylabel('Sales*)

pit.show()

fig=pit.figure()

ax=fig.addsubplot(1,1,1)

ax.scatter(df[*Age*],df['Sales'],s=df['Income'])

#Addedthirdvariableincomeassizeofthebubble

pit.xlabelCAge))

pit.ylabel('Sales))

pit.show()

importmatplotlib.pyplotaspit

importpandasaspd

fig=pit.figure()

ax=fig.addsubplot(1,1,1)

#Variable

ax.boxplot(dfAge'])

pit.title(，BoxfigureofAge*)

pit.show()

vars=[，Age','Sales']

data=df[vars]

pit.show(data,plot(kind='box'))

importmatplotlib.pyplotaspit

importpandasaspd

6

importnumpyasnp

df=pd.read_excel(^G:/2glkx/data/al3-l.xls")

df.head()

var=df.groupby(['Gender']).sum().stack()

temp=var.unstack()

xlist=temp[，Sales']

label_list=temp.index

pit.axis("equal")

pit.pie(x_list)

pit.title(^Pastafatianismexpenses")

pit.show()

frompylabimport*

#makeasquarefigureandaxes

figured,figsize=(6,6))

ax=axes([0.1,0.1,0.8,0.8])

fracs=[60,40]#每一块占得比例，总和为100

explode=(0,0.08)#离开整体的距离，看效果

labels='Male'Female'#对应每一块的标志

pie(fracs,explode=explode,labels=labels,autopct=，%1.If%%',shadow=True,

startangle=90,colors=(〃g〃，〃r"))

title(，RateofMaleandFemale，)#标题

show()

importmatplotlib.pyplotasplot

importpandasaspd

importnumpyasnp

df=pd.readexcel(，，G:/2glkx/data/al3-l.xls〃)

df.head()

var=df.groupby('Gender，).Sales.sum()

#groupedsumofsalesatGenderlevel

fig=pit.figure()

axl=fig.addsubplot(1,1,1)

axl.set_xlabel('Gender，)

axl.set_ylabel('SumofSales')

axl.set_title(，，GenderwiseSumofSales")

var.plot(kind='bar，)

var=df.groupby(['BMI','Gender，]).Sales.sum()

var.unstack().plot(kind=,bar',stacked二True,color=['red','blue'])

importmatplotlib.pyplotasplot

7

importpandasaspd

importnumpyasnp

df=pd.readexcel(，zG:/2glkx/data/al3-l.xls〃)

df.head()

var=df.groupbyCBMT).Sales.sum()

fig=plt.figure()

axl=fig.add_subplot(1,1,1)

axl.set_xlabel(*BMI')

axl.set_ylabel('SumofSales，)

axl.settitle(，zBMIwiseSumofSales")

var.plot(kind='line')

importpandasaspd

importnumpyasnp

data=pd.DataFrame(pd.read_excelCG:\\2glkx\\data\\al3-3.xls'))

data,head()

t=np.array(data[[，year,1])

x=np.array(datatotal*]])

y=np.array(data[new*]])

importpylabaspl

pl.plot(t,x)

pl.plot(t,y)

pl.show()

importpylabaspl

pl.plot(t,x)

pl.plot(t,y)

pl.title(，populationcensus')

pl.xlabel('Time,)

pl.ylabel('Population')

pl.show()

importpandasaspd

importnumpyasnp

data=pd.DataFrame(pd.read_excel('G:\\2glkx\\data\\al3-4.xls'))

data,head()

t=np.array(data[['year'1])

x=np.array(data[['number']])

importpylabaspl

pl.plot(t,x)

pl.title(J1998-2015ofAlistedcompaniesinchina')

pl.xlabel('Time')

pl.ylabel(*companiesnumbers，)

8

pl.show()

importpandasaspd

importnumpyasnp

importpylabaspl

data=pd.DataFrame(pd.read_excelCG:\\2glkx\\data\\al3-4.xls'))

data,head()

t=np.array(data[['year']])

x=np.array(data[[，number511)

importpylabaspl

pl.plot(t,x,'ro')

pl.title(，1998-2015ofAlistedcompaniesinchina，)

pl.xlabel('Time*)

pl.ylabel(*companiesnumbers')

pl.show()

importrandom

importnumpyasnp

importmatplotlibasmpl

importmatplotlib.pyplotaspit

importmatplotlib.datesasmdates

frommpl_toolkits.mplot3dimportAxes3D

mpl.rcParams[，font,size*]=10

fig=pit.figure()

ax=fig.addsubplot(111,projection='3d')

forzin[2011,2012,2013,2014]:

xs=ranged,13)

ys=1000*np.random,rand(12)

color=plt.cm.Set2(random,choice(range(pit.cm.Set2.N)))

ax.bar(xs,ys,zs=z,zdir='y',color=color,alpha=0.8)

ax.xaxis.setmajor_locator(mpl.ticker.FixedLocator(xs))

ax.yaxis.set_major_locator(mpl.ticker.FixedLocator(ys))

ax.setxlabel('Month')

ax.set_ylabel('Year')

ax.setzlabel(JSalesNet[usd]')

pit.show()

frommpl_toolkits.mplot3dimportAxes3D

frommatplotlibimportcm

importmatplotlib.pyplotaspit

importnumpyasnp

n_angles=36

n_radii=8

#Anarrayofradii

9

#Doesnotincluderadiusr=0,thisistoeliminateduplicatepoints

radii=np.1inspace(0.125,1.0,n_radii)

#Anarrayofangles

angles=np.1inspace(0,2*np.pi,n_angles,endpoint=False)

#Repeatallanglesforeachradius

angles=np.repeat(angles,np.newaxis],n_radii,axis=l)

#Convertpolar(radii,angles)coordstocartesian(x,y)coords

#(0,0)isaddedhere.Therearenoduplicatepointsinthe(x,y)plane

x=np.append(0,(radii*np.cos(angles)).flatten())

y=np.append(0,(radii*np.sin(angles)),flatten())

#Pringlesurface

z=np.sin(-x*y)

fig=pit.figure()

ax=fig.gca(projection^Bd))

ax.plottrisurf(x,y,z,cmap=cm.jet,linewidth=0.2)

pit.show()

4Python描述性统计

#两个常用的统计包

importscipy.statsasstats

importnumpyasnp

#我们拿两个数据集来举例

xl=[1,2,2,3,4,5,5,7]

x2=xl+[100]

print('xl的平均值:',sum(xl),len(xl),，=，,np.mean(xl))

print('x2的平均值：'，sum(x2),',len(x2),,=,,np.mean(x2))

xl的平均值：29/8=3.625

x2的平均值：129/9=14.333333333333334

print(,xl的中位数：'，np.median(xl))

print(*x2的中位数：'，np.median(x2))

print(JOnemodeofxl:',stats,mode(xl)[0][0])

#因此我们自定义一个求众数的函数

defmode⑴：

#统计列表中每个元素出现的次数

counts={}

forein1:

ifeincounts:

counts[e]+=1

else:

10

counts[e]=1

#返回出现次数最多的元素

maxcount=0

modes={}

for(key,value)incounts,items():

ifvalue>maxcount:

maxcount=value

modes={key}

elifvalue==maxcount:

modes.add(key)

ifmaxcount>1orlen(l)==1:

returnlist(modes)

return'Nomode,

printCAllofthemodesofxl/,mode(xl))

importscipy.statsasstats

importnumpyasnp

#获取收益率数据并计算出mode

start='2014-01-01'

end=J2015-01-0T

pricing=D.history_data(，000002.SZA',fields=[,close)],start_date=start,

enddate=end)['close，]

returns=pricing.pct_change()[1:]

print('收益率众数：'，stats,mode(returns))

#由于所有的收益率都是不同的，所以我们使用频率分布来变相计算mode

hist,bins=np.histogram(returns,20)#将数据分成20个bin

maxfreq=max(hist)

#找出哪个bin里面出现的数据点次数最大，这个bin就当做计算出来的mode

print(，Modeofbins:',[(bins[i],bins[i+l])fori,jinenumerate(hist)ifj==

maxfreq])

#使用Scipy包中的gmean函数来计算几何平均值

printCxl几何平均值:',stats,gmean(xl))

printCx2几何平均值:',stats,gmean(x2))

#在每个元素上增加1来计算几何平均值

importscipy.statsasstats

importnumpyasnp

ratios=returns+np.ones(len(returns))

11

RG=stats,gmean(ratios)-1

print('收益率的几何平均值:'，R_G)

T=len(returns)

init_price=pricing[0]

final_price=pricing[T]

print('最初价格:',init_price)

print()最终价格:'，final_price)

print('通过几何平均收益率计算的最终价格:',initprice*(1+RG)**T)

print('xl的调和平均值:'，stats,hmean(xl))

print(,x2的调和平均值:，stats,hmean(x2))

importnumpyasnp

np.random.seed(121)

#生成20个小于100的随机整数

X=np.random.randint(100,size=20)

#Sortthem

X=np.sort(X)

print(JX:%s'%(X))

mu=np.mean(X)

printCX的平均值:',mu)

printCRangeofX:%s，%(np.ptp(X)))

abs_dispersion=[np.abs(mu-x)forxinX]

MAI)=np.sum(abs_dispersion)/len(abs_dispersion)

print('X的平均绝对偏差:‘，MAD)

print('X的方差:',np.var(X))

printCX的标准差:，np.std(X))

k=1.25#随便举的一个k值

dist=k*np.std(X)

1=[xforxinXifabs(x-mu)<=dist]

print('k值'，k,'在k倍标准差距离内的样本为:‘，1)

print('验证'，float(len(1))/len(X),>>J,1-l/k**2)

#没有现成的计算下偏方差的函数，因此我们手动计算：

lows=[eforeinXife<=mu]

semivar=np.sum((lows-mu)**2)/len(lows)

print('X的下偏方差:',semivar)

print('X的下偏标准差:，np.sqrt(semivar))

12

B=19#目标为19

lowsB=[eforeinXife<=B]

semivar_B=sum(map(lambdax:(x-B)**2,lows_B))/len(lows_B)

print('X的目标下偏方差：‘，semivarB)

print(JX的目标下偏标准差:',np.sqrt(semivar_B))

importmatplotlib.pyplotaspit

importnumpyasnp

importscipy.statsasstats

xs=np.1inspace(-6,6,300)

normal=stats,norm,pdf(xs)

pit.plot(xs,normal);

#产生数据

xs2=np.1inspace(stats,lognorm,ppf(0.01,.7,loc=~.1),stats,lognorm,ppf(0.99,.7,

loc=-.1),150)

#偏度>0

lognormal=stats,lognorm,pdf(xs2,.7)

pit.plot(xs2,lognormal,label.二'Skew>0,)

#偏度<0

pit.plot(xs2,lognormal[::-1],label=，Skew<O')

pit.legend();

#注：本程序在Bigquant环境中运行。

start='2016-01-01'

end='2018-01-01'

pricing=D.history_data(，000300.SHA',start_date=start,end_date=end,)['close']

returns=pricing,pctchange()[1:]

printCSkew:J,stats,skew(returns))

print('Mean:,np.mean(returns))

print(，Median:,np.median(returns))

pit.hist(returns,30);

fromscipyimportstats

frompandas,coreimportdatetools

importstatsmodels.apiassm#统计相关的库

importnumpyasnp

importpandasaspd

importmatplotlib.pyplotaspit

importtushareasts#财经数据接口包tushare

IndexData=ts.get_k_data(code=，hs300,,start='2016-01-01',end='2018-08-0T)

13

IndexData.index=pd.todatetime(IndexData.date)

close=IndexData.close

returns=(close-close,shift(1))/close,shift(l)

returns二returns,dropna()

print(JSkew/,stats,skew(returns))

print('Mean:,np.mean(returns))

printMedian/,np.median(returns))

pit.hist(returns,30)

pit.plot(xs,stats,laplace,pdf(xs),label='Leptokurtic，)

print(J尖峰的超额峰度：'，(stats,laplace,stats(moments=，k')))

pit.plot(xs,normal,label=,Mesokurtic(normal)')

print('正态分布超额峰度：'，(stats,norm.stats(moments=，k')))

pit.plot(xs,stats,cosine,pdf(xs),label=，Platykurtic，)

printC平峰超额峰度:',(stats,cosine,stats(moments='k')))

pit.legendO;

fromstatsmodels.stats,stattoolsimportjarque_bera

,pvalue,=jarquebera(returns)

ifpvalue>0.05:

print('沪深300收益率数据服从正态分布.')

else:

print('沪深300收益率数据并不服从正态分布.')

fields=[ffs_roe_0,]

start_date='2017-04-21‘

end_date=)2017-04-2f

instruments=D.instruments(startdate,enddate)

roe=D.features(instruments,start_date,end_date,fields=fields)['fs_roe_0']

1描述性统计

printC均值：',roe.mean())

printC标准差：',roe.std())

roe.describe()

2绘制直方图

roe.hist(bins=100)

1固定比例法

roe=D.features(instruments,start_date,end_date,fields=fields)fs_roe_0,]

roe[roe>=roe.quantile(0.99)]=roe.quantile(0.99)

roe[roe<=roe.quantile(0.01)]=roe.quantile(0.01)

print('均值：'，roe.mean())

print('标准差：'，roe.std())

14

roe.hist(bins=100)

2均值标准差方法

通常把三倍标准差之外的值都视为异常值，然后将这些异常值重新赋值

roe=D.features(instruments,start_date,end_date,fields二fields)['fs_roe_0']

roe[roe>=roe.mean()+3*roe.std()]=roe.mean()+3*roe.std()

roe[roe<=roe.mean()-3*roe.std()]=roe.mean()-3*roe.std()

print()均值：'，roe.mean())

print('标准差：'，roe.std())

roe.hist(bins=100)

3MAD方法

roe=D.features(instruments,start_date,end_date,fields=fields)['fs_roe_0']

roe=roe.dropna()

median=np.median(list(roe))

MAD=np.mean(abs(roe)-median)

roe=roe[abs(roe-median)/MAD<=6]#剔除偏离中位数6倍以上的数据

print(J均值：',roe.mean())

pr^int('标准差：'，roe.std())

roe.hist(bins=100)

fromstatsmodels.stats,stattoolsimportmedcouple

roe=D.features(instruments,startdate,enddate,fields=fields)fsroe_0,]

roe=roe.dropna()

defboxplot(data):

#mc可以使用statsmodels包中的medcouple函数直接进行计算

me=medcouple(data)

data,sort()

ql=data[int(0.25*len(data))]

q3=data[int(0.75*len(data))]

iqr=q3-ql

ifme>=0:

1=ql-1.5*np.exp(-3.5*me)*iqr

u=q3+1.5*np.exp(4*me)*iqr

else:

1=ql-1.5*np.exp(-4*me)*iqr

u=q3+5*np・exp(3.5*me)*iqr

data=pd.Series(data)

data[data<1]=1

data[data>u]=u

returndata

print(，均值'，boxplot(list(roe)).mean())

print(，标准差'，boxplot(list(roe)).std())

boxplot(list(roe)),hist(bins=100)

15

5Python参数估计

importnumpyasnp

x=[1,1,0,1,0,0,1,0,1,1,1,0,1,1,0,1,0,0,1,0,1,0,1,0,0,1,1,0,1,1,0,1]

theta=np.mean(x)

h=theta/(1-theta)

print('h=',h)

5.3Python单正态总体均值区间估计

importnumpyasnp

importscipy.statsasss

n=6;p=0.025;sigma=np.sqrt(0.6)

x=[14.6,15.1,14.9,14.8,15.2,15.1]

xbar=np.mean(x)

low=xbar-ss.norm,ppf(q=1-p)*(sigma/np.sqrt(n))

up=xbar+ss.norm,ppf(q=1-p)*(sigma/np.sqrt(n))

print('low='，low)

print('up=,,up)

2.方差b?未知时〃的置信区间

importnumpyasnp

importscipy.statsasss

fromscipy.statsimportt

n=9;p=0.025;s=np.sqrt(1.47)

x=[99.3,98.7,100.5,101.2,98.3,99.7,99.5,102.1,100.5]

xbar=np.mean(x)

low=xbar-ss.t.ppf(1-p,n-l)*(s/np.sqrt(n))

up=xbar+ss.t.ppf(1-p,n-1)*(s/np.sqrt(n))

print('low='，low)

print('up=，,up)

fromscipyimportstats

importnumpyasnp

x=[99.3,98.7,100.5,101.2,98.3,99.7,99.5,102.1,100.5]

stats,t.interval(0.95,len(x)-l,np.mean(x),stats,sem(x))

(99.04599342616191,100.90956212939363)

5.4Python单正态总体方差区间估计

fromscipy.statsimportchi2

n=16;sq=0.0023;p=0.025

low=((n-l)*sq)/chi2.ppf(1-p,n-l)

up=((n-l)*sq)/chi2.ppf(p,n-1)

print('low=,,low)

16

print('up='，up)

5.5Python双正态总体均值差区间估计

importnumpyasnp

importscipy.statsasss

x=[628,583,510,554,612,523,530,615]

y=[535,433,398,470,567,480,498,560,503,426]

nl=len(x);n2=len(y)

xbar=np.mean(x);ybar=np.mean(y)

sigmaq1=2140;sigmaq2=3250;p=0.025

low=xbar-ybar-ss.norm,ppf(q=1-p)*np.sqrt(sigmaql/nl+sigmaq2/n2)

up=xbar-ybar+ss.norm,ppf(q=1-p)*np.sqrt(sigmaql/nl+sigmaq2/n2)

print('low=，,low)

print('up='，up)

2.两方差都未知时两均值的置信区间

importnumpyasnp

importscipy.statsasss

x=[628,583,510,554,612,523,530,615]

y二[535,433,398,470,567,480,498,560,503,426]

nl=l.0*len(x);n2=l.0*len(y)#转为小数

si二np.var(x);s2=np.var(y)

xbar二np.mean(x);ybar=np.mean(y)

p=0.025

sq=((nl-1)*sl+(n2-l)*s2)/(nl-l+n2-l)

low=xbar-ybar-ss.t.ppf(1-p,nl+n2-2)*np.sqrt(sq*(l/nl+l/n2))

up=xbar-ybar+ss.t.ppf(1-p,nl+n2-2)*np.sqrt(sq*(l/nl+l/n2))

print('low=,,low)

print('up='，up)

5.6Python双正态总体方差比区间估计

importnumpyasnp

fromscipy.statsimportf

x=[20.5,19.8,19.7,20.4,20.1,20.0,19.0,19.91

y=[20.7,19.8,19.5,20.8,20.4,19.6,20.2]

sql=np.var(x);sq2=np.var(y)

nl=8;n2=7;p=0.025

f.ppf(0.025,nl-1,n2-l)

low=sql/sq2*l/f.ppf(1-p,nl-l,n2T)

up=sql/sq2*l/f.ppf(p,nl-l,n2T)

print('low='，low)

print('up='，up)

low=0.142168867371

up=4/p>

17

6Python参数假设检验

6.2Python单个样本t检验

importpandasaspd

importnumpyasnp

#读取数据并创建数据表，名称为data。

data=pd.DataFrame(pd.readexcelCG:\\2glkx\\data\\al6-l.xls'))

#查看数据表前5行的内容

data,head()

#取sale数据

x=np.array(data[sale*]])

mu=np.mean(x)

fromscipyimportstatsasss

printmu,ss.ttestIsamp(a=x,popmean=500)

6.3Python两个独立样本t检验

importnumpyasnp

#读取数据并创建数据表，名称为data。

data=pd.DataFrame(pd.readexcel('G:\\2glkx\\data\\al6-2.xls'))

#查看数据表前5行的内容

x=np.array(data[[*faJ]])

y=np.array(data[fbJ]])

fromscipy.statsimportttest_ind

t,p=ttestind(x,y)

printJt=,,t

print，p=，，p

6.4Python配对样本t检验

importpandasaspd

importnumpyasnp

#读取数据并创建数据表，名称为data。

data二pd.DataFrame(pd.read_excel('G:\\2glkx\\data\\al6-3.xls'))

#查看数据表前5行的内容

x=np.array(data[qian，]])

y=np.array(data[fhou，]])

fromscipy.statsimportttest_rel

18

t,p=ttestrel(x,y)

print't='，t

print5p=，,p

6.5Python单样本方差假设检验

importpandasaspd

importnumpyasnp

#读取数据并创建数据表，名称为data。

data=pd.DataFrame(pd.readexcel(*G:\\2glkx\\data\\al6-4.xls'))

#查看数据表前5行的内容

data,head()

#取收益率数据

importnumpyasnp

x=np.array(data[syV]])

n=len(x)

#计算方差

s2=np.var(x)

#计算卡方值

chisquare=(n-l)*s2/0.01

printchisquare

python实现卡方分布

fromscipyimportstats

obs=[102,102,96,105,95,100]

exp=[100,100,100,100,100,100]

stats,chisquare(obs,f_exp=exp)

6.6Python双样本方差假设检验

importpandasaspd

importnumpyasnp

fromscipyimportstats

fromstatsmodels.formula,apiimportols

fromstatsmodels.stats,anovaimportanova_lm

在目录G:\2glkx\data下建立al6-5.xls数据文件后，取数的命令如下:

#读取数据并创建数据表，名称为data。

df=pd.DataFrame(pd.readexcelCG:\\2glkx\\data\\al6-5.xls'))

#查看数据表前5行的内容

df.head()

19

formula=*returnA"returnB5隔离因变量和自变量（左边因变量，右边自变量）

model=ols(formula,df).fit()#根据公式数据建模，拟合

results=anova_Im(model)#计算F和P

printresults

7Python相关分析

7.2使用模拟数据计算变量之间的相关系数和绘图

#导入包

importnumpyasnp

importstatsmodels.tsa.stattoolsassts

importmatplotlib.pyplotaspit

importpandasaspd

importseabornassns

importstatsmodels.apiassm

(1)生成随机变数并绘制图形

X=np.random.randn(lOOO)

Y=np.random.randn(lOOO)

pit.scatter(X,Y)

pit.show()

print("correlationofXandYis〃)

np.corrcoef(X,Y)[0,1]

X=np.random.randn(lOOO)

Y=X+np.random,normal(0,0.1,1000)

pit.scatter(X,Y)

pit.show()

print(^correlationofXandYis〃)

np.corrcoef(X,Y)[0,1]

7.3使用本地数据计算变量之间的相关系数和绘图

importpandasaspd

importnumpyasnp

#读取数据并创建数据表，名称为data。

data二pd.DataFrame(pd.read_excel('G:\\2glkx\\data\\al7-l.xls'))

#查看数据表前5行的内容

data,head()

timeadvsale

013550

20

1250100

2356120

3468180

4570175

#取adv和sale数据

x=np.array(data[fadv*]])

y=np.array(data[sale*]])

importscipy.stats,statsasstats

r=stats.pearsonr(x,y)[0]

Print(r)

7.4使用网上数据计算变量之间的相关系数和绘图

#本程序需在Bigquant平台中运行

#计算两只股票的日收益率

#中国铁建数据

Stockl=

D.historydata([，，601186.SHA〃],startdate=，2016-12-0T,enddate='2017-05-0T,fie

Ids=fclose*!)Vclose，].pct_change()[1:]

#中国中铁数据

Stock2=

D.historydata([z，601390.SHA，Z],startdate=，2016-12-0T,enddate=，2017-05-0T,fie

Ids=['close'])['close'].pct_change()[1:]

pit.scatter(Stockl,Stock2)

pit.xlabel(，z601186.SHAdailyreturn")

pit.ylabel(，z601390.SHAdailyreturn")

pit.show()

print(^thecorrlationfortwostocksis:〃)

Stock2.corr(Stockl)

Stockl=

D.history_data([