Continuous Data¶
Quickstart¶
To run a continuous dataset:
import pybmds
# create a continuous dataset
dataset = pybmds.ContinuousDataset(
doses=[0, 25, 50, 75, 100],
ns=[20, 20, 20, 20, 20],
means=[6, 8, 13, 25, 30],
stdevs=[4, 4.3, 3.8, 4.4, 3.7],
)
# create a session
session = pybmds.Session(dataset=dataset)
session.add_default_models()
# execute the session; recommend best-fitting model
session.execute_and_recommend()
# save excel report
df = session.to_df()
df.to_excel("output/report.xlsx")
# save to a word report
report = session.to_docx()
report.save("output/report.docx")
Continuous datasets¶
A continuous dataset requires a list of doses, mean responses, standard deviations, and the total number of subjects. All lists must be the same size, and there should be one item in each list for each dose-group.
You can also add optional attributes, such as name, dose_name, dose_units, response_name, response_units, etc.
For example:
dataset = pybmds.ContinuousDataset(
name="Body Weight from ChemX Exposure",
dose_name="Concentration",
dose_units="ppm",
response_units="kg",
doses=[0, 25, 50, 75, 100],
ns=[20, 20, 20, 20, 20],
means=[6, 8, 13, 25, 30],
stdevs=[4, 4.3, 3.8, 4.4, 3.7],
)
dataset.plot()
Single model fit¶
You can fit a specific model to the dataset and plot/print the results. The printed results will include the BMD, BMDL, BMDU, p-value, AIC, etc.
For example, to fit the Hill model, run the following code and plot the results:
from pybmds.models import continuous
model = continuous.Hill(dataset)
model.execute()
model.plot()
To view an output report:
print(model.text())
Hill Model
══════════════════════════════
Version: pybmds 24.1 (bmdscore 24.1)
Input Summary:
╒══════════════════════════════╤════════════════════════════╕
│ BMR │ 1.0 Standard Deviation │
│ Distribution │ Normal + Constant variance │
│ Modeling Direction │ Up (↑) │
│ Confidence Level (one sided) │ 0.95 │
│ Modeling Approach │ MLE │
╘══════════════════════════════╧════════════════════════════╛
Parameter Settings:
╒═════════════╤═══════════╤═══════╤═══════╕
│ Parameter │ Initial │ Min │ Max │
╞═════════════╪═══════════╪═══════╪═══════╡
│ g │ 0 │ -100 │ 100 │
│ v │ 0 │ -100 │ 100 │
│ k │ 0 │ 0 │ 5 │
│ n │ 1 │ 1 │ 18 │
│ alpha │ 0 │ -18 │ 18 │
╘═════════════╧═══════════╧═══════╧═══════╛
Modeling Summary:
╒════════════════╤═════════════╕
│ BMD │ 44.3603 │
│ BMDL │ 37.4378 │
│ BMDU │ 49.4299 │
│ AIC │ 570.3 │
│ Log-Likelihood │ -280.15 │
│ P-Value │ 0.163921 │
│ Model d.f. │ 1 │
╘════════════════╧═════════════╛
Model Parameters:
╒════════════╤════════════╤════════════╤═════════════╕
│ Variable │ Estimate │ On Bound │ Std Error │
╞════════════╪════════════╪════════════╪═════════════╡
│ g │ 6.80589 │ no │ 0.689269 │
│ v │ 25.6632 │ no │ 2.65373 │
│ k │ 62.7872 │ no │ 3.59334 │
│ n │ 4.87539 │ no │ 1.08357 │
│ alpha │ 15.881 │ no │ 35.6668 │
╘════════════╧════════════╧════════════╧═════════════╛
Goodness of Fit:
╒════════╤═════╤═══════════════╤═════════════════════╤═══════════════════╕
│ Dose │ N │ Sample Mean │ Model Fitted Mean │ Scaled Residual │
╞════════╪═════╪═══════════════╪═════════════════════╪═══════════════════╡
│ 0 │ 20 │ 6 │ 6.80589 │ -0.904381 │
│ 25 │ 20 │ 8 │ 7.09076 │ 1.02037 │
│ 50 │ 20 │ 13 │ 13.1658 │ -0.186064 │
│ 75 │ 20 │ 25 │ 24.8734 │ 0.14202 │
│ 100 │ 20 │ 30 │ 30.0641 │ -0.0719426 │
╘════════╧═════╧═══════════════╧═════════════════════╧═══════════════════╛
╒════════╤═════╤═════════════╤═══════════════════╕
│ Dose │ N │ Sample SD │ Model Fitted SD │
╞════════╪═════╪═════════════╪═══════════════════╡
│ 0 │ 20 │ 4 │ 3.98509 │
│ 25 │ 20 │ 4.3 │ 3.98509 │
│ 50 │ 20 │ 3.8 │ 3.98509 │
│ 75 │ 20 │ 4.4 │ 3.98509 │
│ 100 │ 20 │ 3.7 │ 3.98509 │
╘════════╧═════╧═════════════╧═══════════════════╛
Likelihoods:
╒═════════╤══════════════════╤════════════╤═════════╕
│ Model │ Log-Likelihood │ # Params │ AIC │
╞═════════╪══════════════════╪════════════╪═════════╡
│ A1 │ -279.181 │ 6 │ 570.362 │
│ A2 │ -278.726 │ 10 │ 577.452 │
│ A3 │ -279.181 │ 6 │ 570.362 │
│ fitted │ -280.15 │ 5 │ 570.3 │
│ reduced │ -374.79 │ 2 │ 753.579 │
╘═════════╧══════════════════╧════════════╧═════════╛
Tests of Mean and Variance Fits:
╒════════╤══════════════════════════════╤═════════════╤═══════════╕
│ Name │ -2 * Log(Likelihood Ratio) │ Test d.f. │ P-Value │
╞════════╪══════════════════════════════╪═════════════╪═══════════╡
│ Test 1 │ 192.127 │ 8 │ 0 │
│ Test 2 │ 0.909828 │ 4 │ 0.923147 │
│ Test 3 │ 0.909828 │ 4 │ 0.923147 │
│ Test 4 │ 1.93767 │ 1 │ 0.163921 │
╘════════╧══════════════════════════════╧═════════════╧═══════════╛
Test 1: Test the null hypothesis that responses and variances don't differ among dose levels
(A2 vs R). If this test fails to reject the null hypothesis (p-value > 0.05), there may not be
a dose-response.
Test 2: Test the null hypothesis that variances are homogenous (A1 vs A2). If this test fails to
reject the null hypothesis (p-value > 0.05), the simpler constant variance model may be appropriate.
Test 3: Test the null hypothesis that the variances are adequately modeled (A3 vs A2). If this test
fails to reject the null hypothesis (p-value > 0.05), it may be inferred that the variances have
been modeled appropriately.
Test 4: Test the null hypothesis that the model for the mean fits the data (Fitted vs A3). If this
test fails to reject the null hypothesis (p-value > 0.1), the user has support for use of the
selected model.
You can change the range and initial value for any parameter in the model by following the same steps above.
Change input settings¶
The default settings when running a continuous analysis are a BMR of 1 standard deviation from the BMD, a 95% confidence interval, and a normally distributed and constant variance between dose groups.
If you fit a single model to your dataset, you can change these settings:
model = continuous.Hill(dataset)
print(model.settings.tbl())
print(model.priors_tbl())
╒══════════════════════════════╤════════════════════════════╕
│ BMR │ 1.0 Standard Deviation │
│ Distribution │ Normal + Constant variance │
│ Modeling Direction │ Up (↑) │
│ Confidence Level (one sided) │ 0.95 │
│ Modeling Approach │ MLE │
│ Degree │ 0 │
╘══════════════════════════════╧════════════════════════════╛
╒═════════════╤═══════════╤═══════╤═══════╕
│ Parameter │ Initial │ Min │ Max │
╞═════════════╪═══════════╪═══════╪═══════╡
│ g │ 0 │ -100 │ 100 │
│ v │ 0 │ -100 │ 100 │
│ k │ 0 │ 0 │ 5 │
│ n │ 1 │ 1 │ 18 │
│ alpha │ 0 │ -18 │ 18 │
╘═════════════╧═══════════╧═══════╧═══════╛
You can change individual model settings:
model = continuous.Hill(
dataset,
settings={
"bmr": 0.1,
"bmr_type": pybmds.ContinuousRiskType.RelativeDeviation,
"disttype": pybmds.ContinuousDistType.normal_ncv,
},
)
model.settings.priors.update("k", max_value=3)
print(model.settings.tbl())
print(model.priors_tbl())
╒══════════════════════════════╤═══════════════════════════════╕
│ BMR │ 10% Relative Deviation │
│ Distribution │ Normal + Nonconstant variance │
│ Modeling Direction │ Up (↑) │
│ Confidence Level (one sided) │ 0.95 │
│ Modeling Approach │ MLE │
│ Degree │ 0 │
╘══════════════════════════════╧═══════════════════════════════╛
╒═════════════╤═══════════╤═══════╤═══════╕
│ Parameter │ Initial │ Min │ Max │
╞═════════════╪═══════════╪═══════╪═══════╡
│ g │ 0 │ -100 │ 100 │
│ v │ 0 │ -100 │ 100 │
│ k │ 0 │ 0 │ 3 │
│ n │ 1 │ 1 │ 18 │
│ rho │ 0 │ -18 │ 18 │
│ alpha │ 0 │ -18 │ 18 │
╘═════════════╧═══════════╧═══════╧═══════╛
For the Polynomial model, degrees can be increased to a maximum of the lesser of n-1 or 8, where n is the number of dose groups.
continuous.Polynomial(dataset, settings={"degree": 2})
continuous.Polynomial(dataset, settings={"degree": 3})
<pybmds.models.continuous.Polynomial at 0x7fd5956e6570>
Change parameter settings¶
To see a preview of the initial parameter settings:
model = continuous.Hill(dataset)
print(model.settings.tbl())
print(model.priors_tbl())
╒══════════════════════════════╤════════════════════════════╕
│ BMR │ 1.0 Standard Deviation │
│ Distribution │ Normal + Constant variance │
│ Modeling Direction │ Up (↑) │
│ Confidence Level (one sided) │ 0.95 │
│ Modeling Approach │ MLE │
│ Degree │ 0 │
╘══════════════════════════════╧════════════════════════════╛
╒═════════════╤═══════════╤═══════╤═══════╕
│ Parameter │ Initial │ Min │ Max │
╞═════════════╪═══════════╪═══════╪═══════╡
│ g │ 0 │ -100 │ 100 │
│ v │ 0 │ -100 │ 100 │
│ k │ 0 │ 0 │ 5 │
│ n │ 1 │ 1 │ 18 │
│ alpha │ 0 │ -18 │ 18 │
╘═════════════╧═══════════╧═══════╧═══════╛
You can also change the initial parameter settings shown above for any run of a single continuous model. Continuing with the Hill model example, you can set the power parameter n to be equal to 1. You do this by changing the initial, minimum, and maximum values of that parameter:
model = continuous.Hill(dataset)
model.settings.priors.update("n", initial_value=1, min_value=1, max_value=1)
print(model.priors_tbl())
╒═════════════╤═══════════╤═══════╤═══════╕
│ Parameter │ Initial │ Min │ Max │
╞═════════════╪═══════════╪═══════╪═══════╡
│ g │ 0 │ -100 │ 100 │
│ v │ 0 │ -100 │ 100 │
│ k │ 0 │ 0 │ 5 │
│ n │ 1 │ 1 │ 1 │
│ alpha │ 0 │ -18 │ 18 │
╘═════════════╧═══════════╧═══════╧═══════╛
Multiple model fit (sessions) and model recommendation¶
To run all the default models, save the results, and save the plot of the fit of the recommended model with the data:
# create a continuous dataset
dataset = pybmds.ContinuousDataset(
doses=[0, 25, 50, 75, 100],
ns=[20, 20, 20, 20, 20],
means=[6, 8, 13, 25, 30],
stdevs=[4, 4.3, 3.8, 4.4, 3.7],
)
session = pybmds.Session(dataset=dataset)
session.add_default_models()
# execute the session
session.execute()
# recommend a best-fitting model
session.recommend()
# print recommended model and plot recommended model with dataset
model_index = session.recommender.results.recommended_model_index
if model_index:
model = session.models[model_index]
model.plot()
You can also plot all models:
session.plot()
session.plot(colorize=False)
Change model settings¶
To change model settings to all default models added to a session:
session = pybmds.Session(dataset=dataset)
session.add_default_models(
settings={
"disttype": pybmds.ContinuousDistType.normal_ncv,
"bmr_type": pybmds.ContinuousRiskType.AbsoluteDeviation,
"bmr": 2,
"alpha": 0.1,
}
)
Run subset of models and select best fit¶
You can select a set of models and find the best fit, rather than using all of the default continuous models.
For example, to model average on the Hill, Linear, and Exponential 3:
session = pybmds.Session(dataset=dataset)
session.add_model(pybmds.Models.Linear)
session.add_model(pybmds.Models.ExponentialM3)
session.add_model(pybmds.Models.Hill)
session.execute()
Model recommendation and selection¶
The pybmds package may recommend a best-fitting model based on a decision tree, but expert judgment may be required for model selection. To run model recommendation and view a recommended model, if one is recommended:
session.recommend()
if session.recommended_model is not None:
display(session.recommended_model.plot())
print(session.recommended_model.text())
Hill Model
══════════════════════════════
Version: pybmds 24.1 (bmdscore 24.1)
Input Summary:
╒══════════════════════════════╤════════════════════════════╕
│ BMR │ 1.0 Standard Deviation │
│ Distribution │ Normal + Constant variance │
│ Modeling Direction │ Up (↑) │
│ Confidence Level (one sided) │ 0.95 │
│ Modeling Approach │ MLE │
╘══════════════════════════════╧════════════════════════════╛
Parameter Settings:
╒═════════════╤═══════════╤═══════╤═══════╕
│ Parameter │ Initial │ Min │ Max │
╞═════════════╪═══════════╪═══════╪═══════╡
│ g │ 0 │ -100 │ 100 │
│ v │ 0 │ -100 │ 100 │
│ k │ 0 │ 0 │ 5 │
│ n │ 1 │ 1 │ 18 │
│ alpha │ 0 │ -18 │ 18 │
╘═════════════╧═══════════╧═══════╧═══════╛
Modeling Summary:
╒════════════════╤═════════════╕
│ BMD │ 44.3603 │
│ BMDL │ 37.4378 │
│ BMDU │ 49.4299 │
│ AIC │ 570.3 │
│ Log-Likelihood │ -280.15 │
│ P-Value │ 0.163921 │
│ Model d.f. │ 1 │
╘════════════════╧═════════════╛
Model Parameters:
╒════════════╤════════════╤════════════╤═════════════╕
│ Variable │ Estimate │ On Bound │ Std Error │
╞════════════╪════════════╪════════════╪═════════════╡
│ g │ 6.80589 │ no │ 0.689269 │
│ v │ 25.6632 │ no │ 2.65373 │
│ k │ 62.7872 │ no │ 3.59334 │
│ n │ 4.87539 │ no │ 1.08357 │
│ alpha │ 15.881 │ no │ 35.6668 │
╘════════════╧════════════╧════════════╧═════════════╛
Goodness of Fit:
╒════════╤═════╤═══════════════╤═════════════════════╤═══════════════════╕
│ Dose │ N │ Sample Mean │ Model Fitted Mean │ Scaled Residual │
╞════════╪═════╪═══════════════╪═════════════════════╪═══════════════════╡
│ 0 │ 20 │ 6 │ 6.80589 │ -0.904381 │
│ 25 │ 20 │ 8 │ 7.09076 │ 1.02037 │
│ 50 │ 20 │ 13 │ 13.1658 │ -0.186064 │
│ 75 │ 20 │ 25 │ 24.8734 │ 0.14202 │
│ 100 │ 20 │ 30 │ 30.0641 │ -0.0719426 │
╘════════╧═════╧═══════════════╧═════════════════════╧═══════════════════╛
╒════════╤═════╤═════════════╤═══════════════════╕
│ Dose │ N │ Sample SD │ Model Fitted SD │
╞════════╪═════╪═════════════╪═══════════════════╡
│ 0 │ 20 │ 4 │ 3.98509 │
│ 25 │ 20 │ 4.3 │ 3.98509 │
│ 50 │ 20 │ 3.8 │ 3.98509 │
│ 75 │ 20 │ 4.4 │ 3.98509 │
│ 100 │ 20 │ 3.7 │ 3.98509 │
╘════════╧═════╧═════════════╧═══════════════════╛
Likelihoods:
╒═════════╤══════════════════╤════════════╤═════════╕
│ Model │ Log-Likelihood │ # Params │ AIC │
╞═════════╪══════════════════╪════════════╪═════════╡
│ A1 │ -279.181 │ 6 │ 570.362 │
│ A2 │ -278.726 │ 10 │ 577.452 │
│ A3 │ -279.181 │ 6 │ 570.362 │
│ fitted │ -280.15 │ 5 │ 570.3 │
│ reduced │ -374.79 │ 2 │ 753.579 │
╘═════════╧══════════════════╧════════════╧═════════╛
Tests of Mean and Variance Fits:
╒════════╤══════════════════════════════╤═════════════╤═══════════╕
│ Name │ -2 * Log(Likelihood Ratio) │ Test d.f. │ P-Value │
╞════════╪══════════════════════════════╪═════════════╪═══════════╡
│ Test 1 │ 192.127 │ 8 │ 0 │
│ Test 2 │ 0.909828 │ 4 │ 0.923147 │
│ Test 3 │ 0.909828 │ 4 │ 0.923147 │
│ Test 4 │ 1.93767 │ 1 │ 0.163921 │
╘════════╧══════════════════════════════╧═════════════╧═══════════╛
Test 1: Test the null hypothesis that responses and variances don't differ among dose levels
(A2 vs R). If this test fails to reject the null hypothesis (p-value > 0.05), there may not be
a dose-response.
Test 2: Test the null hypothesis that variances are homogenous (A1 vs A2). If this test fails to
reject the null hypothesis (p-value > 0.05), the simpler constant variance model may be appropriate.
Test 3: Test the null hypothesis that the variances are adequately modeled (A3 vs A2). If this test
fails to reject the null hypothesis (p-value > 0.05), it may be inferred that the variances have
been modeled appropriately.
Test 4: Test the null hypothesis that the model for the mean fits the data (Fitted vs A3). If this
test fails to reject the null hypothesis (p-value > 0.1), the user has support for use of the
selected model.
You can select a best-fitting model and output reports generated will indicate this selection. This is a manual action. The example below selects the recommended model, but any model in the session could be selected.
session.select(model=session.recommended_model, notes="Lowest AIC; recommended model")