| |
| import matplotlib.pyplot as plt |
| import matplotlib as mpl |
| mpl.rcParams['font.size'] = 20 |
| import numpy as np |
| import pandas as pd |
| from scipy.integrate import odeint |
| from scipy.optimize import curve_fit |
| from scipy.signal import find_peaks, peak_widths |
| def ave(data): |
| return data.mean(axis=1) |
| def std(data): |
| return data.std(axis=1) |
| def sem(data): |
| return data.std(axis=1)/np.sqrt(data.shape[-1]) |
| def calc_stats(sample_fn, blank_fn, dt=10.0/60.0, convfac=1): |
| blank = np.loadtxt(blank_fn) |
| sample = np.loadtxt(sample_fn) |
| if blank.shape[0] == sample.shape[0]: |
| t = np.arange(0, blank.shape[0])*dt |
| sample_ave = (ave(sample) - ave(blank))*convfac |
| sample_err = np.sqrt(sem(sample)**2+sem(blank)**2) |
| return (t, sample_ave, sample_err) |
| else: |
| print('Error: blank and sample do not have the same number of data/time points') |
| return None |
| def remove_baseline_from_data(sample_fn, blank_fn): |
| blank = np.loadtxt(blank_fn) |
| sample = np.loadtxt(sample_fn) |
| new_sample_data = np.zeros_like(sample) |
| baseline = ave(blank) |
| if blank.shape[0] == sample.shape[0]: |
| for col in range(sample.shape[-1]): |
| new_sample_data[:,col] = sample[:,col] - baseline |
| return new_sample_data |
| else: |
| print('Error: blank and sample do not have the same number of data/time points') |
| return None |
| def smooth(x,window_len=11,window='hanning'): |
| if window_len<3: |
| return x |
| ss=np.r_[x[window_len-1:0:-1],x,x[-2:-window_len-1:-1]] |
| if window == 'flat': |
| w=np.ones(window_len,'d') |
| else: |
| w=eval('np.'+window+'(window_len)') |
| y=np.convolve(w/w.sum(),ss,mode='valid') |
| return y |
| def calc_deriv(t, ys, deriv_order=1, sw=11): |
| dt = t[1]-t[0] |
| dys = [] |
| for y in ys: |
| if deriv_order == 1: |
| dys.append(np.diff(smooth(y, sw))/dt) |
| elif deriv_order==2: |
| dys.append(np.diff(smooth(np.diff(smooth(y, sw))/dt, sw))/dt) |
| else: |
| print('Higher order derivatives are not supported') |
| dys.append(np.diff(smooth(y, sw))/dt) |
| dys = np.array(dys) |
| return (t, dys.mean(axis=0)[:len(t)], dys.std(axis=0)[:len(t)]/np.sqrt(len(ys))) |
| def get_peak_info(t, y): |
| res = [] |
| peaks, _ = find_peaks(y, height=np.max(y)/2.0) |
| for p in peaks: |
| res.append([t[p], y[p]]) |
| return np.array(res) |
| plt.rcParams['font.sans-serif'] = ['SimHei'] |
| plt.rcParams['axes.unicode_minus'] = False |
| |
| |
| plt.figure(figsize=(5,4)) |
| t1, s1, e1 = calc_stats('s1_od.txt', 'blank_od.txt') |
| plt.errorbar(t1, s1, e1, fmt='o-',capsize=0.4, ms=1.2, mew=0.2, |
| elinewidth=0.2, color='tab:blue', mfc='none', label='大肠杆菌') |
| plt.xlabel('时间 (hr)') |
| plt.ylabel('光密度 (a.u.)') |
| plt.legend() |
| plt.tight_layout() |
| plt.savefig('./od值.pdf') |
| |
| t,s,e = calc_stats('s1_od.txt', 'blank_od.txt') |
| sample = remove_baseline_from_data('s1_od.txt', 'blank_od.txt') |
| ys = [sample[:,col] for col in range(sample.shape[-1])] |
| t1, d1, e1 = calc_deriv(t, ys, deriv_order=1, sw=11) |
| p1 = get_peak_info(t1, d1)[0] |
| scfac = 1 |
| plt.figure(figsize=(5,4)) |
| plt.errorbar(t1, d1/scfac, e1/scfac, fmt='o-',ms=1.2, mew=0.2, capsize=0.4, |
| elinewidth=0.2, color='tab:blue', mfc='none',label='大肠杆菌') |
| plt.xlabel('时间 (hr)') |
| plt.ylabel('d光密度/dt (a.u.)') |
| ylim = plt.ylim() |
| plt.plot([p1[0],p1[0]],ylim, '--', color='green', zorder=99) |
| plt.ylim(ylim) |
| plt.legend() |
| plt.tight_layout() |
| plt.savefig('./od值一阶导数.pdf') |
| |
| |
| plt.figure(figsize=(5,4)) |
| scfac = 1e5 |
| t1, s1, e1 = calc_stats('s1_fl.txt', 'blank_fl.txt') |
| plt.errorbar(t1, s1/scfac, e1/scfac, fmt='o-',ms=1.2, mew=0.2, capsize=0.4, |
| elinewidth=0.2, color='tab:red', mfc='none',label='大肠杆菌EGFP') |
| plt.xlabel('时间 (hr)') |
| plt.ylabel('荧光强度 (*10^5a.u.)') |
| plt.legend() |
| plt.tight_layout() |
| plt.savefig('./fl值.pdf') |
| |
| t,s,e = calc_stats('s1_fl.txt', 'blank_fl.txt') |
| sample = remove_baseline_from_data('s1_fl.txt', 'blank_fl.txt') |
| ys = [sample[:,col] for col in range(sample.shape[-1])] |
| t1, d1, e1 = calc_deriv(t, ys, deriv_order=2, sw=11) |
| p1 = get_peak_info(t1, d1)[0] |
| scfac = 1e5 |
| plt.figure(figsize=(5,4)) |
| plt.errorbar(t1, d1/scfac, e1/scfac, fmt='o-',ms=1.2, mew=0.2, capsize=0.4, |
| elinewidth=0.2, color='tab:blue', mfc='none',label='大肠杆菌EGFP') |
| plt.xlabel('时间 (hr)') |
| plt.ylabel('d荧光强度/dt (*10^5a.u.)') |
| ylim = plt.ylim() |
| plt.plot([p1[0],p1[0]],ylim, '--', color='black', zorder=99) |
| plt.ylim(ylim) |
| plt.legend() |
| plt.tight_layout() |
| plt.savefig('./fl的一阶导数.pdf') |
| |
| |
| import numpy as np |
| import pandas as pd |
| import matplotlib.pyplot as plt |
| from scipy import optimize as op |
| import math |
| plt.rcParams['font.sans-serif'] = ['SimHei'] |
| plt.rcParams['axes.unicode_minus'] = False |
| |
| |
| def function(x, R0, Nm): |
| return Nm/(1+2.5*Nm*float(math.e)**(-R0*x)-float(math.e)**(-R0*x)) |
| |
| |
| x_group = np.arange(0, 171/6, 1/6) |
| y_group = np.loadtxt('mean_od.txt') |
| |
| |
| R0, Nm = op.curve_fit(function, x_group, y_group)[0] |
| print(R0, Nm) |
| |
| |
| x = x_group |
| y = Nm/(1+2.5*Nm*float(math.e)**(-R0*x)-float(math.e)**(-R0*x)) |
| plt.scatter(x_group, y_group, marker='o',label='实验值') |
| plt.plot(x, y,color='red',label='拟合曲线') |
| plt.xlabel('时间(hr)') |
| plt.ylabel('菌群浓度(10^8 CFU/ml)') |
| plt.legend() |
| plt.show() |
| |
| |
| import numpy as np |
| import pandas as pd |
| import matplotlib.pyplot as plt |
| from scipy import optimize as op |
| import math |
| plt.rcParams['font.sans-serif'] = ['SimHei'] |
| plt.rcParams['axes.unicode_minus'] = False |
| |
| |
| def function(x, R0, Nm): |
| return Nm/(1+2.5*Nm*float(math.e)**(-R0*x)-float(math.e)**(-R0*x)) |
| |
| |
| x_group = np.arange(0, 171/6, 1/6) |
| y_group = np.loadtxt('mean_od.txt') |
| |
| |
| R0, Nm = op.curve_fit(function, x_group, y_group)[0] |
| print(R0, Nm) |
| |
| |
| x = x_group |
| y = Nm/(1+2.5*Nm*float(math.e)**(-R0*x)-float(math.e)**(-R0*x)) |
| plt.scatter(x_group, y_group, marker='o',label='实验值') |
| plt.plot(x, y,color='red',label='拟合曲线') |
| plt.xlabel('时间(hr)') |
| plt.ylabel('菌群浓度(10^8 CFU/ml)') |
| plt.legend() |
| plt.show() |
| |
| import math |
| import numpy as np |
| import matplotlib.pyplot as plt |
| import math |
| |
| def f(x, y): |
| return np.log(y*math.e**(y/1460000000))+(1/8.5)*(250/1460000000+np.log(2))*x-np.log(500)-500/1460000000 |
| |
| x = np.linspace(0, 100, 1001) |
| y = np.linspace(0, 500, 1001) |
| X, Y = np.meshgrid(x, y) |
| Z = f(X, Y) |
| |
| fig, ax = plt.subplots(figsize=(8, 8)) |
| ax.contour(X, Y, Z, levels=[0], colors=('red',), linewidths=1.5) |
| |
| plt.xlabel('Time(hr)',loc='center',font='Microsoft YaHei',fontsize='20') |
| plt.ylabel('Intestinal drug concentration(ng/mL)',loc='center',font='Microsoft YaHei',fontsize='20') |
| fig.savefig(r'C:\Users\Lenovo\Desktop\文献\文章\单次给药后药物浓度.pdf') |
| import numpy as np |
| import matplotlib.pyplot as plt |
| import math |
| from scipy import optimize |
| x=np.linspace(0,51,100) |
| y=np.linspace(0,1000,100) |
| k=1460000000 |
| a=math.e**(np.log(3500)+875/k) |
| def f4(t): |
| return np.array(1460000000*np.log(937.5)+937.5+3*t/2+1460000000*np.log(8)-1460000000*np.log(t)) |
| t = optimize.fsolve(f4, [2, 6]) |
| p = t[1] |
| c1=968.75 |
| t1=4 |
| def f4(b): |
| return np.array(k*np.log(c1)+c1+t1*(b/2+k*np.log(2))-k*np.log(b)-b) |
| b = optimize.fsolve(f4, [2, 6]) |
| m = b[1] |
| c2=984.375 |
| t2=5 |
| def f4(d): |
| return np.array(k*np.log(c2)+c2+t2*(d/2+k*np.log(2))-k*np.log(d)-d) |
| e = optimize.fsolve(f4, [2, 6]) |
| n = e[1] |
| def f(x,y): |
| if 0<=x and x<8.5: |
| return np.log(y*math.e**(y/1460000000))+(1/8.5)*(250/1460000000+np.log(2))*x-np.log(500)-500/1460000000 |
| if 8.5<=x and x<17: |
| return np.log((y)*math.e**((y)/1460000000))+(1/8.5)*(750/1460000000+np.log(2))*x-np.log(1500)-1500/1460000000 |
| if 17<=x and x<25.5: |
| return np.log((y)*math.e**((y)/1460000000))+(1/8.5)*(a/2/1460000000+np.log(2))*x-np.log(a)-a/1460000000 |
| if 25.5<=x and x<34: |
| return np.log((y)*math.e**((y)/1460000000))+(1/8.5)*(p/2/1460000000+np.log(2))*x-np.log(p)-p/1460000000 |
| if 34<=x and x<42.5: |
| return np.log((y)*math.e**((y)/1460000000))+(1/8.5)*(m/2/1460000000+np.log(2))*x-np.log(m)-m/1460000000 |
| if 42.5<=x and x<51: |
| return np.log((y)*math.e**((y)/1460000000))+(1/8.5)*(n/2/1460000000+np.log(2))*x-np.log(n)-n/1460000000 |
| else: |
| return 0 |
| f_vectorized = np.vectorize(f) |
| X, Y = np.meshgrid(x, y) |
| Z = f_vectorized(X, Y) |
| fig, ax = plt.subplots(figsize=(8, 8)) |
| ax.contour(X, Y, Z, levels=[0], colors=('red',), linewidths=1.5) |
| plt.xlabel('Time(hr)',loc='center',font='Microsoft YaHei',fontsize='20') |
| plt.ylabel('Intestinal drug concentration(ng/mL)',loc='center',font='Microsoft YaHei',fontsize='20') |
| |
| plt.show() |
| fig.savefig(r'C:\Users\Lenovo\Desktop\文献\文章\多次给药后药物浓度.pdf') |
| |
| import numpy as np |
| import pandas as pd |
| import matplotlib.pyplot as plt |
| from scipy import optimize as op |
| import math |
| plt.rcParams['font.sans-serif'] = ['SimHei'] |
| plt.rcParams['axes.unicode_minus'] = False |
| |
| |
| x_group = [0,2,4,6,8,10,12] |
| y_group = [120800,17360,1415.333333,496.6666667,184,119.6666667,72.33333333] |
| |
| |
| h = int(x_group[1] - x_group[0]) |
| |
| |
| dy = np.gradient(y_group, h) |
| print(dy) |
| |
| def function(y, r0, Nm, R): |
| return r0*(1-y/Nm)*y-R*y |
| r0, Nm, R = op.curve_fit(function, y_group, dy)[0] |
| print(r0, Nm, R) |
| |
| |
| y = [] |
| for i in range(0,7): |
| y.append(r0*(1-float(y_group[i])/Nm)*float(y_group[i])-R*float(y_group[i])) |
| |
| |
| plt.scatter(y_group, dy, marker='o',label='Experimental Value', color='orange') |
| plt.plot(y_group, y,color='purple',label='Fitted Curve') |
| plt.grid(True, linewidth=0.5) |
| plt.xlabel('Flora Concentration(cfu/ml)') |
| plt.ylabel('d(Flora Concentration)') |
| plt.tick_params(axis='x', labelsize=10) |
| plt.tick_params(axis='y', labelsize=10) |
| plt.grid(True, linewidth=0.5) |
| plt.rcParams.update({'font.size': 10}) |
| plt.legend() |
| plt.savefig('./涂板拟合.pdf') |
| |