추적 값에 안정성 (원치 않는) 기간이있는 이유는 무엇입니까?

Question

나는 pymc3에 맞추려고 노력하고있는 상당히 간단한 테스트 데이터 세트를 가지고있다. traceplot에 의해 생성

결과는 본질적으로 this. 다시 애벌레 다음에 750 반복에 대한 평평한 선으로 이어 100 반복을위한 표준 '애벌레',이 같은 모든 매개 변수 모양의 흔적 같은 것을 보인다.

초기 100 회 반복은 25,000 회 ADVI 반복 및 10,000 회 반복 반복 후에 발생합니다. 이 금액을 변경하면 나는 이러한 원치 않는 안정성을 무작위로 갖거나 갖지 않을 것입니다.

아무도 내가 이것을 어떻게 막을 수 있는지에 대한 조언이 있다면 궁금합니다. 그 원인은 무엇입니까?

감사합니다.

전체 코드는 다음과 같습니다. 간단히 말해서, 나는 y = a (j) * sin (phase) + b (j) * sin (phase)의 값을 갖는 'phase'(-pi -> pi)를 생성하고있다. a 및 b는 각 피사체 j에 대해 무작위로 추출되지만 서로 관련되어 있습니다. 그런 다음 본질적으로이 동일한 모델에 맞추려고합니다.

편집 : Here is a similar example, running for 25,000 iterations. Something goes wrong around iteration 20,000.

#!/usr/bin/env python3 
# -*- coding: utf-8 -*- 
import matplotlib.pyplot as plt 
import numpy as np 
import pymc3 as pm 
%matplotlib inline 

np.random.seed(0) 


n_draw = 2000 
n_tune = 10000 
n_init = 25000 
init_string = 'advi' 
target_accept = 0.95 

## 
# Generate some test data 
# Just generates: 
# x a vector of phases 
# y a vector corresponding to some sinusoidal function of x 
# subject_idx a vector corresponding to which subject x is 

#9 Subjects 
N_j = 9 

#Each with 276 measurements 
N_i = 276 
sigma_y = 1.0 
mean = [0.1, 0.1] 
cov = [[0.01, 0], [0, 0.01]] # diagonal covariance 

x_sub = np.zeros((N_j,N_i)) 
y_sub = np.zeros((N_j,N_i)) 
y_true_sub = np.zeros((N_j,N_i)) 
ab_sub = np.zeros((N_j,2)) 
tuning_sub = np.zeros((N_j,1)) 
sub_ix_sub = np.zeros((N_j,N_i)) 
for j in range(0,N_j): 
    aj,bj = np.random.multivariate_normal(mean, cov) 
    #aj = np.abs(aj) 
    #bj = np.abs(bj) 
    xj = np.random.uniform(-1,1,size = (N_i,1))*np.pi 
    xj = np.sort(xj)#for convenience 
    yj_true = aj*np.sin(xj) + bj*np.cos(xj) 
    yj = yj_true + np.random.normal(scale=sigma_y, size=(N_i,1)) 

    x_sub[j,:] = xj.ravel() 
    y_sub[j,:] = yj.ravel() 
    y_true_sub[j,:] = yj_true.ravel() 

    ab_sub[j,:] = [aj,bj] 
    tuning_sub[j,:] = np.sqrt(((aj**2)+(bj**2))) 
    sub_ix_sub[j,:] = [j]*N_i 

x = np.ravel(x_sub) 
y = np.ravel(y_sub) 
subject_idx = np.ravel(sub_ix_sub) 
subject_idx = np.asarray(subject_idx, dtype=int) 

## 
# Fit model 
hb1_model = pm.Model() 
with hb1_model: 
# Hyperpriors 
    hb1_mu_a = pm.Normal('hb1_mu_a', mu=0., sd=100) 
    hb1_sigma_a = pm.HalfCauchy('hb1_sigma_a', 4) 
    hb1_mu_b = pm.Normal('hb1_mu_b', mu=0., sd=100) 
    hb1_sigma_b = pm.HalfCauchy('hb1_sigma_b', 4) 

    # We fit a mixture of a sine and cosine with these two coeffieicents 
    # allowed to be different for each subject 
    hb1_aj = pm.Normal('hb1_aj', mu=hb1_mu_a, sd=hb1_sigma_a, shape = N_j) 
    hb1_bj = pm.Normal('hb1_bj', mu=hb1_mu_b, sd=hb1_sigma_b, shape = N_j) 

    # Model error 
    hb1_eps = pm.HalfCauchy('hb1_eps', 5) 

    hb1_linear = hb1_aj[subject_idx]*pm.math.sin(x) + hb1_bj[subject_idx]*pm.math.cos(x) 
    hb1_linear_like = pm.Normal('y', mu = hb1_linear, sd=hb1_eps, observed=y) 

with hb1_model: 
    hb1_trace = pm.sample(draws=n_draw, tune = n_tune, 
         init = init_string, n_init = n_init, 
        target_accept = target_accept) 

pm.traceplot(hb1_trace) 

Answer 1

해밀턴 몬테 카를로에 노트와 문학에 설명 된대로 내가, 예를 들어, here 및 here를 참조의 divergencies보고 싶은데.

with model: 
    np.savetxt('diverging.csv', hb1_trace['diverging']) 

더러운 해결책으로, 아마도 target_accept을 늘려보십시오.

행운을 빈다.

Answer 2

부분적으로 내 자신의 질문에 대답하기 : 문제가있을 수 있습니다처럼 잠시 동안 함께 연주 한 후,이 때문에 내가 확실하지 오전 0으로가는 hyperprior 표준 편차를 찾습니다 왜 알고리즘이 거기서 붙어 있어야만 하는가? (작은 표준 편차를 테스트하는 것은 드문 일이 아니다 ...). 어떤 경우

은 (그들이 완전히 제거하지 않지만) 문제를 완화하는 것 두 가지 솔루션은 다음과 같습니다

1) 표준 편차의 정의에 오프셋을 추가합니다. 예컨대 :

offset = 1e-2 
hb1_sigma_a = offset + pm.HalfCauchy('hb1_sigma_a', 4) 

2) 대신 전에 SD 용 HalfCauchy 또는 HalfNormal를 사용하는 가능성이 0이되도록 설정 로그 정규 분포를 사용한다.