Python商业数据分析代码

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

importpandasaspd

#读取本地数据

df=pd.readcsvC,/2glkx/data/al2-l.csv*)

#读取网上数据

importpandasaspd

data_url=

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

/Python/Numorical-Dcscriptions-of-the-Data/data.csv”

df=pd.read_csv(data_url)

print(df.head())

print(df.tai1())

print(df.columns)

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

#Extractingrownamesortheindex

print(df.index)

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

#Transposedata

print(df.T)

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

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

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

print(df.describe0)

fromscipyimportstatsasss

#Performonesamplet-testusing1500asthetruemean

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

TtestlsampResult(statistical.1281738488299586,pvalue=0.26270472069109496)

print(ss.ttcst_lsamp(a=df,popmean=15000))

importmatplotlib.pyplotaspit

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

1

importseabornassns

#Dotheboxplot

pit.show(sns.boxplot(df))

defadd_2int(x,y):

returnx+y

print(add_2int(2,2))

2Python商业数据的存取

importpandasaspd

importnumpyasnp

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=dict(zip(a,b))

d

p=pd.DataFramc(d)

P

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

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

importpandasaspd

importnumpyasnp

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

df.headO

importpandasaspd

importnumpyasnp

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

data=pd.DataFrame(pd.readexcel(JG:\\2glkx\\data\\al2-2.xls,))

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

data,head()

importtushareasts

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

#直接保存

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

#选择数据保存

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

importpandasaspd

importnumpyasnp

df=pd.readcsv(*F:/2glkx/data/000875.csv')

2

df.head()

importtushareasts

importos

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

forcodeinf000875,,*600848,,'000981']:

df=ts.gethistdata(code)

ifos.path,exists(filename):

itdf.to_csv(filename,mode='a',header二None)

df.to_csv(filoname,modo=,a')

else:

df.to_csv(filename)

importpandasaspd

importnumpyasnp

df=pd.rcadexccl('F:/2glkx/data/000875.xls')

df.head()

importtushareasts#需先安装tushare程序包

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

importpandasaspd

importnumpyasnp象中

data=pd.DataFrame()

datal=ts.gct_hist_data(，600000,)

并需要修改上面的时间

datal=datalfclose,]

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

dataf600000']=datal

data2=ts.got_hist_data(，0009801)

data2=data2fclose*]

data2=data2[:：-l]

data['000980']=data2

data3=ts.get_hist_data(，000981,)

data3=data3['close']

data3=data3[:：-l]

dataf00098T]=data3

data,info0#查看数据情况

data=dala.dropna0

data,info()

data,head0

data,tail()

3

data=dataEf600000),‘000981']]

data,head()

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

data.iloc[:2,:3]

importtushareasts

importpandasaspd

pd.setoption(*expand_frame_repr，,False)#显示所有列

ts.settokenCyourtoken*)#获取token号，需要先注册

pro=_api0

stockdata=pro.daily(ts_code=,000001.SZ*,startdate=，20100101，,

end_date=，20190101')

stockdata,head()

importpandasdatareader.dataasweb

importdatetime

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

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

stock=web.DataReader(，，600797.SS^,“yahoo”,start,end)

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

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

importpandasaspd

frompandasimportSeries,DataFramc

importnumpyasnp

importmatplotlib.pyplotaspit

frompandas_datareaderimportdata,wh

fromdatetimeimportdatetime

end=datetime.no\v()

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

df=data.DataReadcrC*600797.SS*,'yahoo',start,end)

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

pit.show()

importnumpyasnp

importpandasaspd

importpandasdatareader.dataasweb

importdatetime

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

4

dfcsvsave=

web.DataReader(“600018.SS","yahoo”,datetime.datetime(2019,1,1),datetime,date.to

day())

print(df_csvsave)

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

True)

importtushareasts

importpandasaspd

pd.set_option(，expand_frame_repr,,False)#显示所有列

ts.set_token(，yourtoken')#获取token号

pro=ts.proapi()

code」ist=f000001.SZ\*600000.SH1,'()00002.ST]

stockdata=pd.DataFrameO

forcodeincode_list:

print(code)

df=pro.daily(ts_codc=codc,start_date=，2018010T,end_date=,20180104))

stock_data=stock_datci.cippend(df,ignore_index=True)

print(stockdata)

grouped=stock_data.groupby(*ts_code')

print(grouped)

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

importmatplotlib.pyplotaspit

importpandasaspd

importnumpyasnp

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

#或者df=pd.read_excelCG:\\2glkx\\data\\al3-l.xls')

df.headO

fig=plt.figureO

ax=fig.addsubplot(1,L1)

ax.hist(df[JAge*],bins=7)

pit.show()

fig=plt.figureO

ax=fig.addsubplot(1,L1)

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

pit.titleCAgedistribution')

5

pit.show()

fig=plt.figure()

ax=fig.add_subplot(l,1.1)

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

plt.titleC1Agedistribution')

pit.xlabel(,Age*)

pit.ylabel(JttEmployee*)

pit.show()

fig=plt.figure()

ax二fig.add_subplot(l,1,1)

ax.scatter(dffAge*],df['Sales'])

#Youcanalsoaddmorevariablesheretorepresentcolorandsize.

pit.title(*Age&SalesScatterofEmployee1)

^Variable

pit.xlabel(*Age*)

pit.ylabel(*Sales,)

pit.show()

fig=plt.figure()

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

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

#Addedthirdvariableincomeassizeofthebubble

pit.xlabel(JAge*)

pit.ylabel('Sales*)

pit.show()

importmatplotlib.pyplotaspit

importpandasaspd

fig二pit.figure()

ax=fig.addsubplot(1,L1)

^Variable

ax.boxplot(dffAge1])

pit.BoxfigureofAge')

pit.show()

vars=[*Age'Sales']

data=df[vars]

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

importmatplotlib.pyplotaspit

importpandasaspd

6

importnumpyasnp

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

df.head()

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

temp=var.unstack()

x_list=temp[，Sales*]

label_list=temp.index

pit.axis("equal")

pit.pie(xlist)

pit.titie(，zPastafatianismexpenses")

pit.show()

frompylabimport*

#makeasquarefigureandaxes

figure(1,figsize=(6,6.))

ax=axcs([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=，%l.If%%,,shadow=True,

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

title(*RateofMaleandFemale*)#标题

show()

importmatplotlib.pyplotasplot

importpandasaspd

importnumpyasnp

df=pd.readexcel("G:/2glkx/data/a13T.xls")

df.headO

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

^groupedsumofsalesatGenderlevel

fig=pit.figure0

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

axl.set_xlabelCGender")

axl.set_ylabel('SumofSales')

axl.settitie(^GenderwiseSumofSales")

var.plot(kind='bar*)

var=df.groupby(fBMT,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(*BMT)

axl.sctylabcl(*SumofSales*)

axl.set_title(，，BMIwiseSumofSales")

var.plot(kind』line')

importpandasaspd

importnumpyasnp

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

data,head()

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

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

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.titleC5populationcensus')

pl.xlabel('Time')

pl.ylabel('PopulatioiV)

pl.show()

importpandasaspd

importnumpyasnp

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

data,head()

t=np.array(data[fyear,]])

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

importpylabaspl

pl.plot(t,x)

pl.title(>1998-2015ofAlistedcompaniesinchina))

pl.xlabelCTime*)

pl.ylabel('companiesnumbers*)

8

pl.show()

importpandasaspd

importnumpyasnp

importpylabaspl

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

data,head()

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

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

importpylabaspl

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

pl.title(>1998-2015ofAlistedcompaniesinchina,)

pl.xlabel('Time')

pl.ylabelCcompaniesnumbers,)

pl.show()

importrandom

importnumpyasnp

importmatplotlibasmpl

importmatplotlib.pyplotaspit

importmatplotlib.datesasmdates

frommpl_toolkits.mplol3dimportAxes3I)

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

fig=pit.figure()

ax=fig.add_subplot(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.setmajorlocator(mpl.ticker.FixedLocator(xs))

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

ax.setxlabel('Month')

ax.set_ylabel('Year*)

ax.setzlabel(*SalesNet[usd]*)

pit.show()

frommpltoolkits.mplot3dimportAxes3D

frommatplotlibimportcm

importmatplotlib.pyplotaspit

importnumpyasnp

nangles=36

n_radii=8

#Anarrayofradii

9

#Doosnotincluderadiusr=0,thisistoeliminateduplicatepoints

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

#Anarrayofangles

angles=np.1inspace(0.2*np.pi,n_cingles,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(tingles)).flatten())

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

#Pringlesurface

z=np.sin(-x*y)

fig=pit.figureO

ax=fig.gca(projection^3d1)

ax.plot_trisurf(x,y,z,cmap=cm.jet,1inewidth=O.2)

pit.show()

4Python描述性统计

#两个常用的统计包

importscipy.statsasstats

importnumpyasnp

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

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的中位数：'，rp.median(xl))

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

print(,Onemodeofxl,stats,mode(xl)[0][0])

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

defmode(l):

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

counts={}

forein1:

ifeincounts:

counts[e]+=1

else:

10

counts[e]=1

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

maxcount=0

modes-{)

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

ifvalue>maxcount:

maxcount=value

modes={key}

clifvalue==maxcount:

modes,add(key)

ifmaxcount>1orlen(1)==1:

returnlist(modes)

return'Nomode,

print(*Allofthemodesofxl,mode(xl))

importscipy.statsasstats

importnumpyasnp

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

start=*2014-01-0r

end='2015-01-01'

pricing=D.historydata('000002.SZA',fields=[,close*],start_date=start,

end_date=end)['close*]

returns=pricing.pct_changc()[1:]

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

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

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

maxfreq=max(hist)

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

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

maxfreq])

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

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

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

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

importscipy.statsasstats

importnumpyasnp

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

11

R_G=stats,gmcan(ratios)-1

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

T=len(returns)

initprice=pricingfO.

final_price=pricing[T]

print('最初价格:，,initprice)

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

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

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

print(1x2的调和平均值:'，stats.hmean(x2))

importnumpyasnp

np.random.seed(121)

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

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

#Sortthem

X=np.sort(X)

printCX:%s'%(X))

mu=np.mean(X)

print('X的平均值:mu)

print(*RcingeofX:%s'%(np.ptp(X)))

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

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

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

printCX的方差：'，np.var(X))

print('X的标准差：'，np.std(X))

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

dist=k*np.std(X)

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

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

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(1X的目标下偏方差：'，semivar_B)

print('X的目标下偏标准差：’，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>O')

#偏度<0

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

pit.legendO;

并注：本程序在Bigquanl环境中运行。

start='2016-01-01'

end='2018-01-01'

pricing=D.historydata(1000300.SHA',startdate=start,enddate=end,)['close*]

returns=pricing.pct_chango()[1:]

print(*Skew:),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

importtushareasIs#财经数据接口包tushare

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

13

IndcxData.index=pd.to_datetimc(IndexData.date)

close=IndexDcita.close

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

returns=returns.dropnaO

print('Skew:*,stats,skew(returns))

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

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

pit.hist(returns,30)

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

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

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

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

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

print('平峰超额峰度：'，(stats,cosine,stats(moments=，k*)))

pit.legend();

fromstatsmodels.stats,stattoolsimportjarque_bera

pvalue,=jarquebera(returns)

ifpvalue>0.05:

printC沪深300收益率数据服从正态分布.’)

else:

printC沪深300收益率数据并不服从正态分布.’)

fields=[,fs_roc_0,]

start_date=>2017-04-21,

enddate=J2017-04-21

instruments=I),instruments(start_date,enddate)

roe=D.features(instruments,startdate,enddate,fields二fields)['fsroe01]

1描述性统计

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

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

roe.describe()

2绘制直方图

roe.hist(bins=100)

1固定比例法

roe=D.features(instruments,startdatc,enddate,fields=fields)[*fspoeOl]

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=0.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.stdO

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

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

roo.hist(bins=100)

3MAI)方法

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

roe=roe.dropnaO

median=np.median(list(roe))

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

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

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

print('标准差：'，roe.$td())

roe.hist(bins=100)

fromstatsmodels.stats,stattoolsimportmedcouple

roe=D.features(instruments,start_date,enddate,fields=fields)fs_roe_01]

roe=roe.dropnaO

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+1.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())

printC标准差',boxplot(list(roe)).std())

boxplot(1ist(roe)).hist(bins=100)

15

5Python参数估计

importnumpyasnp

x=[l,1,0,1,0,0,1,0,1,LI,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/(l-theta)

printCh=',h)

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

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.方差。2未知时〃的置信区间

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(l-p,n-l)*(s/np.sqrt(n))

up=xbar+ss.t.ppf(1-p,n-l)*(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-l)

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)

sigmaql=2140;sigmaq2=3250;p=0.025

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

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

print(*low=,,low)

printup=，,up)

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

importnumpyasnp

importscipy.statsusss

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)#转为小数

sl=np.var(x);s2=np.var(y)

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

p=0.025

sq=((nl-1)*s1+(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)

printup=,,up)

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

importnumpyasnp

fromscipy.statsimportf

x=[20.5,19.8,19.7,20.4;20.1,20.0,19.0,19.93

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,n2-l)

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

printClow=5,low)

print('up=',up)

low=0.142168867371

up=4/p>

17

6Python参数假设检验

6.2Python单个样本t检验

importpandasaspd

importnumpyasnp

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

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

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

data,head()

#取sale数据

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

mu=np.mean(x)

fromscipyimportstatsasss

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

6.3Python两个独立样本t检验

importnumpyasnp

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

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

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

x=np.array(datafa*]])

y=np.array(data[['fb']])

fromscipy.statsimportttestind

t,p=ttest_ind(x,y)

print't=',t

print，p=，,p

6.4Python配对样本t检险

importpandasaspd

importnumpyasnp

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

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

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

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

y=np.array(datahou，]])

fromscipy.statsimportttest_rel

18

t,p=ttest_rel(x,y)

print't='，t

print，p=，,p

6.5Python单样本方差假设检脸

importpandasaspd

importnumpyasnp

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

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

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

data,head0

#取收益率数据

importnumpyasnp

x=np.array(data[[,syl*]])

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,1D0,100,100,100]

stats.chisquare(obs,fexp=exp)

6.6Python双样本方差假设检验

importpandasaspd

importnumpyasnp

fromscipyimportstats

fromstatsmodels.formula.apiimportols

fromstatsmodels.stats.anovaimportanovalm

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

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

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

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

df.head()

19

formula='rcturnA"returnB，《隔离因变量和自变量（左边因变量，右边自变展）

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(lODO)

Y=np.random.randn(lODO)

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,head0

timeadvsale

013550

20

1250100

2356120

3468180

4570175

#取adv和sale数据

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

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

importscipy.stats,statsasstats

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

Print(r)

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

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

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

#中国铁建数据

Stockl=

D.history_data(["601186.SHA"],start_date=，2016-12-0T,end_dato=,2017-05-01,,fie

Ids=fclose*])Vclose1].pct_change()[1:]

#中国中铁数据

Stock2=

D.history_data([z，60139D.SHA"],start_date=2016-12-01',end_date=,2017-05-01,,fie

Ids=fclose*])close*].pct_change()[1:]

pit.scatter(Stockl,Stock2)

pit.xlabel(z，601186.SHAdailyreturn*)

plI.ylabel(“601390.SHAdailyreturn*)

pit.show()

print(^thecorrlationfortwostocksis:")

Stock2.corr(