<jupyter_start><jupyter_text>Developing latitudinal models of $T_\mathrm{opt}, CT_\textrm{min}$ and $CT_\textrm{max}$<jupyter_code>plankton_o = pd.read_csv("../../data/Phytoplankton_temperature_growth_rate_dataset_2016_01_29/traits_derived_2016_01_29.csv", engine='python')<jupyter_output><empty_output><jupyter_text>**Filter** the data by fit quality (perhaps) and marine:<jupyter_code># plankton = plankton_o[(plankton_o.minqual == "good") &
#                     (plankton_o.maxqual == "good") &
#                     (plankton_o.curvequal == "good")]
plankton = plankton_o[plankton_o.habitat == 'marine']<jupyter_output><empty_output><jupyter_text>**Drop NAs**<jupyter_code>plankton = plankton.dropna(
    axis=0, 
    subset=[
        'isolation.latitude',
        'isolation.longitude',
        'mu.g.opt.list', 
        'tmin', 
        'tmax'
    ]
)<jupyter_output><empty_output><jupyter_text>Relevant parameter names/descriptions (from `../../data/Phytoplankton_temperature_growth_rate_dataset_2016_01_29/Dataset explanation.doc`)```23) mu.wlist = estimated thermal niche width (parameter ‘ɷ’ in the thermal reaction norm model)24) mu.alist = estimate of parameter ‘a’ in the thermal reaction norm model25) mu.blist = estimate of parameter ‘b’ in the thermal reaction norm model26) mu.slist = variance parameter for the maximum likelihood model fit. 27) mu.c.opt.list = estimate of parameter ‘z’ in the thermal reaction norm model28) mu.g.opt.val.list = estimated specific growth rate (per day) when temperature is at ‘z’ (i.e. mu.c.opt.list)29) mu.g.opt.list = estimated optimum temperature for growth 30) mu.g.opt.val.list = estimated maximum specific growth rate (per day) based on the thermal reaction norm model fit31) mu.n = number of points (i.e. number of growth rate measurements) in the curve 32) emp.max.growth = maximum specific growth rate (per day) measured during the growth assays. These were used for maximum growth rate analyses, but results did not differ significantly if estimated maximum specific growth rate based on the thermal reaction norm model fit (i.e. mu.g.opt.val.list) was used instead.33) tmin = Tmin, or minimum persistence temperature, estimated from the thermal reaction norm model fit34) tmax = Tmax, or maximum persistence temperature, estimated from the thermal reaction norm model fit35) minqual = quality of Tmin estimate (quality control criteria found in supplementary info).36) maxqual = quality of Tmax estimate (quality control criteria found in supplementary info).37) curvequal = quality of niche width estimate (quality control criteria found in supplementary info).38) abs.curveskew = Estimated absolute skewness of the thermal reaction norm39) rel.curveskew = Estimated relative skewness of the thermal reaction norm```  $T_\mathrm{opt}$<jupyter_code>topt_colname = 'mu.g.opt.list' ## 29) mu.g.opt.list = estimated optimum temperature for growth 
sns.lmplot(
    x = 'isolation.latitude',
    y = topt_colname,
    size = 6, 
    order = 2, 
    data = plankton
)
plt.title("$T_\mathrm{opt}$ by Latitude ($N = " + str(len(plankton)) + ")$")
plt.xlabel("Degrees Latitude")
plt.ylabel("T [deg C]")<jupyter_output>/home/ec2-user/miniconda3/envs/notebook/lib/python3.7/site-packages/seaborn/regression.py:546: UserWarning: The `size` paramter has been renamed to `height`; please update your code.
  warnings.warn(msg, UserWarning)<jupyter_text>`statsmodels` fit<jupyter_code>lat_column = 'isolation.latitude'

topt = plankton[topt_colname]

order2 = PolynomialFeatures(2).fit_transform(plankton[lat_column].to_numpy().reshape(-1,1))

topt_model = sm.OLS(topt, order2).fit()

topt_model.summary(title="Regular Latitude (Second order fit) Results")<jupyter_output><empty_output><jupyter_text>$T_\mathrm{min}$<jupyter_code>sns.lmplot(
    x = 'isolation.latitude',
    y = 'tmin',
    size = 6, 
    order = 2, 
    data = plankton
)
plt.title("$T_\mathrm{min}$ by Latitude ($N = " + str(len(plankton)) + ")$")
plt.xlabel("Degrees Latitude")
plt.ylabel("T [deg C]")<jupyter_output>/home/ec2-user/miniconda3/envs/notebook/lib/python3.7/site-packages/seaborn/regression.py:546: UserWarning: The `size` paramter has been renamed to `height`; please update your code.
  warnings.warn(msg, UserWarning)<jupyter_text>**Something off about that. might need a filter here**.<jupyter_code>plankton_filtered = plankton[
    plankton.curvequal == 'good'
]
sns.lmplot(
    x = 'isolation.latitude',
    y = 'tmin',
    size = 6, 
    order = 2, 
    data = plankton_filtered
)
plt.title("$T_\mathrm{min}$ by Latitude ($N = " + str(len(plankton_filtered)) + ")$")
plt.xlabel("Degrees Latitude")
plt.ylabel("T [deg C]")<jupyter_output><empty_output><jupyter_text>**Perhaps more reasonable**. compare `statsmodels` fits **1) No Data Filter**<jupyter_code>lat_column = 'isolation.latitude'

tmin = plankton['tmin']

order2 = PolynomialFeatures(2).fit_transform(plankton[lat_column].to_numpy().reshape(-1,1))

tmin_model_nofilter = sm.OLS(tmin, order2).fit()

tmin_model_nofilter.summary(title="Regular Latitude (Second order fit) Results")<jupyter_output><empty_output><jupyter_text>**2) With filtered data**<jupyter_code>lat_column = 'isolation.latitude'

tmin = plankton_filtered['tmin']

order2 = PolynomialFeatures(2).fit_transform(plankton_filtered[lat_column].to_numpy().reshape(-1,1))

tmin_model_filter = sm.OLS(tmin, order2).fit()

tmin_model_filter.summary(title="Regular Latitude (Second order fit) Results")<jupyter_output><empty_output><jupyter_text>**Compare Fits:**<jupyter_code>pd.merge(
    tmin_model_filter.params.rename("Filtered"),
    tmin_model_nofilter.params.rename("Not Filtered"),
    left_index=True, right_index=True
)<jupyter_output><empty_output><jupyter_text>**Compare $R^2$**<jupyter_code>print(
    f"Filtered R2: {tmin_model_filter.rsquared:.3f}",
    f"Non-Filtered R2: {tmin_model_nofilter.rsquared:.3f}"
)<jupyter_output>Filtered R2: 0.410 Non-Filtered R2: 0.084<jupyter_text>$T_\mathrm{max}$<jupyter_code>sns.lmplot(
    x = 'isolation.latitude',
    y = 'tmax',
    size = 6, 
    order = 2, 
    data = plankton
)
plt.title("$T_\mathrm{max}$ by Latitude ($N = " + str(len(plankton)) + ")$")
plt.xlabel("Degrees Latitude")
plt.ylabel("T [deg C]")
sns.lmplot(
    x = 'isolation.latitude',
    y = 'tmax',
    size = 6, 
    order = 2, 
    data = plankton_filtered
)
plt.title("$T_\mathrm{max}$ by Latitude ($N = " + str(len(plankton_filtered)) + ")$")
plt.xlabel("Degrees Latitude")
plt.ylabel("T [deg C]")<jupyter_output><empty_output><jupyter_text>Compare `statsmodels` fits**1) No Data Filter**<jupyter_code>lat_column = 'isolation.latitude'

tmax = plankton['tmax']

order2 = PolynomialFeatures(2).fit_transform(plankton[lat_column].to_numpy().reshape(-1,1))

tmax_model_nofilter = sm.OLS(tmax, order2).fit()

tmax_model_nofilter.summary(title="Regular Latitude (Second order fit) Results")<jupyter_output><empty_output><jupyter_text>**2) With Data Filter**<jupyter_code>lat_column = 'isolation.latitude'

tmax = plankton_filtered['tmax']

order2 = PolynomialFeatures(2).fit_transform(plankton_filtered[lat_column].to_numpy().reshape(-1,1))

tmax_model_filter = sm.OLS(tmax, order2).fit()

tmax_model_filter.summary(title="Regular Latitude (Second order fit) Results")<jupyter_output><empty_output><jupyter_text>**Compare Fits:**<jupyter_code>pd.merge(
    tmax_model_filter.params.rename("Filtered"),
    tmax_model_nofilter.params.rename("Not Filtered"),
    left_index=True, right_index=True
)<jupyter_output><empty_output><jupyter_text>Similar. **Compare R^2**:<jupyter_code>print(
    f"Filtered R2: {tmax_model_filter.rsquared:.3f}",
    f"Non-Filtered R2: {tmax_model_nofilter.rsquared:.3f}"
)<jupyter_output>Filtered R2: 0.481 Non-Filtered R2: 0.197<jupyter_text>Develop Generalized TPC Class<jupyter_code>class GeneralizedTPC(object):
    def __init__(self, toptModel, tminModel, tmaxModel):
        self.toptModel = deepcopy(toptModel)
        self.tminModel = deepcopy(tminModel)
        self.tmaxModel = deepcopy(tmaxModel)
        
    def getTPCParameters(self, latitude):
        import numpy as np
        from sklearn.preprocessing import PolynomialFeatures
        latitude = PolynomialFeatures(2).fit_transform(np.array(latitude).reshape(-1, 1))
        _topt = self.toptModel.predict(latitude).item()
        _tmin = self.tminModel.predict(latitude).item()
        _tmax = self.tmaxModel.predict(latitude).item()
        return (_topt, _tmin, _tmax)
    
    def getLatitudinalTPC(self, latitude):
        _topt, _tmin, _tmax = self.getTPCParameters(latitude)
        return(partial(
            tpc.tpc_general, 
            Topt = _topt, 
            CTmin = _tmin, 
            CTmax = _tmax
        ))<jupyter_output><empty_output><jupyter_text>Instantiate + Save to FileWe "ossify" the general class by instantiating it with the models derived in this notebook, above.<jupyter_code>## "ossify"
gtpc = GeneralizedTPC(topt_model, tmin_model_filter, tmax_model_filter)<jupyter_output><empty_output><jupyter_text>To save this to a file we're forced to use `cloudpickle` since it's the only available serialization package with support for interactively-defined user classes.<jupyter_code>with open("gtpc_modeled.pkl", 'wb') as f:
    cloudpickle.dump(gtpc, f)<jupyter_output><empty_output><jupyter_text>---  Approach 2: Decomposition/Multivariate OrdinationReference: https://github.com/HuckleyLab/ThermalStress/blob/master/ToptCTmax_analysis.R<jupyter_code>axes = plankton[[topt_colname, 'tmin', 'tmax']]
pca = PCA()
X = pca.fit_transform(axes)
pca.components_.T
axes.columns
weights = pd.DataFrame(
    pca.components_.T, 
    columns=[f'PC{i+1}' for i in range(pca.components_.T.shape[1])],
    index=axes.columns)
weights
loadings = pca.components_.T * np.sqrt(pca.explained_variance_)
loading_matrix = pd.DataFrame(
    loadings, 
    columns=[f'PC{i+1}' for i in range(pca.components_.T.shape[1])],
    index=axes.columns)
loading_matrix<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13138905592303066, 'fraction_numerical': 0.015233513730206455, 'mean_word_length': 5.240775484677924, 'pattern_counts': {'":': 0, '<': 77, '<?xml version=': 0, '>': 77, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '20', 'dir': 'jupyter-structured-clean-dedup', 'id': '32110', 'license': 'MIT', 'n_tokens_mistral': 3692, 'n_tokens_neox': 3475, 'n_words': 707, 'path': 'analysis/physiology/TPC_params_latitudinal_generalization.ipynb', 'repo_name': 'HuckleyLab/phyto-mhw', 'size': '9978'}