数据分析(线性模型,协方差,相关矩阵,相关系数,多项式几何,数据平滑)

文件地址:

博客里要用的文件都在这里!

布林带

布林带由三条线组成:

中轨:移动平均线

上轨:中轨+2x5日收盘价标准差 (顶部的压力)

下轨:中轨-2x5日收盘价标准差 (底部的支撑力)

布林带收窄代表稳定的趋势,布林带张开代表有较大的波动空间的趋势。

绘制5日均线的布林带

import numpy as np
import datetime as dt
import matplotlib.pyplot as mp
# 日期转换函数
def dmy2ymd(dmy):
    dmy = str(dmy, encoding='utf-8')
    time = dt.datetime.strptime(dmy, '%d-%m-%Y').date()
    t = time.strftime('%Y-%m-%d')
    return t
dates, opening_prices, highest_prices,\
    lowest_prices, closing_prices = np.loadtxt(
        '../da_data/aapl.csv', delimiter=',',
        usecols=(1, 3, 4, 5, 6),
        dtype='M8[D], f8, f8, f8, f8',
        unpack=True, converters={1: dmy2ymd})
# 绘制收盘价的折线图
mp.figure('AAPL', facecolor='lightgray')
mp.title('AAPL', fontsize=16)
mp.xlabel('Date', fontsize=14)
mp.ylabel('closing price', fontsize=14)
mp.grid(linestyle=":")
import matplotlib.dates as md
# 拿到坐标轴
ax = mp.gca()
# 设置主刻度定位器为周定位器(每周一显示主刻度文本)
ax.xaxis.set_major_locator(
    md.WeekdayLocator(byweekday=md.MO))
ax.xaxis.set_major_formatter(
    md.DateFormatter('%d %b %Y'))
# 设置次刻度定位器为日定位器
ax.xaxis.set_minor_locator(md.DayLocator())
dates = dates.astype(md.datetime.datetime)
mp.plot(dates, closing_prices, color='dodgerblue',
        label='AAPL', linestyle='--',
        linewidth=2)
x = np.linspace(-1,0,5)
kernel = np.exp(x)[::-1]
kernel = kernel / kernel.sum()
ma53 = np.convolve(closing_prices,kernel,'valid')
mp.plot(dates[4:],ma53,color='red',label='MA-53')
# 计算上轨于下轨
stds = np.zeros(ma53.size)
for i in range(stds.size):
    stds[i] = closing_prices[i:i+5].std()
upper = ma53 + 2*stds
lower = ma53-2*stds
mp.plot(dates[4:],upper,color='blue',label= 'UPPER')
mp.plot(dates[4:],lower,color = 'orangered',label='LOWER')
mp.fill_between(dates[4:],upper,lower,lower<upper,color='green',alpha=0.2)
mp.legend()
mp.gcf().autofmt_xdate()
mp.show()

线性模型

什么是线性关系?

1 2 3 4 5

60 65 70 75 ?

线性预测

假设一组数据符合一种线型规律,那么就可以预测未来将会出现的数据。

a	b	c	d	e	f	?

ax + by + cz = d
bx + cy + dz = e
cx + dy + ez = f

a	b	c	d	e	f	g   ?


[ a b c b c d c d e ] × [ x y z ] = [ d e f ] A x B \left[ \begin{array}{ccc} a & b & c\\ b & c & d\\ c & d & e\\ \end{array} \right ] \times \left[ \begin{array}{ccc} x\\ y\\ z\\ \end{array} \right ]= \left[ \begin{array}{ccc} d\\ e\\ f\\ \end{array} \right ] \\ \quad A \quad \quad \quad \quad x \quad \quad \quad B abcbcdcde×xyz=defAxB

x = np.linalg.lstsq(A, B)[0] 

根据线性模型的特点可以通过一组历史数据求出线性关系系数x, y, z,从而预测d、e、f下的一个数据是多少。

线性预测需要使用历史数据进行检验,让预测结果可信度更高

案例:使用线性预测,预测下一天的收盘价。

# 整理五元一次方程组    最终获取一组股票走势预测值
N = 5
pred_prices = np.zeros(closing_prices.size - 2 * N + 1)
for i in range(pred_prices.size):
    a = np.zeros((N, N))
    for j in range(N):
        a[j, ] = closing_prices[i + j:i + j + N]
    b = closing_prices[i + N:i + N * 2]
    x = np.linalg.lstsq(a, b)[0]
    pred_prices[i] = b.dot(x)
# 由于预测的是下一天的收盘价,所以想日期数组中追加一个元素,为下一个工作日的日期
dates = dates.astype(md.datetime.datetime)
mp.plot(dates, closing_prices, 'o-', c='lightgray', label='Closing Price')
dates = np.append(dates, dates[-1] + pd.tseries.offsets.BDay())
mp.plot(dates[2 * N:], pred_prices, 'o-',c='orangered', 
        linewidth=3,label='Predicted Price')
mp.legend()
mp.gcf().autofmt_xdate() 
mp.show()
线性拟合

线性拟合可以寻求与一组散点走向趋势规律相适应的线型表达式方程。

有一组散点描述时间序列下的股价:

[x1, y1], [x2, y2], [x3, y3], 
...
[xn, yn]

根据线型 y=kx + b 方程可得:

kx1 + b = y1
kx2 + b = y2
kx3 + b = y3
...
kxn + b = yn

[ x 1 1 x 2 1 x 3 1 x n 1 ] × [ k b ] = [ y 1 y 2 y 3 y n ] \left[ \begin{array}{ccc} x{_1} & 1\\ x{_2} & 1\\ x{_3} & 1 \\ x{_n} & 1 \\ \end{array} \right ] \times \left[ \begin{array}{ccc} k\\ b\\ \end{array} \right ]= \left[ \begin{array}{ccc} y{_1}\\ y{_2}\\ y{_3}\\ y{_n}\\ \end{array} \right ] x1x2x3xn1111×[kb]=y1y2y3yn

样本过多,每两组方程即可求得一组k与b的值。np.linalg.lstsq(a, b) 可以通过最小二乘法求出所有结果中拟合误差最小的k与b的值。

案例:利用线型拟合画出股价的趋势线

  1. 绘制趋势线(趋势可以表示为最高价、最低价、收盘价的均值):
import numpy as np
import datetime as dt
import matplotlib.pyplot as mp
# 日期转换函数


def dmy2ymd(dmy):
    dmy = str(dmy, encoding='utf-8')
    time = dt.datetime.strptime(dmy, '%d-%m-%Y').date()
    t = time.strftime('%Y-%m-%d')
    return t

dates, opening_prices, highest_prices,\
    lowest_prices, closing_prices = np.loadtxt(
        '../da_data/aapl.csv', delimiter=',',
        usecols=(1, 3, 4, 5, 6),
        dtype='M8[D], f8, f8, f8, f8',
        unpack=True, converters={1: dmy2ymd})

# 绘制收盘价的折线图
mp.figure('AAPL', facecolor='lightgray')
mp.title('AAPL', fontsize=16)
mp.xlabel('Date', fontsize=14)
mp.ylabel('closing price', fontsize=14)
mp.grid(linestyle=":")

import matplotlib.dates as md
# 拿到坐标轴
ax = mp.gca()
# 设置主刻度定位器为周定位器(每周一显示主刻度文本)
ax.xaxis.set_major_locator(
    md.WeekdayLocator(byweekday=md.MO))
ax.xaxis.set_major_formatter(
    md.DateFormatter('%d %b %Y'))
# 设置次刻度定位器为日定位器
ax.xaxis.set_minor_locator(md.DayLocator())

dates = dates.astype(md.datetime.datetime)
mp.plot(dates, closing_prices, color='dodgerblue',
        label='AAPL', linestyle='--',
        linewidth=2,alpha = 0.4)

# 求得每天的趋势价格
trend_price = (highest_prices+lowest_prices+closing_prices)/3
mp.scatter(dates,trend_price,marker='o',color='orangered',s=15,label = 'Trend Points')
# 绘制趋势线
days = dates.astype('M8[D]').astype('int32')
A = np.column_stack((days,np.ones_like(days)))
B = trend_price
x = np.linalg.lstsq(A,B)[0]
trend_line = x[0]*days + x[1]
mp.plot(dates,trend_line,color='black',alpha=0.3,label='Trena Line')
  1. 绘制顶部压力线(趋势线+(最高价 - 最低价))
trend_points = (highest_prices + lowest_prices + closing_prices) / 3
spreads = highest_prices - lowest_prices
resistance_points = trend_points + spreads
days = dates.astype(int)
x = np.linalg.lstsq(a, resistance_points)[0]
resistance_line = days * x[0] + x[1]
mp.scatter(dates, resistance_points, c='orangered', alpha=0.5, s=60, zorder=2)
mp.plot(dates, resistance_line, c='orangered', linewidth=3, label='Resistance')
  1. 绘制底部支撑线(趋势线-(最高价 - 最低价))
trend_points = (highest_prices + lowest_prices + closing_prices) / 3
spreads = highest_prices - lowest_prices
support_points = trend_points - spreads
days = dates.astype(int)
x = np.linalg.lstsq(a, support_points)[0]
support_line = days * x[0] + x[1]
mp.scatter(dates, support_points, c='limegreen', alpha=0.5, s=60, zorder=2)
mp.plot(dates, support_line, c='limegreen', linewidth=3, label='Support')

协方差、相关矩阵、相关系数

通过两组统计数据计算而得的协方差可以评估这两组统计数据的相似程度。

样本

A = [a1, a2, ..., an]
B = [b1, b2, ..., bn]

平均值

ave_a = (a1 + a2 +...+ an)/n
ave_b = (b1 + b2 +...+ bn)/n

离差(用样本中的每一个元素减去平均数,求得数据的误差程度):

dev_a = [a1, a2, ..., an] - ave_a
dev_b = [b1, b2, ..., bn] - ave_b

协方差

协方差可以简单反映两组统计样本的相关性,值为正,则为正相关;值为负,则为负相关,绝对值越大相关性越强。

cov_ab = ave(dev_a * dev_b)
cov_ba = ave(dev_b * dev_a)

案例:计算两组数据的协方差,并绘图观察。

import numpy as np
import matplotlib.pyplot as mp

a = np.random.randint(1, 30, 10)
b = np.random.randint(1, 30, 10)
#平均值
ave_a = np.mean(a)
ave_b = np.mean(b)
#离差
dev_a = a - ave_a
dev_b = b - ave_b
#协方差
cov_ab = np.mean(dev_a*dev_b)
cov_ba = np.mean(dev_b*dev_a)
print('a与b数组:', a, b)
print('a与b样本方差:', np.sum(dev_a**2)/(len(dev_a)-1), np.sum(dev_b**2)/(len(dev_b)-1))
print('a与b协方差:',cov_ab, cov_ba)
#绘图,查看两条图线的相关性
mp.figure('COV LINES', facecolor='lightgray')
mp.title('COV LINES', fontsize=16)
mp.xlabel('x', fontsize=14)
mp.ylabel('y', fontsize=14)
x = np.arange(0, 10)
#a,b两条线
mp.plot(x, a, color='dodgerblue', label='Line1')
mp.plot(x, b, color='limegreen', label='Line2')
#a,b两条线的平均线
mp.plot([0, 9], [ave_a, ave_a], color='dodgerblue', linestyle='--', alpha=0.7, linewidth=3)
mp.plot([0, 9], [ave_b, ave_b], color='limegreen', linestyle='--', alpha=0.7, linewidth=3)

mp.grid(linestyle='--', alpha=0.5)
mp.legend()
mp.tight_layout()
mp.show()

相关系数

协方差除去两组统计样本标准差的乘积是一个[-1, 1]之间的数。该结果称为统计样本的相关系数。

# a组样本 与 b组样本做对照后的相关系数
cov_ab/(std_a x std_b)
# b组样本 与 a组样本做对照后的相关系数
cov_ba/(std_b x std_a)
# a样本与a样本作对照   b样本与b样本做对照   二者必然相等
cov_ab/(std_a x std_b)=cov_ba/(std_b x std_a)

通过相关系数可以分析两组数据的相关性:

若相关系数越接近于0,越表示两组样本越不相关。
若相关系数越接近于1,越表示两组样本正相关。
若相关系数越接近于-1,越表示两组样本负相关。

案例:输出案例中两组数据的相关系数。

print('相关系数:', cov_ab/(np.std(a)*np.std(b)), cov_ba/(np.std(a)*np.std(b)))

相关矩阵
[ v a r _ a s t d _ a × s t d _ a c o v _ a b s t d _ a × s t d _ b c o v _ b a s t d _ b × s t d _ a v a r _ b s t d _ b × s t d _ b ] [ a 与 a 的 相 关 系 数 a 与 b 的 相 关 系 数 b 与 a 的 相 关 系 数 b 与 b 的 相 关 系 数 ] \left[ \begin{array}{c} \frac{var\_a}{std\_a \times std\_a} & \frac{cov\_ab}{std\_a \times std\_b} \\ \frac{cov\_ba}{std\_b \times std\_a} & \frac{var\_b}{std\_b \times std\_b}\\ \end{array} \right ] \\ \left[ \begin{array}{c} a与a的相关系数 & a与b的相关系数 \\ b与a的相关系数 & b与b的相关系数 \\ \end{array} \right ] [std_a×std_avar_astd_b×std_acov_bastd_a×std_bcov_abstd_b×std_bvar_b][aabaabbb]
矩阵正对角线上的值都为1。(同组样本自己相比绝对正相关)
[ 1 c o v _ a b s t d _ a × s t d _ b c o v _ b a s t d _ b × s t d _ a 1 ] \left[ \begin{array}{ccc} 1 & \frac{cov\_ab}{std\_a \times std\_b} \\ \frac{cov\_ba}{std\_b \times std\_a} & 1\\ \end{array} \right ] [1std_b×std_acov_bastd_a×std_bcov_ab1]
numpy提供了求得相关矩阵的API:

# 相关矩阵
numpy.corrcoef(a, b)	
# 相关矩阵的分子矩阵 (协方差矩阵)
# [[a方差,ab协方差], [ba协方差, b方差]]
numpy.cov(a, b)			

相关系数的运用
例: 腾讯影片给你推荐你可能喜欢的电影,一种方法就是给你推荐爱好和你相似的用户喜欢的电影。那么怎么知道别人的爱好和你相似呢?通过比较别的用户和你的相关系数即可得到结果。

用户动作悬疑科幻爱情
用户14320
用户28053
用户30359

或者导演啊,演员啊之类的
通过这些去计算相关系数,即可得到用户和你的相似性

多项式拟合

为什么要用多项式拟合?
原因: 当某些情况下线性拟合产生的误差太大,这个时候我们要让线条拐弯,
多项式的一般形式:
y = p 0 x n + p 1 x n − 1 + p 2 x n − 2 + p 3 x n − 3 + . . . + p n y=p_{0}x^n + p_{1}x^{n-1} + p_{2}x^{n-2} + p_{3}x^{n-3} +...+p_{n} y=p0xn+p1xn1+p2xn2+p3xn3+...+pn
多项式拟合的目的是为了找到一组p0-pn,使得拟合方程尽可能的与实际样本数据相符合。

假设拟合得到的多项式如下:
f ( x ) = p 0 x n + p 1 x n − 1 + p 2 x n − 2 + p 3 x n − 3 + . . . + p n f(x)=p_{0}x^n + p_{1}x^{n-1} + p_{2}x^{n-2} + p_{3}x^{n-3} +...+p_{n} f(x)=p0xn+p1xn1+p2xn2+p3xn3+...+pn
则拟合函数与真实结果的差方如下
l o s s = ( y 1 − f ( x 1 ) ) 2 + ( y 2 − f ( x 2 ) ) 2 + . . . + ( y n − f ( x n ) ) 2 loss = (y_1-f(x_1))^2 + (y_2-f(x_2))^2 + ... + (y_n-f(x_n))^2 loss=(y1f(x1))2+(y2f(x2))2+...+(ynf(xn))2
那么多项式拟合的过程即为求取一组p0-pn,使得loss的值最小。
多项式的表示方法很简单只要将其的系数放在数组里面就可以了

X = [x1, x2, ..., xn] - 自变量
Y = [y1, y2, ..., yn] - 实际函数值
Y'= [y1',y2',...,yn'] - 拟合函数值
P = [p0, p1, ..., pn] - 多项式函数中的系数

根据一组样本,并给出最高次幂,求出拟合系数
np.polyfit(X, Y, 最高次幂)->P

多项式运算相关API:

根据拟合系数与自变量求出拟合值, 由此可得拟合曲线坐标样本数据 [X, Y']
np.polyval(P, X)->Y'

多项式函数求导,根据拟合系数求出多项式函数导函数的系数
np.polyder(P)->Q 

已知多项式系数Q 求多项式函数的根(与x轴交点的横坐标)
xs = np.roots(Q)

两个多项式函数的差函数的系数(可以通过差函数的根求取两个曲线的交点)
Q = np.polysub(P1, P2)

案例:求多项式 y = 4x3 + 3x2 - 1000x + 1曲线驻点的坐标。

'''
1. 求出多项式的导函数
2. 求出导函数的根,若导函数的根为实数,则该点则为曲线拐点。
'''
import numpy as np
import matplotlib.pyplot as mp

x = np.linspace(-20, 20, 1000)
y = 4*x**3 + 3*x**2  - 1000*x + 1
Q = np.polyder([4,3,-1000,1])
xs = np.roots(Q)
ys =  4*xs**3 + 3*xs**2  - 1000*xs + 1
mp.plot(x, y)
mp.scatter(xs, ys, s=80, c='orangered')
mp.show()

案例:使用多项式函数拟合两只股票bhp、vale的差价函数:

'''
1. 计算两只股票的差价
2. 利用多项式拟合求出与两只股票差价相近的多项式系数,最高次为4
3. 把该曲线的拐点都标出来。
'''
dates, bhp_closing_prices = np.loadtxt('../../data/bhp.csv', 
                                       delimiter=',',usecols=(1, 6), unpack=True, 
                                       dtype='M8[D], f8', conv erters={1: dmy2ymd})
vale_closing_prices = np.loa dtxt('../../data/vale.csv', delimiter=',',
                                 usecols=(6), unpack=True)
diff_closing_prices = bhp_closing_prices - vale_closing_prices
days = dates.astype(int)
p = np.polyfit(days, diff_closing_prices, 4)
poly_closing_prices = np.polyval(p, days)
q = np.polyder(p)
roots_x = np.roots(q)
roots_y = np.polyval(p, roots_x)
mp.figure('Polynomial Fitting', facecolor='lightgray')
mp.title('Polynomial Fitting', fontsize=20)
mp.xlabel('Date', fontsize=14)
mp.ylabel('Difference Price', fontsize=14)
ax = mp.gca()
ax.xaxis.set_major_locator(md.WeekdayLocator(byweekday=md.MO))
ax.xaxis.set_minor_locator(md.DayLocator())
ax.xaxis.set_major_formatter(md.DateFormatter('%d %b %Y'))
mp.tick_params(labelsize=10)
mp.grid(linestyle=':')
dates = dates.astype(md.datetime.datetime)
mp.plot(dates, poly_closing_prices, c='limegreen',
        linewidth=3, label='Polynomial Fitting')
mp.scatter(dates, diff_closing_prices, c='dodgerblue',
           alpha=0.5, s=60, label='Difference Price')
roots_x = roots_x.astype(int).astype('M8[D]').astype(
    		md.datetime.datetime)
mp.scatter(roots_x, roots_y, marker='^', s=80,
           c='orangered', label='Peek', zorder=4)
mp.legend()
mp.gcf().autofmt_xdate()
mp.show()

拟合需要注意的点:

  1. 拟合的次方不要过低或过高,否者会出现过拟合和欠拟合的现象
  2. 用拟合去判断数值不可以超出拟合的区间,比如用多项式拟合去预测股票的数据就是一种不可取的方法,因为未来的股票不在拟合区间里面

数据平滑

数据的平滑处理通常包含有降噪拟合等操作。降噪的功能意在****去除额外的影响因素**,拟合的目的意在数学模型化,可以通过更多的数学方法识别曲线特征。

案例:绘制两只股票收益率曲线。收益率 =(后一天收盘价-前一天收盘价) / 前一天收盘价

  1. 使用卷积完成数据降噪。
import numpy as np
import datetime as dt
import matplotlib.pyplot as mp
import matplotlib.dates as md
# 日期转换函数
def dmy2ymd(dmy):
    dmy = str(dmy, encoding='utf-8')
    time = dt.datetime.strptime(dmy, '%d-%m-%Y').date()
    t = time.strftime('%Y-%m-%d')
    return t
# 获取数据
dates,bhp_closing_price = np.loadtxt('../da_data/bhp.csv',
                             delimiter=',',
                             usecols=(1,6),
                             dtype='M8[D],f8',
                             converters={1:dmy2ymd},
                             unpack=True)

vale_closing_price = np.loadtxt('../da_data/vale.csv',
                        delimiter=',',
                        usecols=(6),
                        dtype='f8',
                        unpack=True)

# 计算收益率  收益率 =  差价/总价
bhp_returns = np.diff(bhp_closing_price)/bhp_closing_price[:-1]
vale_returns = np.diff(vale_closing_price)/vale_closing_price[:-1]
dates = dates[:-1]
# 绘制这条曲线
mp.figure('BHP VALE RETURNS', facecolor='lightgray')
mp.title('BHP VALE RETURNS', fontsize=20)
mp.xlabel('Date')
mp.ylabel('Price')
mp.grid(linestyle=':')
ax = mp.gca()
ax.xaxis.set_major_locator(md.WeekdayLocator(byweekday=md.MO))
ax.xaxis.set_minor_locator(md.DayLocator())
ax.xaxis.set_major_formatter(md.DateFormatter('%Y %m %d'))
dates = dates.astype(md.datetime.datetime)
# 绘制收益线
mp.plot(dates,bhp_returns,color='dodgerblue',linestyle='--',label='bhp returns')
mp.plot(dates,vale_returns,color='dodgerblue',linestyle='--',label='vale returns')
# 卷积降噪
convolve_core = np.hanning(8) # 得到权重
convolve_core /= convolve_core.sum()
bhp_returns_convolved = np.convolve(bhp_returns,convolve_core,'valid')
vale_returns_convolved = np.convolve(vale_returns,convolve_core,'valid')
# 绘制卷积降噪线
mp.plot(dates[7:],bhp_returns_convolved,color='orangered',alpha=0.5,label='bhp returns convolved')
mp.plot(dates[7:],vale_returns_convolved,color='orangered',alpha=0.5,label='vale returns convolved')

mp.show()
  1. 对处理过的股票收益率做多项式拟合。
days = dates.astype('M8[D]').astype('int32')
bhp_polyfit = np.polyfit(days[7:],bhp_returns_convolved,3)
bhp_polyfit_y = np.polyval(bhp_polyfit,days[7:])
vale_polyfit = np.polyfit(days[7:],vale_returns_convolved,3)
vale_polyfit_y = np.polyval(vale_polyfit,days[7:])
# 绘制拟合后的曲线
mp.plot(dates[7:],bhp_polyfit_y,color='green',linewidth=2,label='bhp poly',alpha= 0.7)
mp.plot(dates[7:],vale_polyfit_y,color='green',linewidth=2,label='vale poly',alpha= 0.7)
  1. 通过获取两个函数的焦点可以分析两只股票的投资收益比。
# 计算交叉点
diff_p = np.polysub(bhp_polyfit,vale_polyfit)
xs = np.roots(diff_p)
# 绘制这些点
roots_y = np.polyval(bhp_polyfit,xs)
xs = xs.astype('M8[D]')
mp.scatter(xs,roots_y,color= 'orangered')
  • 0
    点赞
  • 4
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值