diff --git a/contextualized/analysis/pvals.py b/contextualized/analysis/pvals.py index 23a2519..2d2b523 100644 --- a/contextualized/analysis/pvals.py +++ b/contextualized/analysis/pvals.py @@ -15,7 +15,7 @@ from contextualized.easy.wrappers import SKLearnWrapper -def get_pval_range(num_bootstraps: int) -> list: +def get_possible_pvals(num_bootstraps: int) -> list: """ Get the range of possible p-values based on the number of bootstraps. @@ -30,6 +30,29 @@ def get_pval_range(num_bootstraps: int) -> list: return [min_pval, max_pval] +def _validate_args(n_bootstraps: int, verbose: bool = False) -> None: + """ + Check that the test has a sufficient number of bootstraps. + + Args: + num_bootstraps (int): The number of bootstraps. + + Raises: + ValueError: If the number of bootstraps is less than 2. + """ + if n_bootstraps < 2: + raise ValueError(f"P-values are not well defined without multiple bootstrap samples.") + min_pval, max_pval = get_possible_pvals(n_bootstraps) + if verbose: + print( + "########################################################################################\n" + f"You are testing a model which contains {n_bootstraps} bootstraps.\n" + f"The minimum possible p-value is {min_pval}.\n" + f"To allow for lower p-values, increase the model's n_bootstraps.\n" + "########################################################################################" + ) + + def calc_pval_bootstraps_one_sided(estimates, thresh=0, laplace_smoothing=1): """ Calculate p-values from bootstrapped estimates. @@ -80,18 +103,11 @@ def calc_homogeneous_context_effects_pvals( Returns: np.ndarray: P-values of shape (n_contexts, n_outcomes) testing whether the sign of the direct effect of context on outcomes is consistent across bootstraps. + + Raises: + ValueError: If the model's n_bootstraps is less than 2. """ - if model.n_bootstraps < 2: - print("P values are not well defined without multiple bootstrap samples.") - return None - if verbose: - print( - "########################################################################################\n" - f"This model contains {model.n_bootstraps} bootstraps.\n" - f"The range of possible p-values is {get_pval_range(model.n_bootstraps)}.\n" - f"To get lower p-values, increase the model's n_bootstraps.\n" - "########################################################################################" - ) + _validate_args(model.n_bootstraps, verbose = verbose) _, effects = get_homogeneous_context_effects(model, C, **kwargs) # effects.shape: (n_contexts, n_bootstraps, n_context_vals, n_outcomes) diffs = effects[:, :, -1] - effects[:, :, 0] # Test whether the sign is consistent @@ -126,18 +142,11 @@ def calc_homogeneous_predictor_effects_pvals( Returns: np.ndarray: P-values of shape (n_predictors, n_outcomes) testing whether the sign of the context-invariant predictor effects are consistent across bootstraps. + + Raises: + ValueError: If the model's n_bootstraps is less than 2. """ - if model.n_bootstraps < 2: - print("P values are not well defined without multiple bootstrap samples.") - return None - if verbose: - print( - "########################################################################################\n" - f"This model contains {model.n_bootstraps} bootstraps.\n" - f"The range of possible p-values is {get_pval_range(model.n_bootstraps)}.\n" - f"To get lower p-values, increase the model's n_bootstraps.\n" - "########################################################################################" - ) + _validate_args(model.n_bootstraps, verbose = verbose) _, effects = get_homogeneous_predictor_effects(model, C, **kwargs) # effects.shape: (n_predictors, n_bootstraps, n_outcomes) pvals = np.array( @@ -171,18 +180,11 @@ def calc_heterogeneous_predictor_effects_pvals( Returns: np.ndarray: P-values of shape (n_contexts, n_predictors, n_outcomes) testing whether the context-varying parameter range is consistent across bootstraps. + + Raises: + ValueError: If the model's n_bootstraps is less than 2. """ - if model.n_bootstraps < 2: - print("P values are not well defined without multiple bootstrap samples.") - return None - if verbose: - print( - "########################################################################################\n" - f"This model contains {model.n_bootstraps} bootstraps.\n" - f"The range of possible p-values is {get_pval_range(model.n_bootstraps)}.\n" - f"To get lower p-values, increase the model's n_bootstraps.\n" - "########################################################################################" - ) + _validate_args(model.n_bootstraps, verbose = verbose) _, effects = get_heterogeneous_predictor_effects(model, C, **kwargs) # effects.shape is (n_contexts, n_predictors, n_bootstraps, n_context_vals, n_outcomes) diffs = ( @@ -233,6 +235,9 @@ def test_each_context( Returns: pd.DataFrame: A DataFrame of p-values for each (context, predictor, outcome) combination, describing how much the predictor's effect on the outcome varies across the context. + + Raises: + ValueError: If the model's n_bootstraps is less than 2. """ default_fit_params = { "encoder_type": "mlp", @@ -247,15 +252,7 @@ def test_each_context( "Target": [], "Pvals": [], } - if verbose: - print( - "########################################################################################\n" - f'Running a model with {fit_params["n_bootstraps"]} bootstraps for each of the {len(C.columns)} contexts.\n' - f'The range of possible p-values is {get_pval_range(fit_params["n_bootstraps"])}.\n' - f"To get lower p-values, set n_bootstraps to a higher value.\n" - "########################################################################################" - ) - + _validate_args(fit_params["n_bootstraps"], verbose = verbose) for context in C.columns: context_col = C[[context]].values model = model_constructor(**fit_params) diff --git a/docs/demos/sequential_testing.ipynb b/docs/demos/sequential_testing.ipynb index eff9d1b..80f68cb 100644 --- a/docs/demos/sequential_testing.ipynb +++ b/docs/demos/sequential_testing.ipynb @@ -13,13 +13,13 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import pandas as pd\n", - "from contextualized.analysis.pvals import test_each_context, get_pval_range\n", + "from contextualized.analysis.pvals import test_each_context, get_possible_pvals\n", "from contextualized.easy import ContextualizedRegressor\n", "\n", "import logging\n", @@ -36,7 +36,7 @@ }, { "cell_type": "code", - "execution_count": 2, + "execution_count": 3, "metadata": {}, "outputs": [], "source": [ @@ -72,7 +72,7 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 4, "metadata": {}, "outputs": [ { @@ -81,19 +81,19 @@ "[0.024390243902439025, 0.975609756097561]" ] }, - "execution_count": 3, + "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# getting the range of possible p-values for 40 bootstraps\n", - "get_pval_range(40)" + "get_possible_pvals(40)" ] }, { "cell_type": "code", - "execution_count": 4, + "execution_count": 5, "metadata": {}, "outputs": [ { @@ -129,7 +129,7 @@ " C0\n", " X0\n", " Y\n", - " 0.292683\n", + " 0.170732\n", " \n", " \n", " 1\n", @@ -143,14 +143,14 @@ " C1\n", " X0\n", " Y\n", - " 0.439024\n", + " 0.341463\n", " \n", " \n", " 3\n", " C1\n", " X1\n", " Y\n", - " 0.365854\n", + " 0.390244\n", " \n", " \n", "\n", @@ -158,13 +158,13 @@ ], "text/plain": [ " Context Predictor Target Pvals\n", - "0 C0 X0 Y 0.292683\n", + "0 C0 X0 Y 0.170732\n", "1 C0 X1 Y 0.024390\n", - "2 C1 X0 Y 0.439024\n", - "3 C1 X1 Y 0.365854" + "2 C1 X0 Y 0.341463\n", + "3 C1 X1 Y 0.390244" ] }, - "execution_count": 4, + "execution_count": 5, "metadata": {}, "output_type": "execute_result" } @@ -187,7 +187,7 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": 6, "metadata": {}, "outputs": [], "source": [ @@ -210,37 +210,37 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "%%capture\n", - "pvals = test_each_context(ContextualizedRegressor, C_train_df, X_train_df, Y_train_df, encoder_type=\"mlp\", max_epochs=20, learning_rate=1e-2, n_bootstraps=20)" + "pvals = test_each_context(ContextualizedRegressor, C_train_df, X_train_df, Y_train_df, encoder_type=\"mlp\", max_epochs=20, learning_rate=1e-2, n_bootstraps=40)" ] }, { "cell_type": "code", - "execution_count": 7, + "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ - "[0.047619047619047616, 0.9523809523809523]" + "[0.024390243902439025, 0.975609756097561]" ] }, - "execution_count": 7, + "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ - "get_pval_range(20)" + "get_possible_pvals(40)" ] }, { "cell_type": "code", - "execution_count": 8, + "execution_count": 9, "metadata": {}, "outputs": [ { @@ -276,147 +276,147 @@ " age\n", " bp\n", " 0\n", - " 0.428571\n", + " 0.560976\n", " \n", " \n", " 1\n", " age\n", " s1\n", " 0\n", - " 0.428571\n", + " 0.341463\n", " \n", " \n", " 2\n", " age\n", " s2\n", " 0\n", - " 0.523810\n", + " 0.512195\n", " \n", " \n", " 3\n", " age\n", " s3\n", " 0\n", - " 0.428571\n", + " 0.560976\n", " \n", " \n", " 4\n", " age\n", " s4\n", " 0\n", - " 0.476190\n", + " 0.560976\n", " \n", " \n", " 5\n", " age\n", " s5\n", " 0\n", - " 0.428571\n", + " 0.560976\n", " \n", " \n", " 6\n", " age\n", " s6\n", " 0\n", - " 0.428571\n", + " 0.560976\n", " \n", " \n", " 7\n", " sex\n", " bp\n", " 0\n", - " 0.095238\n", + " 0.170732\n", " \n", " \n", " 8\n", " sex\n", " s1\n", " 0\n", - " 0.666667\n", + " 0.414634\n", " \n", " \n", " 9\n", " sex\n", " s2\n", " 0\n", - " 0.523810\n", + " 0.609756\n", " \n", " \n", " 10\n", " sex\n", " s3\n", " 0\n", - " 0.095238\n", + " 0.170732\n", " \n", " \n", " 11\n", " sex\n", " s4\n", " 0\n", - " 0.190476\n", + " 0.170732\n", " \n", " \n", " 12\n", " sex\n", " s5\n", " 0\n", - " 0.095238\n", + " 0.170732\n", " \n", " \n", " 13\n", " sex\n", " s6\n", " 0\n", - " 0.095238\n", + " 0.170732\n", " \n", " \n", " 14\n", " bmi\n", " bp\n", " 0\n", - " 0.047619\n", + " 0.024390\n", " \n", " \n", " 15\n", " bmi\n", " s1\n", " 0\n", - " 0.380952\n", + " 0.390244\n", " \n", " \n", " 16\n", " bmi\n", " s2\n", " 0\n", - " 0.333333\n", + " 0.219512\n", " \n", " \n", " 17\n", " bmi\n", " s3\n", " 0\n", - " 0.047619\n", + " 0.024390\n", " \n", " \n", " 18\n", " bmi\n", " s4\n", " 0\n", - " 0.095238\n", + " 0.097561\n", " \n", " \n", " 19\n", " bmi\n", " s5\n", " 0\n", - " 0.047619\n", + " 0.024390\n", " \n", " \n", " 20\n", " bmi\n", " s6\n", " 0\n", - " 0.047619\n", + " 0.024390\n", " \n", " \n", "\n", @@ -424,30 +424,30 @@ ], "text/plain": [ " Context Predictor Target Pvals\n", - "0 age bp 0 0.428571\n", - "1 age s1 0 0.428571\n", - "2 age s2 0 0.523810\n", - "3 age s3 0 0.428571\n", - "4 age s4 0 0.476190\n", - "5 age s5 0 0.428571\n", - "6 age s6 0 0.428571\n", - "7 sex bp 0 0.095238\n", - "8 sex s1 0 0.666667\n", - "9 sex s2 0 0.523810\n", - "10 sex s3 0 0.095238\n", - "11 sex s4 0 0.190476\n", - "12 sex s5 0 0.095238\n", - "13 sex s6 0 0.095238\n", - "14 bmi bp 0 0.047619\n", - "15 bmi s1 0 0.380952\n", - "16 bmi s2 0 0.333333\n", - "17 bmi s3 0 0.047619\n", - "18 bmi s4 0 0.095238\n", - "19 bmi s5 0 0.047619\n", - "20 bmi s6 0 0.047619" + "0 age bp 0 0.560976\n", + "1 age s1 0 0.341463\n", + "2 age s2 0 0.512195\n", + "3 age s3 0 0.560976\n", + "4 age s4 0 0.560976\n", + "5 age s5 0 0.560976\n", + "6 age s6 0 0.560976\n", + "7 sex bp 0 0.170732\n", + "8 sex s1 0 0.414634\n", + "9 sex s2 0 0.609756\n", + "10 sex s3 0 0.170732\n", + "11 sex s4 0 0.170732\n", + "12 sex s5 0 0.170732\n", + "13 sex s6 0 0.170732\n", + "14 bmi bp 0 0.024390\n", + "15 bmi s1 0 0.390244\n", + "16 bmi s2 0 0.219512\n", + "17 bmi s3 0 0.024390\n", + "18 bmi s4 0 0.097561\n", + "19 bmi s5 0 0.024390\n", + "20 bmi s6 0 0.024390" ] }, - "execution_count": 8, + "execution_count": 9, "metadata": {}, "output_type": "execute_result" }