Lower bound in parameter vector causes GP to have nan values

I try to explore the likelihood surface of a Gaussian Process with regard to its parameters.
All of this is done with pystan.

If I run the following python/Stan code I sometimes get the following error:

Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Exception: multi_normal_lpdf: Covariance matrix is not symmetric. Covariance matrix[1,2] = -nan, but Covariance matrix[2,1] = -nan (in ‘/tmp/httpstan_3dvn9x_q/model_zg2xtgso.stan’, line 34, column 8 to column 32)
If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.

Interestingly, once I remove the <lower=-3.0> of my theta parameter, it works and starts to throw nonPSD errors.
(These nonPSD errors are the main reason I put a bound on the parameters in the first place)

Is there something I am doing wrong here? Is this intended behaviour and I am just overreacting?
I find it weird that those errors occur only with the bound and don’t feel that this is intended to happen.

STAN_data = 
{'N': 100, 
'D': 3, 
'x': [-4.974874496459961, -4.8743720054626465, -4.773869514465332, -4.673367023468018, -4.572864532470703, -4.472362041473389, -4.371859073638916, -4.271356582641602, -4.170854091644287, -4.070351600646973, -3.969849109649658, -3.8693466186523438, -3.7688441276550293, -3.6683413982391357, -3.5678391456604004, -3.467336654663086, -3.3668341636657715, -3.266331672668457, -3.1658291816711426, -3.065326690673828, -2.9648241996765137, -2.86432147026062, -2.7638189792633057, -2.663316488265991, -2.5628139972686768, -2.4623115062713623, -2.361809015274048, -2.2613062858581543, -2.16080379486084, -2.0603013038635254, -1.9597989320755005, -1.859296441078186, -1.7587939500808716, -1.6582913398742676, -1.5577888488769531, -1.4572863578796387, -1.3567838668823242, -1.2562813758850098, -1.1557788848876953, -1.0552763938903809, -0.9547738432884216, -0.8542713522911072, -0.753768801689148, -0.6532663106918335, -0.552763819694519, -0.4522612690925598, -0.3517587184906006, -0.25125619769096375, -0.1507536917924881, -0.05025118589401245, 0.05025133490562439, 0.15075385570526123, 0.2512563467025757, 0.3517588973045349, 0.4522612392902374, 0.5527637600898743, 0.6532662510871887, 0.753768801689148, 0.8542712926864624, 0.9547737836837769, 1.0552762746810913, 1.1557788848876953, 1.2562813758850098, 1.3567838668823242, 1.4572863578796387, 1.5577888488769531, 1.6582914590835571, 1.7587939500808716, 1.859296441078186, 1.95979905128479, 2.0603015422821045, 2.160804033279419, 2.2613065242767334, 2.361809015274048, 2.4623115062713623, 2.562814235687256, 2.6633167266845703, 2.7638192176818848, 2.86432147026062, 2.9648239612579346, 3.065326452255249, 3.1658289432525635, 3.266331434249878, 3.3668341636657715, 3.467336654663086, 3.5678391456604004, 3.668341636657715, 3.7688441276550293, 3.8693466186523438, 3.969849109649658, 4.070351600646973, 4.170854091644287, 4.271356582641602, 4.371859073638916, 4.4723615646362305, 4.572864055633545, 4.673366546630859, 4.773869037628174, 4.874371528625488, 4.974874019622803], 
'y': [-1.0819181203842163, -0.8043457269668579, 0.3109801113605499, 0.1782567799091339, 0.9251151084899902, 1.159156322479248, 1.0118178129196167, -1.263861060142517, -0.6103272438049316, 0.2958620488643646, 0.2093818038702011, 1.01191246509552, 1.1917515993118286, 0.8286805152893066, -1.3925228118896484, -0.4313831627368927, 0.26310399174690247, 0.265430212020874, 1.068880319595337, 1.2319667339324951, 0.613210916519165, -1.468259572982788, -0.24759787321090698, 0.21203632652759552, 0.320745587348938, 1.1077324151992798, 1.2762025594711304, 0.36927279829978943, -1.4970594644546509, -0.0908873975276947, 0.1682669222354889, 0.39178094267845154, 1.1335960626602173, 1.3149442672729492, 0.10947317630052567, -1.4625505208969116, 0.04500656947493553, 0.13798962533473969, 0.47033414244651794, 1.1364368200302124, 1.3364341259002686, -0.16160064935684204, -1.3814617395401, 0.158389151096344, 0.11528673022985458, 0.5607840418815613, 1.1354470252990723, 1.3368878364562988, -0.4361729025840759, -1.2612437009811401, 0.24097222089767456, 0.11090315878391266, 0.6612271070480347, 1.1329166889190674, 1.3034095764160156, -0.6954591274261475, -1.1072725057601929, 0.2919555902481079, 0.12293160706758499, 0.7663706541061401, 1.132538914680481, 1.2274646759033203, -0.9349056482315063, -0.9334529042243958, 0.3128975033760071, 0.1500050127506256, 0.8675965666770935, 1.1453498601913452, 1.1156777143478394, -1.1452274322509766, -0.7439820766448975, 0.3070670962333679, 0.18555894494056702, 0.9559884667396545, 1.167628526687622, 0.9565341472625732, -1.308119535446167, -0.554887592792511, 0.2829528748989105, 0.2271667718887329, 1.022806167602539, 1.2109051942825317, 0.765887975692749, -1.4273319244384766, -0.3651908040046692, 0.24531379342079163, 0.28091368079185486, 1.0820955038070679, 1.2519733905792236, 0.5359655618667603, -1.4846811294555664, -0.1943015307188034, 0.19299228489398956, 0.34298449754714966, 1.1204502582550049, 1.2861793041229248, 0.2922634482383728, -1.4888073205947876, -0.03836183622479439, 0.16435670852661133], 
't_mu': [-1.7920000553131104, -2.1630001068115234, 0.890999972820282], 
't_sigma': [[3.2660000324249268, 0.0, 0.0], [0.0, 2.447999954223633, 0.0], [0.0, 0.0, 2.194999933242798]]}
functions {
    array[] real softplus(array[] real v){
        array[num_elements(v)] real r;
        for (d in 1:num_elements(v)){
            r[d] = log(1.0 + exp(v[d]));
        }
        return r;
    }
    real softplus(real v){
        return log(1.0 + exp(v));
    }
}
    
data {
    int N;
    int D;
    array[N] real x;
    vector[N] y;
    vector[D] t_mu;
    matrix[D, D] t_sigma;
}
  
parameters {
    vector<lower=-3.0>[D] theta;
}
    
model {
    matrix[N, N] K;
    vector[N] mu;
    theta ~ multi_normal(t_mu, t_sigma);
    K = identity_matrix(dims(x)[1])*softplus(theta[1]) + (softplus(theta[2]) * gp_exp_quad_cov(x, 1.0, softplus(theta[3])));
    mu = zeros_vector(N);
    y ~ multi_normal(mu, K);
}

Thank you!

  • Andreas