首先，导入所需的库：numpy、pandas、matplotlib、math。代码如下：
```python
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import math
```
1. 对数据进行预处理
读取数据，并将数据按照时间顺序进行排序。代码如下：
```python
df = pd.read_csv('data.csv', delimiter='\t')
df = df.sort_values(by=['date'])
```
对于灰色马尔科夫链模型，需要对数据进行一次累加生成新的序列。代码如下：
```python
df['cumulative'] = df['close'].cumsum()
```
2. 对数据进行灰色马尔科夫链建模
根据灰色马尔科夫链模型，需要首先将原始数据序列转化为矩阵。代码如下：
```python
n = len(df)
X0 = np.array(df['cumulative'][:-1])
X1 = np.array(df['cumulative'][1:])
X1 = X1.reshape((n-1, 1))
B = np.ones((n-1, 2))
B[:, 1] = -1 * np.arange(1, n)
Y = df['close'][1:].values
```
接着，可以使用最小二乘法求解出参数a和u，并计算出残差序列e。代码如下：
```python
a, u = np.dot(np.dot(np.linalg.inv(np.dot(B.T, B)), B.T), Y)
e = Y - a*X1[:, 0] - u
```
3. 对模型预测的结果进行检验
首先，可以绘制出原始数据序列和预测序列的图像。代码如下：
```python
Y_predict = np.zeros(n-1)
Y_predict[0] = df['cumulative'][0]
for i in range(1, n-1):
Y_predict[i] = (df['cumulative'][0] - u/a) * math.exp(-a*i) + u/a
df['predict'] = np.concatenate(([df['close'][0]], np.diff(Y_predict)))
plt.plot(df['date'], df['close'], 'b-', label='Original')
plt.plot(df['date'], df['predict'], 'r-', label='Predict')
plt.legend(loc='upper left')
plt.xticks(rotation=45)
plt.show()
```
接着，可以计算出残差序列的均值、标准差和相关系数，并绘制出残差序列的图像。代码如下：
```python
mean_e = np.mean(e)
std_e = np.std(e)
corrcoef_e = np.corrcoef(e[:-1], e[1:])[0][1]
df['e'] = np.concatenate(([0], e))
plt.plot(df['date'], df['e'], 'b-')
plt.xticks(rotation=45)
plt.show()
```
最后，可以使用后验差检验来检验预测精度。代码如下：
```python
delta = np.abs(e) / Y[1:]
C = delta.mean()
P = (np.sum(delta) - delta.max()) / np.sum(delta)
Q = 1 - P
print('C: %.4f' % C)
print('P: %.4f' % P)
print('Q: %.4f' % Q)
```
4. 划分系统状态，检验所得序列是否具有马氏性
首先，需要将残差序列划分为两个状态，即正向和负向。代码如下：
```python
e_mean = np.mean(e)
df['state'] = df['e'].apply(lambda x: 1 if x >= e_mean else -1)
```
接着，可以计算出状态转移概率矩阵，并绘制出状态转移图。代码如下：
```python
P11 = np.sum(df['state'][1:] == 1) / (n-2)
P12 = 1 - P11
P21 = 1 - P11
P22 = np.sum(df['state'][1:] == -1) / (n-2)
P = np.array([[P11, P12], [P21, P22]])
print('P: ')
print(P)
plt.figure(figsize=(4, 4))
plt.imshow(P, cmap='Blues')
plt.xticks([0, 1], ['1', '-1'])
plt.yticks([0, 1], ['1', '-1'])
for i in range(2):
for j in range(2):
plt.text(j, i, '%.2f' % P[i][j], ha='center', va='center', fontsize=18)
plt.show()
```
5. 计算灰色马尔可夫链理论下的状态转移概率矩阵
根据灰色马尔科夫链模型，可以计算出灰色马尔可夫链理论下的状态转移概率矩阵。代码如下：
```python
alpha = 0.5
P_predict = np.zeros((2, 2))
P_predict[0][0] = alpha + (1-alpha)*P[0][0]
P_predict[0][1] = (1-alpha)*P[0][1]
P_predict[1][0] = (1-alpha)*P[1][0]
P_predict[1][1] = alpha + (1-alpha)*P[1][1]
print('P_predict: ')
print(P_predict)
```
6. 对灰色马尔科夫链模型进行预测
根据灰色马尔科夫链模型，可以预测出未来的状态概率分布和预测值。代码如下：
```python
state = np.zeros(n)
state[0] = df['state'][0]
for i in range(1, n):
state[i] = np.random.choice([-1, 1], p=P_predict[int(state[i-1] == 1)])
df['state_predict'] = state
df['predict_gm'] = 0
for i in range(1, n):
if df['state_predict'][i] == 1:
df['predict_gm'][i] = df['predict_gm'][i-1] + abs(df['predict'][i])
else:
df['predict_gm'][i] = df['predict_gm'][i-1] - abs(df['predict'][i])
plt.plot(df['date'], df['predict_gm'], 'r-', label='Predict')
plt.legend(loc='upper left')
plt.xticks(rotation=45)
plt.show()
```
7. 用加权灰色马尔科夫链模型进行建模
根据加权灰色马尔科夫链模型，需要首先确定权重的选择和调整。这里使用指数平均法来确定权重，并设置初始权重为0.5。代码如下：
```python
alpha = 0.5
w = np.zeros(n)
w[0] = 0.5
for i in range(1, n):
w[i] = alpha * w[i-1] + (1-alpha) * (abs(df['predict'][i]) / abs(df['e'][i]))
df['w'] = w
```
接着，根据加权灰色马尔科夫链模型，需要对数据进行二次累加。代码如下：
```python
df['cumulative2'] = df['cumulative'].cumsum()
```
接着，可以将加权灰色马尔科夫链模型转化为灰色马尔科夫链模型，并使用最小二乘法求解出参数a和u，并计算出残差序列e。代码如下：
```python
X0_w = np.array(df['cumulative2'][:-1])
X1_w = np.array(df['cumulative2'][1:])
X1_w = X1_w.reshape((n-1, 1))
Y_w = df['close'][1:].values
B_w = np.ones((n-1, 2))
B_w[:, 1] = -1 * np.arange(1, n)
W = np.diag(df['w'][1:])
a_w, u_w = np.dot(np.dot(np.dot(np.linalg.inv(np.dot(np.dot(B_w.T, W), B_w)), B_w.T), W), Y_w)
e_w = Y_w - a_w*X1_w[:, 0] - u_w
```
8. 计算加权灰色马尔可夫链理论下的状态转移概率矩阵，对加权灰色马尔科夫链模型进行预测，得到未来的预测值
根据加权灰色马尔科夫链模型，可以计算出加权灰色马尔可夫链理论下的状态转移概率矩阵，并预测出未来的预测值。代码如下：
```python
alpha_w = 0.5
P_predict_w = np.zeros((2, 2))
P_predict_w[0][0] = alpha_w + (1-alpha_w)*P[0][0]
P_predict_w[0][1] = (1-alpha_w)*P[0][1]
P_predict_w[1][0] = (1-alpha_w)*P[1][0]
P_predict_w[1][1] = alpha_w + (1-alpha_w)*P[1][1]
print('P_predict_w: ')
print(P_predict_w)
state_w = np.zeros(n)
state_w[0] = df['state'][0]
for i in range(1, n):
state_w[i] = np.random.choice([-1, 1], p=P_predict_w[int(state_w[i-1] == 1)])
df['state_predict_w'] = state_w
df['predict_gm_w'] = 0
for i in range(1, n):
if df['state_predict_w'][i] == 1:
df['predict_gm_w'][i] = df['predict_gm_w'][i-1] + abs(df['predict'][i])
else:
df['predict_gm_w'][i] = df['predict_gm_w'][i-1] - abs(df['predict'][i])
plt.plot(df['date'], df['predict_gm_w'], 'r-', label='Predict')
plt.legend(loc='upper left')
plt.xticks(rotation=45)
plt.show()
```
9. 可视化以上所有的预测结果
绘制出原始数据序列、灰色马尔科夫链模型预测序列和加权灰色马尔科夫链模型预测序列的图像。代码如下：
```python
plt.plot(df['date'], df['close'], 'b-', label='Original')
plt.plot(df['date'], df['predict_gm'], 'r-', label='Predict GM')
plt.plot(df['date'], df['predict_gm_w'], 'g-', label='Predict WGM')
plt.legend(loc='upper left')
plt.xticks(rotation=45)
plt.show()
```
904
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
