Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Bump BioMASS.jl version to 0.5.0 #11

Merged
merged 7 commits into from
Jul 12, 2021
Merged

Bump BioMASS.jl version to 0.5.0 #11

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
1012e14
Code refactoring
himoto Jun 25, 2021
7150ad3
Update patient_model.py
himoto Jun 29, 2021
dac12a2
Update test for extracting response characteristics
himoto Jul 2, 2021
88e4c89
Specify BioMASS.jl version
himoto Jul 4, 2021
19acc3d
Update README.md
himoto Jul 4, 2021
07c6cb8
Remove save_format
himoto Jul 10, 2021
e103b40
Update README.md
himoto Jul 11, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
19 changes: 12 additions & 7 deletions README.md
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,13 +1,15 @@
# Breast cancer [![Actions Status](https://github.com/pasmopy/breast_cancer/workflows/Tests/badge.svg)](https://github.com/pasmopy/breast_cancer/actions)

Workflow for classifying breast cancer subtypes based on intracellular signaling dynamics.
This repository contains analysis code for the following paper:



## Requirements

| Language | Dependent packages |
| ------------- | -------------------------------------------------------------- |
| Python >= 3.7 | See [`requirements.txt`](requirements.txt) |
| Julia >= 1.5 | [BioMASS.jl](https://github.com/biomass-dev/BioMASS.jl) |
| Julia >= 1.5 | [`BioMASS.jl==0.5.0`](https://github.com/biomass-dev/BioMASS.jl) |
| R >= 4.0 | TCGAbiolinks, sva, biomaRt, edgeR, ComplexHeatmap, viridisLite |

## Table of contents
Expand Down Expand Up @@ -171,11 +173,11 @@ Workflow for classifying breast cancer subtypes based on intracellular signaling

breast_cancer_models = []
path_to_models = os.path.join("models", "breast")
for f in os.listdir(path_to_models):
if os.path.isdir(os.path.join(path_to_models, f)) and (
f.startswith("TCGA_") or f.endswith("_BREAST")
for model in os.listdir(path_to_models):
if os.path.isdir(os.path.join(path_to_models, model)) and (
model.startswith("TCGA_") or model.endswith("_BREAST")
):
breast_cancer_models.append(f)
breast_cancer_models.append(model)
# Set optimized parameters
for model in breast_cancer_models:
shutil.copytree(
Expand Down Expand Up @@ -241,11 +243,14 @@ Workflow for classifying breast cancer subtypes based on intracellular signaling
- Calculate sensitivity coefficients by varying the amount of each nonzero species

```python
import os

from pasmopy import PatientModelAnalyses

import models.breast


with open ("selected_tnbc.txt" mode="r") as f:
with open (os.path.join("models", "breast", "selected_tnbc.txt"), mode="r") as f:
TNBC_ID = f.read().splitlines()
analyses = PatientModelAnalyses(
models.breast.__package__,
Expand Down
216 changes: 88 additions & 128 deletions drug_response/drug/database.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -115,14 +115,60 @@ def _get_drug_response(
)
return drug_response_info

def _check_args(self, drug: str, save_format: str) -> Optional[NoReturn]:
def _check_args(self, drug: str) -> Optional[NoReturn]:
if drug not in set(self.drug_response_data["Compound"]):
raise ValueError(
f"{drug} doesn't exist in the CCLE drug response data.\n"
f"Should be one of {', '.join(list(set(self.drug_response_data['Compound'])))}"
)
if save_format not in ["pdf", "png"]:
raise ValueError("save_format must be either 'pdf' of 'png'.")

@staticmethod
def _plot_dose_response_curve(
x: List[float],
population: List[DrugResponse],
color: str,
label: str,
show_individual: bool,
error :str,
):
if show_individual:
for i, _ in enumerate(population):
plt.plot(
x,
np.interp(x, population[i].doses, population[i].activity_data) + 100,
"o",
color=color,
)
y_mean = np.mean(
[
(np.interp(x, population[i].doses, population[i].activity_data) + 100)
for i, _ in enumerate(population)
],
axis=0
)
y_err = np.std(
[
(np.interp(x, population[i].doses, population[i].activity_data) + 100)
for i, _ in enumerate(population)
],
axis=0,
ddof=1
) / (1 if error == "std" else np.sqrt(len(population)))
plt.plot(
x,
y_mean,
"-",
color=color,
label=label,
)
plt.fill_between(
x,
y_mean - y_err,
y_mean + y_err,
lw=0,
color=color,
alpha=0.1,
)

def save_dose_response_curve(
self,
Expand All @@ -132,17 +178,13 @@ def save_dose_response_curve(
config: Optional[dict] = None,
show_individual: bool = False,
error: str = "std",
save_format: str = "pdf",
) -> None:
self._check_args(drug, save_format)
self._check_args(drug)
if error not in ["std", "sem"]:
raise ValueError("error must be either 'std' or 'sem'.")

os.makedirs(
os.path.join(
"dose_response",
f"{self._drug2target(drug)}",
),
os.path.join("dose_response", f"{self._drug2target(drug)}"),
exist_ok=True
)

Expand All @@ -154,107 +196,36 @@ def save_dose_response_curve(
elif level == "low":
bottom30.append(erbb_expression_ratio.index[i])

egfr_high = self._get_drug_response(
top30,
[drug] * len(top30),
)
egfr_low = self._get_drug_response(
bottom30,
[drug] * len(bottom30),
)
egfr_high = self._get_drug_response(top30, [drug] * len(top30))
egfr_low = self._get_drug_response(bottom30, [drug] * len(bottom30))

if config is None:
config = {}
config.setdefault('figure.figsize', (7, 5))
config.setdefault('font.size', 24)
config.setdefault('font.family', 'Arial')
config.setdefault('axes.linewidth', 2.4)
config.setdefault('xtick.major.width', 2.4)
config.setdefault('ytick.major.width', 2.4)
config.setdefault('lines.linewidth', 3)
config.setdefault("figure.figsize", (7, 5))
config.setdefault("font.size", 24)
config.setdefault("font.family", "Arial")
config.setdefault("axes.linewidth", 2.4)
config.setdefault("xtick.major.width", 2.4)
config.setdefault("ytick.major.width", 2.4)
config.setdefault("lines.linewidth", 3)
config.setdefault("lines.markersize", 1)
config.setdefault("savefig.bbox", "tight")
config.setdefault("savefig.format", "pdf")

plt.rcParams.update(config)

plt.gca().spines["right"].set_visible(False)
plt.gca().spines["top"].set_visible(False)

DOSE = [2.50e-03, 8.00e-03, 2.50e-02, 8.00e-02, 2.50e-01, 8.00e-01, 2.53e+00, 8.00e+00]

if show_individual:
for i, _ in enumerate(egfr_high):
plt.plot(
DOSE,
np.interp(DOSE, egfr_high[i].doses, egfr_high[i].activity_data) + 100,
"o",
color="darkmagenta",
)
y_mean = np.mean(
[
(np.interp(DOSE, egfr_high[i].doses, egfr_high[i].activity_data) + 100)
for i, _ in enumerate(egfr_high)
],
axis=0
self._plot_dose_response_curve(
DOSE, egfr_high, color="darkmagenta", label="EGFR high",
show_individual=show_individual, error=error
)
y_err = np.std(
[
(np.interp(DOSE, egfr_high[i].doses, egfr_high[i].activity_data) + 100)
for i, _ in enumerate(egfr_high)
],
axis=0,
ddof=1
) / (1 if error == "std" else np.sqrt(len(egfr_high)))
plt.plot(
DOSE,
y_mean,
"-",
color="darkmagenta",
label="EGFR high",
)
plt.fill_between(
DOSE,
y_mean - y_err,
y_mean + y_err,
lw=0,
color="darkmagenta",
alpha=0.1,
)

if show_individual:
for i, _ in enumerate(egfr_low):
plt.plot(
DOSE,
np.interp(DOSE, egfr_low[i].doses, egfr_low[i].activity_data) + 100,
"o",
color="goldenrod",
)
y_mean = np.mean(
[
(np.interp(DOSE, egfr_low[i].doses, egfr_low[i].activity_data) + 100)
for i, _ in enumerate(egfr_low)
],
axis=0
)
y_err = np.std(
[
(np.interp(DOSE, egfr_low[i].doses, egfr_low[i].activity_data) + 100)
for i, _ in enumerate(egfr_low)
],
axis=0,
ddof=1
) / (1 if error == "std" else np.sqrt(len(egfr_high)))
plt.plot(
DOSE,
y_mean,
"-",
color="goldenrod",
label="EGFR low",
)
plt.fill_between(
DOSE,
y_mean - y_err,
y_mean + y_err,
lw=0,
color="goldenrod",
alpha=0.1,
self._plot_dose_response_curve(
DOSE, egfr_low, color="goldenrod", label="EGFR low",
show_individual=show_individual, error=error
)

plt.xscale("log")
Expand All @@ -270,10 +241,8 @@ def save_dose_response_curve(
os.path.join(
"dose_response",
f"{self._drug2target(drug)}",
f"{self._convert_drug_name(drug)}.{save_format}",
f"{self._convert_drug_name(drug)}",
),
dpi=1200 if save_format == "png" else None,
bbox_inches="tight",
)
plt.close()

Expand All @@ -284,9 +253,8 @@ def save_activity_area(
drug: str,
*,
config: Optional[dict] = None,
save_format: str = "pdf",
) -> None:
self._check_args(drug, save_format)
self._check_args(drug)

os.makedirs(
os.path.join(
Expand All @@ -304,33 +272,30 @@ def save_activity_area(
elif level == "low":
bottom30.append(erbb_expression_ratio.index[i])

egfr_high = self._get_drug_response(
top30,
[drug]*len(top30),
)
egfr_low = self._get_drug_response(
bottom30,
[drug]*len(bottom30),
)

egfr_high = self._get_drug_response(top30, [drug] * len(top30))
egfr_low = self._get_drug_response(bottom30, [drug] * len(bottom30))

activity_area = np.array(
[
[egfr_high[i].act_area for i, _ in enumerate(egfr_high)],
[egfr_low[i].act_area for i, _ in enumerate(egfr_low)],
],
dtype=object,
)

if config is None:
config = {}
config.setdefault('figure.figsize', (3, 5))
config.setdefault('font.family', 'Arial')
config.setdefault('mathtext.it', 'Arial:italic')
config.setdefault('font.size', 24)
config.setdefault('axes.linewidth', 2.4)
config.setdefault('xtick.major.width', 2.4)
config.setdefault('ytick.major.width', 2.4)
config.setdefault('lines.linewidth', 1.8)
config.setdefault('lines.markersize', 11)
config.setdefault("figure.figsize", (3, 5))
config.setdefault("font.family", "Arial")
config.setdefault("mathtext.it", "Arial:italic")
config.setdefault("font.size", 24)
config.setdefault("axes.linewidth", 2.4)
config.setdefault("xtick.major.width", 2.4)
config.setdefault("ytick.major.width", 2.4)
config.setdefault("lines.linewidth", 1.8)
config.setdefault("lines.markersize", 11)
config.setdefault("savefig.bbox", "tight")
config.setdefault("savefig.format", "pdf")

plt.rcParams.update(config)

Expand All @@ -357,10 +322,8 @@ def save_activity_area(
os.path.join(
"activity_area",
f"{self._drug2target(drug)}",
f"{self._convert_drug_name(drug)}.{save_format}",
f"{self._convert_drug_name(drug)}",
),
dpi=1200 if save_format == "png" else None,
bbox_inches="tight",
)
plt.close()

Expand All @@ -371,17 +334,14 @@ def save_all(
*,
config_dose_response_curve: Optional[dict] = None,
config_activity_area: Optional[dict] = None,
save_format: str = "pdf",
) -> None:
self.save_dose_response_curve(
erbb_expression_ratio,
drug,
config=config_dose_response_curve,
save_format=save_format,
)
self.save_activity_area(
erbb_expression_ratio,
drug,
config=config_activity_area,
save_format=save_format,
)
Loading