3

Update: solved! It is producing parameters with the correct signs now, and they do fit the curve. The problem was defining func(a,b,c,x) but curve_fit needs to read x first: func(x,a,b,c). Thanks everyone for all the help! I'll have quantitative analysis when I meet with my boss today :)

Here's some of the new fits: http://imgur.com/NHnzR2A

(I still get a run-time error:

RuntimeWarning: overflow encountered in power
return a*(math.e**(b*(math.e**(c*x))))

)

Can anyone help me figure out what's wrong with this code? I am new to scipy. I am trying to model bacterial growth with the Gompertz equation, but my code produces a curve_fit that's completely wrong. You can see images of my real data, the model equation, and the fit this code produces in this imgur album Thanks!

Fixed code:

#!/usr/bin/python
from numpy import *
from scipy.optimize import curve_fit

values = numpy.asarray(values)  
y = values[:2000//5].astype(numpy.float)
y - y[0] #subtracting blank value
x = numpy.arange(len(y))*5

def Function(x,a,b,c):
  #a = upper asymptote
  #b = negative = x axis displacement
  #c = negative = growth rate
  return a*(math.e**(b*(math.e**(c*x))))

parameters, pcov = curve_fit(Function, x, y, p0=[0.1,-1300,-0.0077])

#Graph data and fit to compare
yaj = Function(  numpy.asarray(x), parameters[0], parameters[1], parameters[2] )
figure(1, figsize=(8.5,11))
subplot(211)
plot(x,y,'g-')
xlim(min(x),max(x))
ylim(min(y),max(y))
subplot(212)
plot(x,yaj,'r-')
xlim(min(x),max(x))
ylim(min(yaj),max(yaj))
savefig('tempgraph.pdf')

return parameters
4

1 回答 1

3

进口： 

import numpy as np
import matplotlib.pyplot as plt
import scipy.optimize as opt

样本值：

values = np.array('0.400    0.400   0.397   0.395   0.396   0.394   0.392   0.390   0.395   0.393   0.392   0.392   0.390   0.388   0.390   0.388   0.385   0.383   0.388   0.387   0.387   0.387   0.385   0.386   0.387   0.379   0.379   0.378   0.375   0.376   0.374   0.373   0.372   0.368   0.373   0.370   0.371   0.370   0.370   0.370   0.367   0.368   0.368   0.365   0.365   0.366   0.364   0.361   0.361   0.356   0.355   0.357   0.354   0.353   0.350   0.351   0.353   0.355   0.350   0.354   0.352   0.351   0.348   0.348   0.347   0.345   0.346   0.343   0.348   0.346   0.344   0.343   0.342   0.341   0.346   0.346   0.345   0.343   0.348   0.345   0.346   0.342   0.344   0.344   0.340   0.341   0.345   0.345   0.343   0.339   0.343   0.344   0.343   0.346   0.344   0.344   0.345   0.347   0.344   0.344   0.338   0.340   0.343   0.340   0.342   0.336   0.334   0.336   0.337   0.338   0.338   0.343   0.342   0.342   0.336   0.334   0.336   0.330   0.325   0.324   0.323   0.319   0.323   0.322   0.318   0.314   0.314   0.319   0.315   0.316   0.313   0.315   0.314   0.314   0.315   0.313   0.308   0.312   0.311   0.310   0.312   0.311'
                  ' 0.311   0.309   0.309   0.316   0.317   0.312   0.309   0.311   0.308   0.310   0.312'.split('\t'), dtype=float)

旧数据准备：

x=[]
y=[]
x_val = 0
for i in values: #values is a list of several thousand experimental data points
  if x_val < 100:
    x.append(float(x_val))
    y.append(float(i))
  x_val += 5
x = np.asarray(x)
y = np.asarray(y)

更轻松的数据准备：

y1 = values[:100//5]
x1 = np.arange(len(y1))*5

检查是否相同：

print np.allclose(y, y1)
print np.allclose(x, x1)

使用 numpy 定义拟合函数：

def function(x, a,b,c):
    #a = upper asymptote
    #b = negative = x axis displacement
    #c = negative = growth rate
    return a*(np.exp(b*(np.exp(c*x))))

使用起点 p0 拟合：

pars, pcov = opt.curve_fit(function, x1, y1, p0=[0.1, -10, 0.1])

画：

yaj = function(x1, *pars)
plt.figure(1, figsize=(8.5, 11))
plt.plot(x1, y1, 'g-', x1, yaj, 'r-')
plt.xlim(min(x1), max(x1))
plt.ylim(min(y1), max(y1))
plt.show()
于 2014-02-21T08:45:57.730 回答