如何最小化基于数据价格 return 的一个分布的 kstest 提供的 python 中的 -p 值?

How can I minimize -p-value in python provided by the kstest of one distribution based on data price return?

我正在尝试从 pandas datareader 库下载股票价格并根据我提供的代码计算(每日、每周、每月等...)return。

下载数据后,我对该数据的分布执行 kstest,并根据提供的 p 值评估它是否类似于双正态分布(两个正态分布之和)。

因为我只为这个分布执行一个 kstest,所以我想利用 Python 中的 "minimize" 库来最大化 p 值(最小化 -p 值),改变均值,这两个分布的标准差和权重。

import warnings
import numpy as np
import pandas as pd
import scipy.stats as st
from scipy.optimize import minimize
import statsmodels as sm
import matplotlib
import matplotlib.pyplot as plt
from pandas_datareader import data
import time
import xlwt
import matplotlib.ticker as mtick
from sklearn import datasets

def Puxa_Preco(ticker,start_date,end_date,lag):    
    dados= data.get_data_yahoo(ticker, start_date, end_date )

    from sklearn import datasets
    data_set =  np.log(dados['Close'])-np.log(dados['Close'] .shift(lag))

    data_set = data_set.fillna(method='ffill') 
    data_set = data_set.dropna() 

    y = pd.DataFrame()
    y=data_set

    x = np.arange(len(y))
    size = len(y)
    print(y)
    return y

def mixnormal_cdf(distribuicao, weight1, mean1, stdv1,weight2, mean2, stdv2):
    """
    CDF of a mixture of two normal distributions.
    """
    return (weight1*st.norm.cdf(distribuicao, mean1, stdv1) +
            weight2*st.norm.cdf(distribuicao, mean2, stdv2))

def Objetivo(X,distribuicao):
    peso_dist_1 = X[0]
    mi1 = X[1]
    sigma1 = X[2]
    peso_dist_2 = 1-X[0]
    mi2 = X[3]
    sigma2 = X[4]

    stat2, pvalue = st.kstest(distribuicao, cdf=mixnormal_cdf,
                                args=(peso_dist_1, mi1, sigma1,peso_dist_2, mi2, sigma2))
    ''' Kolmogorov-Smirnov Test, to test whether or not the data is from a given distribution. The 
        returned p-value indicates the probability that the data is from the given distribution, 
        i.e. a low p-value means the data are likely not from the tested distribution.
        Note that, for this test, it is necessary to specify shape, location, and scale parameters,
        to obtain meaningful results (c,loc,scale). 

        stat2:     the test statistic, e.g. the max distance between the
        cumulated distributions '''

    return -pvalue

ticker = 'PETR4.SA'
start_date  = '2010-01-02'      #yyyy-mm-dd
end_date    = '2015-01-02'

for lag in range(1,503):

    distribuicao = Puxa_Preco(ticker,start_date,end_date,lag)
    n = len(distribuicao)

    ChuteInicial=[0.3,0.0010,0.0010,-0.0030,0.0830]                                      #peso_dist_1, mi1, sigma1, mi2, sigma2
    test = [0.2,0.0020,0.0110,0.8,-0.0020,0.0230]
    Limites = ([0,1],[-50,+50],[0,+50],[0,1],[-50,+50],[0,+50])                              #peso_dist_1, mi1, sigma1, peso_dist_2,mi2, sigma2
    print("------------------------------------------------------------------------------------------------")
    print("Validation Test:")

    print(-Objetivo(test,distribuicao))                                             #the value should be around -0.90 to verify that the objective function it is ok

    solution = minimize(fun=Objetivo,x0=ChuteInicial,args=distribuicao,method='SLSQP',bounds = Limites)             #minimize - p-valor
    print("------------------------------------------------------------------------------------------------")
    print("solution:")
    print(solution)

找到以下解决方案:

         fun: -8.098252265651002e-53
         jac: array([-2.13080032e-35,  0.00000000e+00,  0.00000000e+00, -1.93307671e-34, 7.91878934e-35])
     message: 'Optimization terminated successfully.'
        nfev: 8
         nit: 1
        njev: 1
      status: 6
     success: True
           x: array([ 0.3  ,  0.001,  0.001, -0.003,  0.083])

但我知道正确答案应该类似于 (test) : [0.2,0.0020,0.0110,0.8,-0.0020,0.0230] 产生的 p 值为 0.90

在我看来,这只是 运行 几次模拟,因为它没有改变 p 值,所以它停止了。

有没有一种方法可以确保 "minimize" 只有在找到大于 0.9 的 p 值后才会停止? 有人可以帮我吗?

我尝试使用最小化考虑 Nelder Mead 方法,看起来更准确但甚至不接近应该是答案的 0.9 p 值,我不知道 Nelder Mead 方法是否考虑了我的限制提供。

#solution = minimize(fun=Objetivo,x0=(ChuteInicial),args=distribuicao,method='Nelder-Mead',bounds = Limites,options={'int':1000000})            

通过在定义 cdf 函数时能够最小化 k-s 统计量而不是 p 值和其他修改,我想我能够估计参数。这是我的代码和优化的参数估计。我从本文 (https://www.researchgate.net/publication/250298633_Minimum_Kolmogorov-Smirnov_test_statistic_parameter_estimates)

中得到最小化 k-s 统计量的想法
import warnings
import numpy as np
import pandas as pd
import scipy.stats as st
from scipy.optimize import minimize
import statsmodels as sm
import matplotlib
import matplotlib.pyplot as plt
from pandas_datareader import data
import time
import xlwt
import matplotlib.ticker as mtick
from sklearn import datasets

def Puxa_Preco(ticker,start_date,end_date,lag):    
    dados= data.get_data_yahoo(ticker, start_date, end_date )

    data_set =  np.log(dados['Close'])-np.log(dados['Close'] .shift(lag))

    data_set = data_set.fillna(method='ffill') 
    data_set = data_set.dropna() 

    y = pd.DataFrame()
    y=data_set

    x = np.arange(len(y))
    size = len(y)

    return y

def mixnormal_cdf(distribuicao, weight1, mean1, stdv1,mean2, stdv2):
    ## CDF of a mixture of two normal distributions.

    return (weight1*st.norm.cdf(distribuicao, mean1, stdv1) +
            (1-weight1)*st.norm.cdf(distribuicao, mean2, stdv2))

def Objetivo(X, distribuicao):
    stat2, pvalue = st.kstest(distribuicao, cdf=mixnormal_cdf, args=X)
    return stat2

ticker = 'PETR4.SA'
start_date  = '2010-01-02'      #yyyy-mm-dd
end_date    = '2015-01-02'

for lag in range(1,10):
    distribuicao = Puxa_Preco(ticker,start_date,end_date,lag)

    ChuteInicial=[0.2,0.0020,0.0110,-0.0020,0.0230]                                      #peso_dist_1, mi1, sigma1, mi2, sigma2

    Limites = ([0,1],[-0.1,+0.1],[0,None],[-0.1,0.1],[0,None])                              #peso_dist_1, mi1, sigma1, peso_dist_2,mi2, sigma2

    print("------------------------------------------------------------------------------------------------")
    print("Intial Guess")
    print(" weight1: {:.4f}, mean1: {:.4f}, stdv1: {:.4f}, mean2: {:.4f}, stdv2: {:.4f}".format(*ChuteInicial))
    ks_stat, p_value = st.kstest(distribuicao, cdf=mixnormal_cdf, args=ChuteInicial)
    print("k-s stat: {}, p-value: {}".format(ks_stat, p_value))
    print("------------------------------------------------------------------------------------------------")

    solution = minimize(fun=Objetivo,x0=ChuteInicial,args=(distribuicao), method='SLSQP',bounds = Limites, options={'ftol':1e-12})

    print("Optimized Solution:")
    print("------------------------------------------------------------------------------------------------")
    print(solution.message)
    ks_stat, p_value = st.kstest(distribuicao, cdf=mixnormal_cdf, args=solution.x)
    print("Optimized k-s stat: {}, p-value: {}".format(ks_stat, p_value))
    print("Optimized weight1: {:.4f}, mean1: {:.4f}, stdv1: {:.4f}, mean2: {:.4f}, stdv2: {:.4f}".format(*solution.x))
    print("------------------------------------------------------------------------------------------------")

这导致

 LAG: 1
------------------------------------------------------------------------------------------------
Intial Guess
 weight1: 0.2000, mean1: 0.0020, stdv1: 0.0110, mean2: -0.0020, stdv2: 0.0230
k-s stat: 0.016419395777755863, p-value: 0.9027260234690881
------------------------------------------------------------------------------------------------
Optimized Solution:
------------------------------------------------------------------------------------------------
Optimization terminated successfully.
Optimized k-s stat: 0.014896801186217778, p-value: 0.952748910069967
Optimized weight1: 0.2016, mean1: 0.0016, stdv1: 0.0108, mean2: -0.0017, stdv2: 0.0227
------------------------------------------------------------------------------------------------
LAG: 2
------------------------------------------------------------------------------------------------
Intial Guess
 weight1: 0.2000, mean1: 0.0020, stdv1: 0.0110, mean2: -0.0020, stdv2: 0.0230
k-s stat: 0.0791846286223467, p-value: 5.496852766236095e-07
------------------------------------------------------------------------------------------------
Optimized Solution:
------------------------------------------------------------------------------------------------
Optimization terminated successfully.
Optimized k-s stat: 0.013690695822088705, p-value: 0.978164746056248
Optimized weight1: 0.2291, mean1: -0.0047, stdv1: 0.0125, mean2: -0.0012, stdv2: 0.0351
------------------------------------------------------------------------------------------------
LAG: 3
------------------------------------------------------------------------------------------------
Intial Guess
 weight1: 0.2000, mean1: 0.0020, stdv1: 0.0110, mean2: -0.0020, stdv2: 0.0230
k-s stat: 0.13436209329730314, p-value: 2.533666062868497e-19
------------------------------------------------------------------------------------------------
Optimized Solution:
------------------------------------------------------------------------------------------------
Optimization terminated successfully.
Optimized k-s stat: 0.01622740259609401, p-value: 0.9106251238446841
Optimized weight1: 0.2638, mean1: -0.0058, stdv1: 0.0199, mean2: -0.0020, stdv2: 0.0433
------------------------------------------------------------------------------------------------

我找到了一种方法来确保这一点,但这不是最聪明的方法...

  for lag in range(1,10):
     distribuicao = Puxa_Preco(ticker,start_date,end_date,lag)
     delta1 = 0
     delta2 = 0
     delta3 = 0
     for j in range(1,101): 
        ChuteInicial=[0.05+delta1,0.0010+delta2,0.0010+delta3,-0.0030+delta2,0.0830+delta3]                                      #peso_dist_1, mi1, sigma1, mi2, sigma2

        Limites = ([0,1],[-0.5,+0.5],[0,None],[-0.5,0.5],[0,None])                              #peso_dist_1, mi1, sigma1, peso_dist_2,mi2, sigma2

        print("------------------------------------------------------------------------------------------------")
        print("------------------------------------Lag " + str(lag) + " -------------------------------------")
        print("------------------------------------------------------------------------------------------------")
        print("Intial Guess")
        print(" weight1: {:.4f}, mean1: {:.4f}, stdv1: {:.4f}, mean2: {:.4f}, stdv2: {:.4f}".format(*ChuteInicial))
        ks_stat, p_value = st.kstest(distribuicao, cdf=mixnormal_cdf, args=ChuteInicial)
        print("k-s stat: {}, p-value: {}".format(ks_stat, p_value))
        print("------------------------------------------------------------------------------------------------")

        solution = minimize(fun=Objetivo,x0=ChuteInicial,args=(distribuicao), method='SLSQP',bounds = Limites)             #minimize - p-valor

        print("Optimized Solution:")
        print("------------------------------------------------------------------------------------------------")
        print(solution)
        ks_stat, p_value = st.kstest(distribuicao, cdf=mixnormal_cdf, args=solution.x)
        print("Optimized k-s stat: {}, p-value: {}".format(ks_stat, p_value))
        print("Optimized weight1: {:.4f}, mean1: {:.4f}, stdv1: {:.4f}, mean2: {:.4f}, stdv2: {:.4f}".format(*solution.x))
        print("------------------------------------------------------------------------------------------------")
        if ((p_value <= 0.9) or (ks_stat >= 0.02)):
          delta1 += 0.005
          delta2 += 0.001
          delta3 += 0.001
        else:

          break