add sigmoid, hill

droplet_growth/fit.py

@@ -228,15 +228,20 @@ def lin_exp_fit(curve:list, plot=False, fun=lin_exp_fun, **kwargs):
     return popt
 
 
-def fit_exp_on_baseline(data, x=None, p0=(1,1,0), **kwargs):
+def fit_exp_on_baseline(data, x=None, p0=(1,1,0), plot=False, **kwargs):
     """fits a + b * exp(c * x)
-    kwargs: same as in fit.fit_exponent
+    kwargs: same as in fit.fit_exponent and fit.plot_fit
     returns:
     (a, b, c)
     """
     try:
-        return fit_exponent(data, bins=x, fun=exp_on_baseline,
-                            p0=p0, **kwargs)
+        popt, bins = fit_exponent(data, bins=x, fun=exp_on_baseline,
+                            p0=p0, plot=False, return_bins=True)
+        if plot:
+            curve = data - popt[0]
+            fit = exp_on_baseline(bins, 0, popt[1], popt[2])
+            plot_fit(curve, fit, bins, plot=plot, **kwargs)
+        return (popt)
     except RuntimeError:
         print('No fit')
         return (None, None, None)
@@ -247,15 +252,36 @@ def lag_exponent(x, lag, c):
     return np.exp((x - lag) / c)
 
 
-def sigmoid(t, a, b, c):
+def gompertz(t, a, b, c):
     '''
     https://en.wikipedia.org/wiki/Gompertz_function
     '''
     return a * np.exp(- b * np.exp(- c * t))
 
+def sigmoid(x, a, b):
+    return 1 / (1 + np.exp(a * x + b))
+
+
+def hill(x, n, K):
+    return x ** n/(x ** n + K)
+
+
+def fit_sigmoid(probs, ax, fun=sigmoid, fit_name='sigmoid', p0=(2., -2.)):
+    probs = probs.copy()
+    vector = probs.index
+#     print(vector)
+#     probs.loc[:, 'negative'] = 1 - probs.positive
+    popt, pcov = curve_fit(fun, vector, probs.q, p0=p0)
+    a, b = popt
+    da, db = np.sqrt(np.diag(pcov))
+    ax.plot(vector, probs.q, '.', label='data')
+    ax.plot((x := sorted(vector)), fun(x, *popt), lw=10, alpha=.5, label=f'{fit_name} fit')
+    plt.legend()
+    return popt, (da, db)
+
 
 def fit_exponent(curve, bins=None, fun=exponent, p0=(1., .5,), bounds=(-np.inf, np.inf),
-                 plot=False, plot_init=False, **kwargs):
+                 plot=False, plot_init=False, return_bins=False, **kwargs):
     '''
     Fits exponent to 1D curve
     plot: ['linear', 'log', None]
@@ -263,14 +289,18 @@ def fit_exponent(curve, bins=None, fun=exponent, p0=(1., .5,), bounds=(-np.inf,
     if bins is None:
         bins = np.arange(len(curve))
     popt, _ = curve_fit(f=fun, xdata=bins, ydata=curve, p0=p0, bounds=bounds)
+    
     fit_result = fun(bins, *popt)
     plot_fit(curve, fit_result, bins, plot=plot, **kwargs)
     if plot_init:
         plot_fit(curve, fun(bins, *p0), bins, plot=plot, labels=['init', 'data'], **kwargs)
+    if return_bins:
+        return popt, bins
     return popt
 
 
-def plot_fit(curve=None, fit=None, vector=None, plot='linear', labels=['data', 'fit'], markers=['.', 'o'], legend=True, **kwargs):
+def plot_fit(curve=None, fit=None, vector=None, plot='linear', 
+    labels=['data', 'fit'], markers=['.', '-'], legend=True, **kwargs):
     '''
     plot: ['linear', 'log', None]
     '''
@@ -296,11 +326,11 @@ def add_doubling_time(
     rate_column:str='c',
     frame_rate:float=1.,
     new_column='Doubling time, min',
-    fun = lambda c, rate_value: 1 / c * rate_value
+    convert_rate_time = lambda c, frame_rate: np.log(2) / c * frame_rate
 ):
     '''
     adds a column with doubling time using rate 'c'
     '''
     table = df.copy()
-    table[new_column] = fun(table[rate_column].values, frame_rate) 
+    table[new_column] = convert_rate_time(table[rate_column].values, frame_rate) 
     return table
\ No newline at end of file