In [2]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from scipy.optimize import curve_fit
1. Определение концентрации фермента.¶
In [3]:
protein_volume = 1.5e-6 # л
protein_mass_concentration = 0.5 # мг/мл = г/л
protein_molecular_mass = 6.5e4 # Да = г/моль
cell_volume = 1.25e-3 # л
protein_mass = protein_volume * protein_mass_concentration # 7.5е-7 г
protein_amount = protein_mass / protein_molecular_mass # 1.1538е-14 моль
protein_concentration = protein_amount / cell_volume # 9.23е-9 моль/л
In [4]:
protein_concentration
Out[4]:
Data upload¶
In [5]:
data = pd.read_csv("kinetics_data.csv", sep='\t', index_col=None, header=[0,1])
In [6]:
data.head()
Out[6]:
In [7]:
abts_concentrations = [7, 10, 11.5, 15, 24, 33, 50, 100, 200, 500]
Задание 2.¶
2.1. Линейная аппроксимация начальных скоростей.¶
In [8]:
len(abts_concentrations)
Out[8]:
In [9]:
plt.figure(figsize=(10, 25))
for i, conc in enumerate(abts_concentrations):
plt.subplot(5,2,i+1)
plt.plot(data[f'ABTS, {conc} uM']['time'], data[f'ABTS, {conc} uM']['consumption.S'], 'o')
plt.title(f'ABTS, {conc} uM')
plt.xlabel('time, sec')
plt.ylabel('O2 concentration, uM')
plt.tight_layout()
Видно, что при малых концентрациях ABTS нет красивой линейной зависимости. Посмотрим на первые 5 секунд реакции.
In [10]:
plt.figure(figsize=(10, 25))
for i, conc in enumerate(abts_concentrations):
plt.subplot(5,2,i+1)
plt.plot(data[f'ABTS, {conc} uM'].query('time <= 5')['time'],
data[f'ABTS, {conc} uM'].query('time <= 5')['consumption.S'], 'o')
plt.title(f'ABTS, {conc} uM, initial')
plt.xlabel('time, sec')
plt.ylabel('O2 concentration, uM')
plt.tight_layout()
Всё равно выглядит диковато, но пусть так.
In [11]:
o2_speed = list()
for conc in abts_concentrations:
slope, intercept = np.polyfit(x=data[f'ABTS, {conc} uM'].query('time <= 5')['time'],
y=data[f'ABTS, {conc} uM'].query('time <= 5')['consumption.S'],
deg=1)
o2_speed.append([slope, intercept])
o2_speed = np.array(o2_speed, dtype=float)
In [12]:
o2_speed
Out[12]:
In [13]:
plt.plot(abts_concentrations, o2_speed[:, 0])
plt.xlabel("ABTS concentration, uM")
plt.ylabel("Speed")
plt.show()
Мы видим некоторое расхождение в известной нам зависимости, связано, по-видимому, с особенностями аппроксимации, а именно с большой дисперсией поглощения кислорода при концентрации ABTS 7 uM, 11.5 uM, 15 uM. Но в целом, очертания похожи на правду.
In [14]:
initial_speed = pd.DataFrame({'concentration': abts_concentrations, 'speed': o2_speed[:, 0] * 4}) # домножили на 4, так как ABTS в 4 раза больше расходуется, чем выделяется О2 по уравнению реакции.
In [15]:
initial_speed
Out[15]:
Теперь посмотрим, как зафитились наши кривые в начале.
In [16]:
plt.figure(figsize=(10, 25))
for i, conc in enumerate(abts_concentrations):
plt.subplot(5,2,i+1)
plt.plot(data[f'ABTS, {conc} uM'].query('time <= 5')['time'],
data[f'ABTS, {conc} uM'].query('time <= 5')['consumption.S'],
'o', label='data')
plt.plot(data[f'ABTS, {conc} uM'].query('time <= 5')['time'],
data[f'ABTS, {conc} uM'].query('time <= 5')['time'] * o2_speed[i, 0] + o2_speed[i, 1],
color='red', label='fitted line')
plt.legend()
plt.title(f'ABTS, {conc} uM, initial')
plt.xlabel('time, sec')
plt.ylabel('O2 concentration, uM')
plt.tight_layout()
Чем больше концентрация ABTS, тем более логичной выглядит аппроксимация.
2.2. Аппроксимация экспонентой.¶
In [17]:
def expassoc(x, y0, a1, t1, a2, t2):
return y0 + a1 * (1 - np.exp(-x / t1)) + a2 * (1 - np.exp(-x / t2))
По совету преподавателя введём ограничение для данных — первые 60 секунд.
In [23]:
params_list = list()
params_names = ['y0', 'a1', 't1', 'a2', 't2']
for conc in abts_concentrations:
try:
exp_ind = data[f'ABTS, {conc} uM'].dropna()['time'] < 60
params, pcov = curve_fit(expassoc,
data[f'ABTS, {conc} uM'].dropna()['time'][exp_ind],
data[f'ABTS, {conc} uM'].dropna()['consumption.S'][exp_ind],
bounds=([-np.inf, -np.inf, 0, -np.inf, 0], np.inf))
except RuntimeError:
params = [0 for _ in range(5)]
print(f"{conc} didn't fit params")
params_list.append(dict(zip(params_names, params)))
Удивительно, но 200 и 500 не зафитились, несмотря на то, что установил ограничения на параметры, как было написано в документации. Возможно, проблема в том, что первоначальные значения все равны 1 и Trust Region Reflective algorithm не смог подобрать правильные параметры. В Origin используют алгоритм Levenberg-Marquardt, но в scipy он не умеет работать, если заданы ограничения. Тогда можно попробовать dogleg algorithm или поменять первоначальные значения.
In [90]:
pd.DataFrame(params_list)
Out[90]:
In [87]:
params_list = list()
params_names = ['y0', 'a1', 't1', 'a2', 't2']
for conc in abts_concentrations:
try:
params, pcov = curve_fit(expassoc,
data[f'ABTS, {conc} uM'].dropna()['time'],
data[f'ABTS, {conc} uM'].dropna()['consumption.S'],
bounds=([-np.inf, -np.inf, 0, -np.inf, 0], np.inf),
method='dogbox')
except RuntimeError:
params = [0 for _ in range(5)]
print(f"{conc} didn't fit params")
params_list.append(dict(zip(params_names, params)))
Ха, dogbox вообще упал. Тогда возвращается к дефолту, меняем первоначальные параметры. Допустим, y0 = -1, t1 = t2 = 10, A1 = A2 = 5 исходя из предыдущего результата.
In [92]:
params_list = list()
params_names = ['y0', 'a1', 't1', 'a2', 't2']
for conc in abts_concentrations:
try:
params, pcov = curve_fit(expassoc,
data[f'ABTS, {conc} uM'].dropna()['time'],
data[f'ABTS, {conc} uM'].dropna()['consumption.S'],
bounds=([-np.inf, -np.inf, 0, -np.inf, 0], np.inf),
p0=[-1, 5, 10, 5, 10])
except RuntimeError:
params = [0 for _ in range(5)]
print(f"{conc} didn't fit params")
params_list.append(dict(zip(params_names, params)))
Ну, мы зафитили 200, но зато не смогли 50. Что же это такое :(
In [93]:
pd.DataFrame(params_list)
Out[93]:
Давайте попробуем с дефолтными первоначальными значениями увеличить число вычислений функции. По умолчанию, это 100 n, где n -- количество точек. Сделаем это 1000 n.
In [160]:
params_list = list()
params_names = ['y0', 'a1', 't1', 'a2', 't2']
for conc in abts_concentrations:
try:
params, pcov = curve_fit(expassoc,
data[f'ABTS, {conc} uM'].dropna()['time'],
data[f'ABTS, {conc} uM'].dropna()['consumption.S'],
bounds=([-np.inf, -np.inf, 0, -np.inf, 0], np.inf),
max_nfev=1000 * data[f'ABTS, {conc} uM'].dropna()['time'].shape[0])
except RuntimeError:
params = [0 for _ in range(5)]
print(f"{conc} didn't fit params")
params_list.append(dict(zip(params_names, params)))
In [161]:
exp_fit = pd.DataFrame(params_list)
Ура, всё зафитилось, давайте нарисуем.
In [162]:
plt.figure(figsize=(10, 25))
for i, conc in enumerate(abts_concentrations):
plt.subplot(5,2,i+1)
plt.plot(data[f'ABTS, {conc} uM']['time'], data[f'ABTS, {conc} uM']['consumption.S'], 'o', label='data')
plt.plot(data[f'ABTS, {conc} uM']['time'],
expassoc(data[f'ABTS, {conc} uM']['time'], **exp_fit.loc[i, :]), color='red', label='fitted line')
plt.legend()
plt.title(f'ABTS, {conc} uM')
plt.xlabel('time, sec')
plt.ylabel('O2 concentration, uM')
plt.tight_layout()
In [163]:
exp_fit
Out[163]:
Пока дисперсия точек большая, кривая фитится весьма красиво. Но при больших концентрациях становится видно, что сумма двух экспонент не очень-то правильно описывает эксперимент.
In [164]:
der_speed = pd.DataFrame({'concentration': abts_concentrations, 'speed': (exp_fit['a1'] / exp_fit['t1'] + exp_fit['a2'] / exp_fit['t2']) * 4})
In [165]:
plt.plot(der_speed['concentration'], der_speed['speed'])
Out[165]:
Зависимость вообще получается какая-то не очень красивая. По-видимому, это связано с тем, что сумма экспонент не очень хорошо фитирует кривые. А ещё при концентрации в 15 uM скорость отрицательная. Это может потом вредно сказать на фитинге данных по Михаэлису-Ментен.
Может, удастся лучше зафитить, если в качестве первоначального приближения задать какие-то числа?
In [166]:
params_list = list()
params_names = ['y0', 'a1', 't1', 'a2', 't2']
for conc in abts_concentrations:
try:
params, pcov = curve_fit(expassoc,
data[f'ABTS, {conc} uM'].dropna()['time'],
data[f'ABTS, {conc} uM'].dropna()['consumption.S'],
bounds=([-np.inf, -np.inf, 0, -np.inf, 0], np.inf),
max_nfev=1000 * data[f'ABTS, {conc} uM'].dropna()['time'].shape[0],
p0=[-1, 5, 10, 5, 10])
except RuntimeError:
params = [0 for _ in range(5)]
print(f"{conc} didn't fit params")
params_list.append(dict(zip(params_names, params)))
In [167]:
exp_fit = pd.DataFrame(params_list)
In [168]:
plt.figure(figsize=(10, 25))
for i, conc in enumerate(abts_concentrations):
plt.subplot(5,2,i+1)
plt.plot(data[f'ABTS, {conc} uM']['time'], data[f'ABTS, {conc} uM']['consumption.S'], 'o', label='data')
plt.plot(data[f'ABTS, {conc} uM']['time'],
expassoc(data[f'ABTS, {conc} uM']['time'], **exp_fit.loc[i, :]), color='red', label='fitted line')
plt.legend()
plt.title(f'ABTS, {conc} uM')
plt.xlabel('time, sec')
plt.ylabel('O2 concentration, uM')
plt.tight_layout()
На первый взгляд, зафитил так же.
In [169]:
der_speed = pd.DataFrame({'concentration': abts_concentrations,
'speed': (exp_fit['a1'] / exp_fit['t1'] + exp_fit['a2'] / exp_fit['t2']) * 4})
In [202]:
plt.plot(der_speed['concentration'], der_speed['speed'])
plt.xlabel("ABTS concentration, uM")
plt.ylabel("Speed")
plt.show()
Но на деле нет отрицательной скорости.
In [173]:
der_speed
Out[173]:
In [174]:
plt.plot(initial_speed['speed'], der_speed['speed'], 'o')
plt.xlabel("Speed from initial speed approximation")
plt.ylabel("Speed from exponential approximation")
plt.show()
Скорости несильно похожи друг на друга. Наверняка это связано с плохим качеством одной из аппроксимаций. Почему-то кажется, что проблемы именно в начальных скоростях.
Посмотрим на оценки MSE для каждого из методов для каждой концентрации, чтобы были основания для оценки алгоритмов.
In [175]:
mse_asses = list()
for i, conc in enumerate(abts_concentrations):
mse_linear = ((data[f'ABTS, {conc} uM'].query('time <= 5')['time'] * o2_speed[i, 0] + o2_speed[i, 1] - data[f'ABTS, {conc} uM'].query('time <= 5')['consumption.S']) ** 2).mean()
mse_expassoc = ((expassoc(data[f'ABTS, {conc} uM']['time'], **exp_fit.loc[i, :]) - data[f'ABTS, {conc} uM']['consumption.S']) ** 2).mean()
mse_expassoc_init = ((expassoc(data[f'ABTS, {conc} uM'].query('time <=5')['time'], **exp_fit.loc[i, :]) - data[f'ABTS, {conc} uM'].query('time <=5')['consumption.S']) ** 2).mean()
mse_asses.append({'concentration': conc, 'initial speed': mse_linear, 'expassoc': mse_expassoc, 'expassoc_init': mse_expassoc_init})
mse_pd = pd.DataFrame(mse_asses)
In [176]:
plt.plot(mse_pd['concentration'], mse_pd['initial speed'], 'o', label='initial speed', alpha=.5)
plt.plot(mse_pd['concentration'], mse_pd['expassoc'], 'o', label='expassoc all', color='red', alpha=.5)
plt.plot(mse_pd['concentration'], mse_pd['expassoc_init'], 'o', label='expassoc initial', color='green', alpha=.5)
plt.legend(loc='upper right')
plt.xlabel('concentration')
plt.ylabel('MSE')
plt.show()
При малых концентрациях ABTS MSE обоих методов примерно одинакова, но затем линейная аппроксимация остаётся стабильно хорошей, а вот expassoc уходит вверх. Это связано с тем, что expassoc не смогла повторить контуры кривых при высоких концентрациях при большом времени. Но здесь у нас была линейная аппроксимация начальных скоростей до 5 секунд, а expassoc фитировала весь датасет. Если посмотреть на expassoc на первых 5 секундах, то ситуация становится хуже.
2.* Пересчитываем expassoc¶
По совету преподавателя, пересчитываем expassoc по первым 60 секундам.
In [22]:
params_list = list()
params_names = ['y0', 'a1', 't1', 'a2', 't2']
for conc in abts_concentrations:
try:
exp_ind = data[f'ABTS, {conc} uM'].dropna()['time'] < 60
params, pcov = curve_fit(expassoc,
data[f'ABTS, {conc} uM'].dropna()['time'][exp_ind],
data[f'ABTS, {conc} uM'].dropna()['consumption.S'][exp_ind],
bounds=([-np.inf, -np.inf, 0, -np.inf, 0], np.inf))
except RuntimeError:
params = [0 for _ in range(5)]
print(f"{conc} didn't fit params")
params_list.append(dict(zip(params_names, params)))
Не зафитились 24, 33, 100 мкМ, что очень грустно. Посмотрим на последний успешный фиттинг по всему времени (с увеличенным количеством итераций и начальными параметрами).
In [30]:
params_list = list()
params_names = ['y0', 'a1', 't1', 'a2', 't2']
for conc in abts_concentrations:
try:
exp_ind = data[f'ABTS, {conc} uM'].dropna()['time'] < 60
params, pcov = curve_fit(expassoc,
data[f'ABTS, {conc} uM'].dropna()['time'][exp_ind],
data[f'ABTS, {conc} uM'].dropna()['consumption.S'][exp_ind],
bounds=([-np.inf, -np.inf, 0, -np.inf, 0], np.inf),
max_nfev=1000 * np.sum(exp_ind),
p0=[-1, 5, 10, 5, 10])
except RuntimeError:
params = [0 for _ in range(5)]
print(f"{conc} didn't fit params")
params_list.append(dict(zip(params_names, params)))
Ну тут всё зафитилось. Посмотрим, как.
In [31]:
exp_fit = pd.DataFrame(params_list)
plt.figure(figsize=(10, 25))
for i, conc in enumerate(abts_concentrations):
plt.subplot(5,2,i+1)
plt.plot(data[f'ABTS, {conc} uM']['time'], data[f'ABTS, {conc} uM']['consumption.S'], 'o', label='data')
plt.plot(data[f'ABTS, {conc} uM']['time'],
expassoc(data[f'ABTS, {conc} uM']['time'], **exp_fit.loc[i, :]), color='red', label='fitted line')
plt.legend()
plt.title(f'ABTS, {conc} uM')
plt.xlabel('time, sec')
plt.ylabel('O2 concentration, uM')
plt.tight_layout()
В целом, начальные участки описаны хорошо.
In [32]:
exp_fit
Out[32]:
Но параметры неоднородны.
In [33]:
der_speed = pd.DataFrame({'concentration': abts_concentrations,
'speed': (exp_fit['a1'] / exp_fit['t1'] + exp_fit['a2'] / exp_fit['t2']) * 4})
plt.plot(der_speed['concentration'], der_speed['speed'])
plt.xlabel("ABTS concentration, uM")
plt.ylabel("Speed")
plt.show()
In [38]:
der_speed
Out[38]:
А тут вообще дичь какая-то со скоростями.
In [34]:
params_list = list()
params_names = ['y0', 'a1', 't1', 'a2', 't2']
for conc in abts_concentrations:
try:
exp_ind = data[f'ABTS, {conc} uM'].dropna()['time'] < 60
params, pcov = curve_fit(expassoc,
data[f'ABTS, {conc} uM'].dropna()['time'][exp_ind],
data[f'ABTS, {conc} uM'].dropna()['consumption.S'][exp_ind],
bounds=([-np.inf, -np.inf, 0, -np.inf, 0], np.inf),
max_nfev=1000 * np.sum(exp_ind))
except RuntimeError:
params = [0 for _ in range(5)]
print(f"{conc} didn't fit params")
params_list.append(dict(zip(params_names, params)))
Ну тут всё зафитилось. Посмотрим, как.
In [36]:
exp_fit = pd.DataFrame(params_list)
plt.figure(figsize=(10, 25))
for i, conc in enumerate(abts_concentrations):
plt.subplot(5,2,i+1)
plt.plot(data[f'ABTS, {conc} uM']['time'], data[f'ABTS, {conc} uM']['consumption.S'], 'o', label='data')
plt.plot(data[f'ABTS, {conc} uM']['time'],
expassoc(data[f'ABTS, {conc} uM']['time'], **exp_fit.loc[i, :]), color='red', label='fitted line')
plt.legend()
plt.title(f'ABTS, {conc} uM')
plt.xlabel('time, sec')
plt.ylabel('O2 concentration, uM')
plt.tight_layout()
In [37]:
der_speed = pd.DataFrame({'concentration': abts_concentrations,
'speed': (exp_fit['a1'] / exp_fit['t1'] + exp_fit['a2'] / exp_fit['t2']) * 4})
plt.plot(der_speed['concentration'], der_speed['speed'])
plt.xlabel("ABTS concentration, uM")
plt.ylabel("Speed")
plt.show()
In [39]:
der_speed
Out[39]:
Какая-то дичь со скоростями. Из графиков видно, что для больших концентраций правый хвост плохо аппроксимируется. Может, стоит увеличить время подсчёта так, чтобы захватывали линейный кусок.
In [60]:
right_bounds = {7: 60,
10: 60,
11.5: 60,
15: 70,
24: 70,
33: 80,
50: 80,
100: 80,
200: 150,
500: 250
}
In [61]:
params_list = list()
params_names = ['y0', 'a1', 't1', 'a2', 't2']
for conc in abts_concentrations:
try:
exp_ind = data[f'ABTS, {conc} uM'].dropna()['time'] < right_bounds[conc]
params, pcov = curve_fit(expassoc,
data[f'ABTS, {conc} uM'].dropna()['time'][exp_ind],
data[f'ABTS, {conc} uM'].dropna()['consumption.S'][exp_ind],
bounds=([-np.inf, -np.inf, 0, -np.inf, 0], np.inf),
max_nfev=1000 * np.sum(exp_ind))
except RuntimeError:
params = [0 for _ in range(5)]
print(f"{conc} didn't fit params")
params_list.append(dict(zip(params_names, params)))
In [62]:
exp_fit = pd.DataFrame(params_list)
plt.figure(figsize=(10, 25))
for i, conc in enumerate(abts_concentrations):
plt.subplot(5,2,i+1)
plt.plot(data[f'ABTS, {conc} uM']['time'], data[f'ABTS, {conc} uM']['consumption.S'], 'o', label='data')
plt.plot(data[f'ABTS, {conc} uM']['time'],
expassoc(data[f'ABTS, {conc} uM']['time'], **exp_fit.loc[i, :]), color='red', label='fitted line')
plt.legend()
plt.title(f'ABTS, {conc} uM')
plt.xlabel('time, sec')
plt.ylabel('O2 concentration, uM')
plt.tight_layout()
In [63]:
der_speed = pd.DataFrame({'concentration': abts_concentrations,
'speed': (exp_fit['a1'] / exp_fit['t1'] + exp_fit['a2'] / exp_fit['t2']) * 4})
plt.plot(der_speed['concentration'], der_speed['speed'])
plt.xlabel("ABTS concentration, uM")
plt.ylabel("Speed")
plt.show()
100 мкМ нормально не зафитировалась, но хотя бы скорости более-менее приличные теперь.
Задание 3. Вычисление Km и Vm.¶
Согласно модели Михаэлиса-Ментен Vm = k2 * E0, поэтому мы можем посчитать заодно и k2.
3.1. Прямые координаты.¶
In [64]:
def hyperbl(x, p1, p2):
return p1 * x / (p2 + x)
In [65]:
hyperbl_init, cov = curve_fit(hyperbl, initial_speed['concentration'], initial_speed['speed'])
hyperbl_exp, cov = curve_fit(hyperbl, der_speed['concentration'], der_speed['speed'])
In [83]:
print("Hyperbolic fitting V(S)")
print("Type\tKm\tVm\tk2")
print(f"Init\t{hyperbl_init[1]:.4f}\t{hyperbl_init[0]:.4f}\t{hyperbl_init[0]/protein_concentration / 1e6:.4f}")
print(f"Exp\t{hyperbl_exp[1]:.4f}\t{hyperbl_exp[0]:.4f}\t{hyperbl_exp[0]/protein_concentration / 1e6:.4f}")
In [67]:
plt.figure(figsize=(10, 5))
plt.subplot(1,2,1)
plt.plot(initial_speed['concentration'], initial_speed['speed'], 'o', label='data')
plt.plot(initial_speed['concentration'], hyperbl(initial_speed['concentration'], *hyperbl_init), label='fit', color='red')
plt.xlabel('S')
plt.ylabel('V')
plt.legend(loc='lower right')
plt.title('Hyperbolic fitting of initial speed data')
plt.subplot(1,2,2)
plt.plot(der_speed['concentration'], der_speed['speed'], 'o', label='data')
plt.plot(der_speed['concentration'], hyperbl(der_speed['concentration'], *hyperbl_exp), label='fit', color='red')
plt.xlabel('S')
plt.ylabel('V')
plt.legend(loc='lower right')
plt.title('Hyperbolic fitting of expassoc data')
plt.tight_layout()
3.2. Координаты Хейнса.¶
In [78]:
B_init, A_init = np.polyfit(initial_speed['concentration'], initial_speed['concentration'] / initial_speed['speed'], deg=1)
B_exp, A_exp = np.polyfit(der_speed['concentration'], der_speed['concentration'] / der_speed['speed'], deg=1)
print("Hance fitting V/S(S)")
print("Type\tKm\tVm\tk2")
print(f"Init\t{A_init / B_init:.4f}\t{1 / B_init:.4f}\t{1 / B_init / protein_concentration / 1e6:.4f}")
print(f"Exp\t{A_exp / B_exp:.4f}\t{1 / B_exp:.4f}\t{1 / B_exp / protein_concentration / 1e6:.4f}")
In [69]:
plt.figure(figsize=(10, 5))
plt.subplot(1,2,1)
plt.plot(initial_speed['concentration'], initial_speed['concentration'] / initial_speed['speed'], 'o', label='data')
plt.plot(initial_speed['concentration'], A_init + B_init * initial_speed['concentration'], label='fit', color='red')
plt.xlabel('S')
plt.ylabel('S/V')
plt.legend(loc='lower right')
plt.title('Hance fitting of initial speed data')
plt.subplot(1,2,2)
plt.plot(der_speed['concentration'], der_speed['concentration'] / der_speed['speed'], 'o', label='data')
plt.plot(der_speed['concentration'], A_exp + B_exp * der_speed['concentration'], label='fit', color='red')
plt.xlabel('S')
plt.ylabel('S/V')
plt.legend(loc='lower right')
plt.title('Hance fitting of expassoc data')
plt.tight_layout()
Для expassoc получили отрицательную Km, что в природе быть не может. Это связано с тем, скорее всего, что мы видим выброс для 500 uM. Надо его убрать и снова зафитить.
In [79]:
B_init, A_init = np.polyfit(initial_speed['concentration'], initial_speed['concentration'] / initial_speed['speed'], deg=1)
B_exp, A_exp = np.polyfit(der_speed['concentration'][:-1], der_speed['concentration'][:-1] / der_speed['speed'][:-1], deg=1)
print("Hance fitting V/S(S)")
print("Type\tKm\tVm\tk2")
print(f"Init\t{A_init / B_init:.4f}\t{1 / B_init:.4f}\t{1 / B_init / protein_concentration / 1e6:.4f}")
print(f"Exp\t{A_exp / B_exp:.4f}\t{1 / B_exp:.4f}\t{1 / B_exp / protein_concentration / 1e6:.4f}")
Уже данные неотрицательные, но значения сильно отличаются от начальных скоростей. Всё же возьмём их, чтобы показать проблематичность экспоненциальных данных.
In [71]:
plt.figure(figsize=(10, 5))
plt.subplot(1,2,1)
plt.plot(initial_speed['concentration'], initial_speed['concentration'] / initial_speed['speed'], 'o', label='data')
plt.plot(initial_speed['concentration'], A_init + B_init * initial_speed['concentration'], label='fit', color='red')
plt.xlabel('S')
plt.ylabel('S/V')
plt.legend(loc='lower right')
plt.title('Hance fitting of initial speed data')
plt.subplot(1,2,2)
plt.plot(der_speed['concentration'], der_speed['concentration'] / der_speed['speed'], 'o', label='data')
plt.plot(der_speed['concentration'], A_exp + B_exp * der_speed['concentration'], label='fit', color='red')
plt.xlabel('S')
plt.ylabel('S/V')
plt.legend(loc='lower right')
plt.title('Hance fitting of expassoc data')
plt.tight_layout()
По-видимому, проблема всё же в 100 и 200 uM.
3.3. Координаты Иди-Хоффсти.¶
In [80]:
B_init, A_init = np.polyfit(initial_speed['speed'] / initial_speed['concentration'], initial_speed['speed'], deg=1)
B_exp, A_exp = np.polyfit(der_speed['speed'] / der_speed['concentration'], der_speed['speed'], deg=1)
print("IF fitting V(V/S)")
print("Type\tKm\tVm\tk2")
print(f"Init\t{-B_init:.4f}\t{A_init:.4f}\t{A_init / protein_concentration / 1e6:.4f}")
print(f"Exp\t{-B_exp:.4f}\t{A_exp:.4f}\t{A_exp / protein_concentration / 1e6:.4f}")
In [73]:
plt.figure(figsize=(10, 5))
plt.subplot(1,2,1)
plt.plot(initial_speed['speed'] / initial_speed['concentration'], initial_speed['speed'], 'o', label='data')
plt.plot(initial_speed['speed'] / initial_speed['concentration'],
A_init + B_init * initial_speed['speed'] / initial_speed['concentration'], label='fit', color='red')
plt.xlabel('V / S')
plt.ylabel('V')
plt.legend(loc='lower left')
plt.title('IF fitting of initial speed data')
plt.subplot(1,2,2)
plt.plot(der_speed['speed'] / der_speed['concentration'], der_speed['speed'], 'o', label='data')
plt.plot(der_speed['speed'] / der_speed['concentration'],
A_exp + B_exp * der_speed['speed'] / der_speed['concentration'], label='fit', color='red')
plt.xlabel('V / S')
plt.ylabel('V')
plt.legend(loc='lower left')
plt.title('IF fitting of expassoc data')
plt.tight_layout()
Данные начальных скоростей выглядят в этих координатах более качественно линейно фитируемыми.
3.4. Координаты Лайнуивера-Берка.¶
In [82]:
B_init, A_init = np.polyfit(1 / initial_speed['concentration'], 1 / initial_speed['speed'], deg=1)
B_exp, A_exp = np.polyfit(1 / der_speed['concentration'], 1 / der_speed['speed'], deg=1)
print("IF fitting V(V/S)")
print("Type\tKm\tVm\tk2")
print(f"Init\t{B_init / A_init:.4f}\t{1 / A_init:.4f}\t{1 / A_init / protein_concentration / 1e6:.4f}")
print(f"Exp\t{B_exp / A_exp:.4f}\t{1 / A_exp:.4f}\t{1 / A_exp / protein_concentration / 1e6:.4f}")
In [75]:
plt.figure(figsize=(10, 5))
plt.subplot(1,2,1)
plt.plot(1 / initial_speed['concentration'], 1 / initial_speed['speed'], 'o', label='data')
plt.plot(1 / initial_speed['concentration'],
A_init + B_init / initial_speed['concentration'], label='fit', color='red')
plt.xlabel('1 / S')
plt.ylabel('1 / V')
plt.legend(loc='lower right')
plt.title('LB fitting of initial speed data')
plt.subplot(1,2,2)
plt.plot(1 / der_speed['concentration'], 1 / der_speed['speed'], 'o', label='data')
plt.plot(1 / der_speed['concentration'],
A_exp + B_exp / der_speed['concentration'], label='fit', color='red')
plt.xlabel('1 / S')
plt.ylabel('1 / V')
plt.legend(loc='lower right')
plt.title('LB fitting of expassoc data')
plt.tight_layout()
И вновь данные начальных скоростей оказались более аккуратными.
Задание 5. Сравнение кинетических параметров.¶
In [87]:
kinet_params = pd.read_csv('kinetics_params.tsv', sep='\t', index_col=None, header=[0])
In [88]:
kinet_melt = pd.melt(kinet_params, id_vars=['Coordinates', 'Data source'], value_vars=['Km', 'Vm', 'k2'])
In [89]:
plt.figure(figsize=(15,5))
plt.subplot(1,3,1)
sns.scatterplot(data=kinet_params, x='Coordinates', y='Km', hue='Data source')
plt.title('Km')
plt.subplot(1,3,2)
sns.scatterplot(data=kinet_params, x='Coordinates', y='Vm', hue='Data source')
plt.title('Vm')
plt.subplot(1,3,3)
sns.scatterplot(data=kinet_params, x='Coordinates', y='k2', hue='Data source')
plt.title('k2')
plt.tight_layout()
Видно, что кинетические параметры из данных начальных скоростей согласованы друг с другом, а экспоненты -- нет.
Доп. Фитируем первоначальные данные с гиперболой.¶
In [90]:
params_list = list()
params_names = ['p1', 'p2']
for conc in abts_concentrations:
try:
params, pcov = curve_fit(hyperbl,
data[f'ABTS, {conc} uM'].dropna()['time'],
data[f'ABTS, {conc} uM'].dropna()['consumption.S'],
bounds=(0, np.inf))
except RuntimeError:
params = [0 for _ in range(5)]
print(f"{conc} didn't fit params")
params_list.append(dict(zip(params_names, params)))
In [91]:
hyper_fit = pd.DataFrame(params_list)
In [92]:
plt.figure(figsize=(10, 25))
for i, conc in enumerate(abts_concentrations):
plt.subplot(5,2,i+1)
plt.plot(data[f'ABTS, {conc} uM']['time'], data[f'ABTS, {conc} uM']['consumption.S'], 'o', label='data')
plt.plot(data[f'ABTS, {conc} uM']['time'],
hyperbl(data[f'ABTS, {conc} uM']['time'], **hyper_fit.loc[i, :]), color='red', label='fitted line')
plt.legend()
plt.title(f'ABTS, {conc} uM')
plt.xlabel('time, sec')
plt.ylabel('O2 concentration, uM')
plt.tight_layout()
Зафитировала хуже, чем expassoc.
In [93]:
hyper_speed = pd.DataFrame({'concentration': abts_concentrations, 'speed': hyper_fit['p1'] / hyper_fit['p2'] * 4})
In [94]:
plt.plot(hyper_speed['concentration'], hyper_speed['speed'])
plt.xlabel("ABTS, uM")
plt.ylabel("Speed")
plt.show()
Сравним все 3 метода.
In [95]:
plt.plot(initial_speed['concentration'], initial_speed['speed'], label='Initial')
plt.plot(der_speed['concentration'], der_speed['speed'], color='red', label='Expassoc')
plt.plot(hyper_speed['concentration'], hyper_speed['speed'], color='green', label='Hyperbl')
plt.xlabel("ABTS, uM")
plt.ylabel("Speed")
plt.legend()
plt.show()
Гипербола ещё хуже, чем expassoc, поэтому дальше и не будем считать из неё кинетические параметры.