Can't get the trace to converge in my simple polynomial fitting program. I will include an example of the 'a' parameter below. They only thing I can think of is that because I'm using logarithmic data, and the log-likelihood is causing a problem? I'm quite new to this, and I'm unsure. Hence I'm directing my question here.
This is the code I'm running:
y=vc
x=vt
y=np.asarray(y)
x=np.asarray(x)
z = np.polyfit(x, y, 2) # the traditional chi-square fit
print 'The chi-square result: ', z
#The chi-square result: [ -0.63802271 5.83847593 -13.24021273]
#priors
sig = pymc.Uniform('sig', 0, 100, value=1)
a = pymc.Uniform('a', -100, 100,value=0)
b = pymc.Uniform('b', -100, 100,value=0)
c = pymc.Uniform('c', -100, 100,value=0)
#model
@pymc.deterministic(plot=False)
def mod_quadratic(x=x, a=a, b=b, c=c):
return a*x**2 + b*x + c
#likelihood
y = pymc.Normal('y', mu=mod_quadratic, tau=1.0/sig**2, value=y, observed=True)
I have logarithmic x and y values as the input, and I fit a simple x**2 polynomial.
I keep the input file modest:
import pymc,v
# load the model file
import numpy as np
import matplotlib.pyplot as plt
R = pymc.MCMC(v) # build the model
R.sample(iter=100000,burn=20000,thin=100) # populate and run it
print 'a ', R.a.value # print outputs
print 'b ', R.b.value
print 'c ', R.c.value
print 'a ', R.a.stats()
print 'b ', R.b.stats()
print 'c ', R.c.stats()
The trace does not converge. If I use non-logarithmic data the trace converges, however the fit is not acceptable as a polynomial does not fit the non-logarithmic data.
sig trace window b trace window
To save space I have only include the b trace, but they are show the same non-convergance.
Does anyone see my flaw? I'm sure I'm doing something wrong here..
Thanks
