importnumpy asnpchem_data = [( 0.00080, 99),( 0.00800, 91),( 0.08000, 89),( 0.40000, 89),( 0.80000, 79),( 1.60000, 61),( 4.00000, 39),( 8.00000, 25),( 80.00000, 4)] importpandas aspdchem_df = pd.DataFrame(chem_data)chem_df.columns = [ 'concentration', 'activity']chem_df[ 'concentration_log'] = chem_df[ 'concentration'].apply( lambdax:np.log10(x))

# df.set_index('concentration', inplace=True)参数问题

给出数据，化学品的IC50 值是多少, 以及其周围的不确定性？

说明

数据defplot_chemical_data(log=True):fig = plt.figure(figsize=( 10, 6)) ax = fig.add_subplot( 1, 1, 1)

iflog: ax.scatter(x=chem_df[ 'concentration_log'], y=chem_df[ 'activity']) ax.set_xlabel( 'log10(concentration (mM))', fontsize= 20) else: ax.scatter(x=chem_df[ 'concentration'], y=chem_df[ 'activity']) ax.set_xlabel( 'concentration (mM)', fontsize= 20) ax.set_xticklabels([int(i) fori inax.get_xticks()], fontsize= 18) ax.set_yticklabels([int(i) fori inax.get_yticks()], fontsize= 18) plt.hlines(y= 50, xmin=min(ax.get_xlim()), xmax=max(ax.get_xlim()), linestyles= '--',) returnfigfig = plot_chemical_data(log= True)plt.show()

withpm.Model() asic50_model: beta = pm.HalfNormal( 'beta', sd= 100** 2) ic50_log10 = pm.Flat( 'IC50_log10') # Flat prior# MATH WITH DISTRIBUTION OBJECTS!measurements = beta / ( 1+ np.exp(chem_df[ 'concentration_log'].values - ic50_log10)) y_like = pm.Normal( 'y_like', mu=measurements, observed=chem_df[ 'activity'])

# Deterministic transformations.ic50 = pm.Deterministic( 'IC50', np.power( 10, ic50_log10)) MCMC Inference Button (TM)withic50_model: step = pm.Metropolis() ic50_trace = pm.sample( 10000, step=step) pm.traceplot(ic50_trace[ 2000:], varnames=[ 'IC50_log10', 'IC50']) # live: sample from step 2000 onwards.

plt.show()

结果pm.plot_posterior(ic50_trace[ 4000:], varnames=[ 'IC50'], color= '#87ceeb', point_estimate= 'mean')plt.show()

该化学物质的IC50在约 [2mM，2.4mM]（95％HPD）。 这是一种不好的化学物质。

问题类型2：实验组之间的比较

实验组和对照组的不同

例子1：药物IQ问题

药物治疗是否影响 IQ Scores

drug = [ 99., 110., 107., 104., 省略]placebo = [ 95., 105., 103., 99., 省略]

defECDF(data):x = np.sort(data) y = np.cumsum(x) / np.sum(x)

returnx, y

defplot_drug():fig = plt.figure() ax = fig.add_subplot( 1, 1, 1) x_drug, y_drug = ECDF(drug) ax.plot(x_drug, y_drug, label= 'drug, n={0}'.format(len(drug))) x_placebo, y_placebo = ECDF(placebo) ax.plot(x_placebo, y_placebo, label= 'placebo, n={0}'.format(len(placebo))) ax.legend() ax.set_xlabel( 'IQ Score') ax.set_ylabel( 'Cumulative Frequency') ax.hlines( 0.5, ax.get_xlim()[ 0], ax.get_xlim()[ 1], linestyle= '--')

returnfig fromscipy.stats importttest_indttest_ind(drug, placebo) Ttest_indResult(statistic=2.2806701634329549, pvalue=0.025011500508647616)

实验

随机将参与者分配给两个实验组：

+drug vs. -drug

测量每个参与者的 IQ Scores

fig = plot_drug()plt.show()