From 8608d41eca69dc9508f8ae9a78e791ac2c9279f6 Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 5 Nov 2025 11:25:50 -0800 Subject: [PATCH 01/15] add plot_MIP_with_mask_outlines --- pytheranostics/misc_tools/tools.py | 49 ++++++++++++++++++++++++++++++ 1 file changed, 49 insertions(+) diff --git a/pytheranostics/misc_tools/tools.py b/pytheranostics/misc_tools/tools.py index 44b447a..c355428 100644 --- a/pytheranostics/misc_tools/tools.py +++ b/pytheranostics/misc_tools/tools.py @@ -389,3 +389,52 @@ def plot_MIP(SPECT=None, vmax=300000, figsize=(10, 5), ax=None): plt.yticks([]) return fig, ax + + +def plot_MIP_with_mask_outlines(ax, SPECT, masks=None, vmax=300000): + """Plot Maximum Intensity Projection (MIP) of SPECT data with masks outlines. + + Parameters + ---------- + ax : _type_ + _description_ + SPECT : _type_ + _description_ + masks : _type_, optional + _description_, by default None + vmax : int, optional + _description_, by default 300000 + """ + plt.sca(ax) + spect_mip = SPECT.max(axis=0) + plt.imshow(spect_mip.T, cmap="Greys", interpolation="Gaussian", vmax=vmax, vmin=0) + + if masks is not None: + for organ, mask in masks.items(): + organ_lower = organ.lower() + if "kidney" in organ_lower: + color = "lime" + elif "gland" in organ_lower: + color = "red" + elif "lesion" in organ_lower: + color = "m" + else: + continue + + mip_mask = mask.max(axis=0) + if mip_mask.shape != spect_mip.shape: + mip_mask = mip_mask.T + + plt.contour( + numpy.transpose(mip_mask, (1, 0)), + levels=[0.5], + colors=[color], + linewidths=1.5, + alpha=0.5, + ) + + plt.xlim(30, 100) + plt.ylim(0, 234) + plt.axis("off") + plt.xticks([]) + plt.yticks([]) From 1aaa6d4e35e0d37cdde31218e45f50b6255d9ccf Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 5 Nov 2025 11:32:52 -0800 Subject: [PATCH 02/15] improve formatting of variables --- pytheranostics/dosimetry/organ_s_dosimetry.py | 130 +++++++++++------- 1 file changed, 77 insertions(+), 53 deletions(-) diff --git a/pytheranostics/dosimetry/organ_s_dosimetry.py b/pytheranostics/dosimetry/organ_s_dosimetry.py index 0670872..4088f7c 100644 --- a/pytheranostics/dosimetry/organ_s_dosimetry.py +++ b/pytheranostics/dosimetry/organ_s_dosimetry.py @@ -117,22 +117,27 @@ def apply_sphere_method(self, df: pandas.DataFrame) -> pandas.DataFrame: def calculate_ttb(self): """Compute Total Tumor Burden (TTB) metrics and append to results_lesions.""" metrics = { - "Mass_g": self.results_lesions["Mass_g"].sum(), - "Volume_CT_mL": self.results_lesions["Volume_CT_mL"].sum(), - "TIA_h": self.results_lesions["TIA_h"].sum(), + "Mass_g": self.results_dosimetry_lesions["Mass_g"].sum(), + "Volume_CT_mL": self.results_dosimetry_lesions["Volume_CT_mL"].sum(), + "TIA_h": self.results_dosimetry_lesions["TIA_h"].sum(), "AD_Gy": ( - (self.results_lesions["Mass_g"] * self.results_lesions["AD_Gy"]).sum() + ( + self.results_dosimetry_lesions["Mass_g"] + * self.results_dosimetry_lesions["AD_Gy"] + ).sum() ) / ( - self.results_lesions["Mass_g"].sum() - if self.results_lesions["Mass_g"].sum() > 0 + self.results_dosimetry_lesions["Mass_g"].sum() + if self.results_dosimetry_lesions["Mass_g"].sum() > 0 else 0 ), } TTB = pandas.DataFrame(metrics, index=["TTB"]) - self.results_lesions = pandas.concat([self.results_lesions, TTB], axis=0) - return self.results_lesions + self.results_dosimetry_lesions = pandas.concat( + [self.results_dosimetry_lesions, TTB], axis=0 + ) + return self.results_dosimetry_lesions def prepare_data(self) -> None: """ @@ -223,6 +228,7 @@ def prepare_data(self) -> None: self.results_fitting.loc["Red Marrow"][ "Volume_CT_mL" ] = 1170 # TODO volume hardcoded, think about alternatives + self.results_fitting.loc["RemainderOfBody"]["Volume_CT_mL"] = ( self.config["PatientWeight_g"] - self.results_fitting.loc[ @@ -241,9 +247,9 @@ def prepare_data(self) -> None: self.results_fitting_organs = self.results_fitting[ ~lesion_mask ].copy() # all non-lesion entries - self.results_fitting_lesions = self.results_fitting[ - lesion_mask - ].copy() # only lesion entries + self.results_fitting_lesions = self.results[ + self.results.index.str.contains("Lesion") + ] if "TotalTumorBurden" in self.results_fitting.index: self.results_fitting.drop("TotalTumorBurden", axis=0, inplace=True) @@ -251,13 +257,13 @@ def prepare_data(self) -> None: if output_type == "Export": fmt = organ_conf["Output"]["ExportFormat"] if fmt.lower() == "olinda": - self.results_fitting = self.results_fitting.rename( + self.results_fitting_organs = self.results_fitting_organs.rename( index={"RemainderOfBody": "Total Body"} ) - self.results_fitting.loc["Total Body"]["Volume_CT_mL"] = ( + self.results_fitting_organs.loc["Total Body"]["Volume_CT_mL"] = ( self.config["PatientWeight_g"] - - self.results_fitting.loc[ - ~self.results_fitting.index.isin( + - self.results_fitting_organs.loc[ + ~self.results_fitting_organs.index.isin( ["Total Body", "RemainderOfBody"] ), "Volume_CT_mL", @@ -327,7 +333,7 @@ def process_dosimetry(self) -> None: "LesionDosimetry" ): self.results_dosimetry_lesions = self.apply_sphere_method( - self.results_fitting_lesions.index.str.contains("Lesion") + self.results_fitting_lesions ) self.results_dosimetry_lesions = self.calculate_ttb() if "Yes" in self.config["OrganLevel"]["AdditionalOptions"].get( @@ -342,12 +348,12 @@ def calculate_absorbed_dose(self) -> pandas.DataFrame: """Calculate absorbed dose per target organ based on model and disintegration data.""" model_files = { "Female": { - "beta": f'177Lu_S_values_female_{self.config["OrganLevel"]["Calculation"]["SValueSource"].lower()}_BETA.csv', - "gamma": f'177Lu_S_values_female_{self.config["OrganLevel"]["Calculation"]["SValueSource"].lower()}_GAMMA.csv', + "gamma": f'{self.config["Radionuclide"]}_S_values_female_{self.config["OrganLevel"]["Calculation"]["SValueSource"].lower()}_GAMMA.csv', + "beta": f'{self.config["Radionuclide"]}_S_values_female_{self.config["OrganLevel"]["Calculation"]["SValueSource"].lower()}_BETA.csv', }, "Male": { - "beta": f'177Lu_S_values_male_{self.config["OrganLevel"]["Calculation"]["SValueSource"].lower()}_BETA.csv', - "gamma": f'177Lu_S_values_male_{self.config["OrganLevel"]["Calculation"]["SValueSource"].lower()}_GAMMA.csv', + "gamma": f'{self.config["Radionuclide"]}_S_values_male_{self.config["OrganLevel"]["Calculation"]["SValueSource"].lower()}_GAMMA.csv', + "beta": f'{self.config["Radionuclide"]}_S_values_male_{self.config["OrganLevel"]["Calculation"]["SValueSource"].lower()}_BETA.csv', }, } @@ -359,7 +365,7 @@ def calculate_absorbed_dose(self) -> pandas.DataFrame: ) print("Source organs available in the model:", svalues_beta.columns.tolist()) - print("Source organs present :", self.results_fitting.index.tolist()) + print("Source organs present :", self.results_fitting_organs.index.tolist()) self.source_organs_missing = set(svalues_beta.columns) - set( self.results_fitting.index @@ -391,7 +397,7 @@ def calculate_absorbed_dose(self) -> pandas.DataFrame: ) # Apply mass scaling - dose_df = self.perform_mass_scaling(dose_df, self.config["Gender"]) + dose_df = self.perform_mass_scaling(dose_df) # Calculate absorbed dose in Gy for injected activity dose_df["AD_total[Gy/GBq]"] = ( @@ -403,7 +409,10 @@ def calculate_absorbed_dose(self) -> pandas.DataFrame: dose_df["AD_total[Gy]"] = dose_df["AD_total[Gy/GBq]"] / 1000 * injected_activity dose_df = dose_df.reset_index(drop=True) - self.df_ad = dose_df.copy() + self.results_dosimetry_organs = dose_df.copy() + self.results_dosimetry_organs = self.results_dosimetry_organs.set_index( + "Target organ" + ) return dose_df def hollow_organ_correction(self, df: pandas.DataFrame) -> pandas.DataFrame: @@ -476,9 +485,12 @@ def redistribute_ROB_into_source_organs_missing( def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: """Multiply S-values by TIA to compute dose matrix for radiation type.""" # Path to organ masses - masses_path = path.join(MASSES_PATH, "ICRP_mass_male.csv") + masses_path = path.join( + MASSES_PATH, f"ICRP_mass_{self.config['Gender'].lower()}.csv" + ) self.organ_masses = pandas.read_csv(masses_path, index_col=0) + print(f"Applying S-values for {radiation_type} radiation...") # Handle remainder of the body # Redistribute ROB TIA into missing source organs if needed - approach consistent with MIRDcalc software if "RemainderOfBody" in tia_df.index: @@ -545,38 +557,47 @@ def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: dose_df = s_values_subset.multiply(tia_series, axis=1) # ROB + + print("Calculating Remainder of Body contributions...") dose_df["RemainderOfBody"] = 0 - total_body_mass = self.config["PatientWeight_g"] + total_body_mass = self.organ_masses.loc["Total Body", "Mass_g"] ROB_mass = tia_df.loc["RemainderOfBody", "Volume_CT_mL"] + if radiation_type == "beta": organs = s_values_subset.index.difference( tia_series.index.difference(["Red Marrow", "Osteogenic Cells"]) ) if radiation_type == "gamma": organs = s_values_subset.index - # adjust source organs so that they are different for gamma and beta radiation (in beta it is total - source organs + bone + skeleton) - + self.contribution_from_sources_dict = {radiation_type: {}} for target_organ in organs: contribution_from_sources = 0 + + self.contribution_from_sources_dict[radiation_type][ + target_organ + ] = {"sources": {}} + for source_organ in s_values.columns.difference(tia_series.index): + if radiation_type == "beta": + if target_organ == "Red Marrow" and source_organ in [ + "Trabecular Bone" + ]: + continue + if target_organ == "Osteogenic Cells" and source_organ in [ + "Cortical Bone", + "Trabecular Bone", + "Red Marrow", + ]: + continue - if target_organ.split()[0] == source_organ.split()[0]: - continue - if target_organ == "Osteogenic Cells" and source_organ in [ - "Cortical Bone", - "Trabecular Bone", - "Red Marrow", - ]: + if source_organ == "Total Body": continue - if target_organ == "Red Marrow" and source_organ in [ - "Trabecular Bone" - ]: + if target_organ.split()[0] == source_organ.split()[0]: continue if target_organ == "Total Body": continue - if source_organ == "Total Body": - continue + print(f" From source organ: {source_organ}") source_organ_mass = self.organ_masses.loc[ source_organ, "Mass_g" ] @@ -587,6 +608,11 @@ def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: contribution_from_sources += s_value_source_to_target * ( source_organ_mass / ROB_mass ) + self.contribution_from_sources_dict[radiation_type][ + target_organ + ]["sources"][ + f"contribution_from_{source_organ}" + ] = s_value_source_to_target s_value_ROB_to_target = ( s_values.loc[target_organ, "Total Body"] @@ -603,10 +629,13 @@ def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: return dose_df def perform_mass_scaling( - self, df: pandas.DataFrame, gender: str + self, + df: pandas.DataFrame, ) -> pandas.DataFrame: """Apply mass scaling to absorbed dose calculations based on patient-specific organ masses.""" - masses_path = path.join(MASSES_PATH, f"ICRP_mass_{gender.lower()}_target.csv") + masses_path = path.join( + MASSES_PATH, f"ICRP_mass_{self.config['Gender'].lower()}_target.csv" + ) model_masses_df = pandas.read_csv(masses_path, index_col=0) print("Performing mass scaling...") @@ -649,14 +678,9 @@ def create_Olinda_file(self, dirname: str, savefile: bool = False) -> None: this_dir = path.dirname(__file__) TEMPLATE_PATH = path.join(this_dir, "olindaTemplates") - if self.config["Gender"] == "Male": - template = pandas.read_csv(path.join(TEMPLATE_PATH, "adult_male.cas")) - elif self.config["Gender"] == "Female": - template = pandas.read_csv(path.join(TEMPLATE_PATH, "adult_female.cas")) - else: - print( - "Ensure that you correctly wrote patient gender in config file. Olinda supports: Male and Female." - ) + template = pandas.read_csv( + path.join(TEMPLATE_PATH, f"adult_{self.config['Gender'].lower()}.cas") + ) template.columns = ["Data"] match = re.match(r"([a-zA-Z]+)([0-9]+)", self.config["Radionuclide"]) @@ -666,13 +690,13 @@ def create_Olinda_file(self, dirname: str, savefile: bool = False) -> None: ind = template[template["Data"] == "[BEGIN NUCLIDES]"].index template.loc[ind[0] + 1, "Data"] = formatted_radionuclide + "|" - for organ in self.results_fitting.index: + for organ in self.results_fitting_organs.index: indices = template[template["Data"].str.contains(organ)].index source_organ = template.iloc[indices[0]].str.split("|")[0][0] mass_phantom = template.iloc[indices[0]].str.split("|")[0][1] - kinetic_data = self.results_fitting.loc[organ]["TIA_h"] - mass_data = round(self.results_fitting.loc[organ]["Volume_CT_mL"], 1) + kinetic_data = self.results_fitting_organs.loc[organ]["TIA_h"] + mass_data = round(self.results_fitting_organs.loc[organ]["Volume_CT_mL"], 1) # Update the template DataFrame template.iloc[indices[0]] = ( @@ -755,7 +779,7 @@ def read_Olinda_results(self, olinda_results_path: str) -> None: df_ad["Total"].fillna(0, inplace=True) df_ad["AD[Gy]"] = df_ad["Total"] * float(self.config["InjectedActivity"]) / 1000 df_ad = df_ad.rename(columns={"Total": "AD[Gy/GBq]"}) - self.df_ad = df_ad + self.results_dosimetry_organs = df_ad def compute_dose(self): """Compute Time Integrated Activity.""" From 4623ba7b476cfa22ff88982f9f6de26b67220c43 Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 5 Nov 2025 11:42:28 -0800 Subject: [PATCH 03/15] rename variables for consistency --- pytheranostics/dosimetry/base_dosimetry.py | 153 ++++++++++----------- 1 file changed, 73 insertions(+), 80 deletions(-) diff --git a/pytheranostics/dosimetry/base_dosimetry.py b/pytheranostics/dosimetry/base_dosimetry.py index 3380a44..2691c8c 100644 --- a/pytheranostics/dosimetry/base_dosimetry.py +++ b/pytheranostics/dosimetry/base_dosimetry.py @@ -122,9 +122,9 @@ def __init__( # DataFrame storing results self.results = self.initialize() - self.results_lesions = pandas.DataFrame() - self.results_salivaryglands = pandas.DataFrame() - self.df_ad = pandas.DataFrame() + self.results_dosimetry_lesions = pandas.DataFrame() + self.results_dosimetry_salivaryglands = pandas.DataFrame() + self.results_dosimetry_organs = pandas.DataFrame() # Sanity Checks: self.sanity_checks(metric="Volume_CT_mL") @@ -459,23 +459,7 @@ def compute_tia(self) -> None: def smart_fit_selection( self, region_data: pandas.Series, region: str ) -> lmfit.model.ModelResult: - """Select the best fit based on Akaike Information Criterion. - - If enabled in self.config, iterates through different valid fits orders and select best fit based on Akaike Information Criterion. - If a fit order is specified by the user, then the method will just perform fit following user's selected order and configuration. - - Parameters - ---------- - region_data : pandas.Series - Series containing Time and Activity. - region : str - Region of Interest - - Returns - ------- - lmfit.model.ModelResult - The best fit model based on Akaike Information Criterion. - """ + """Select the best fit based on Akaike Information Criterion.""" # If fit_order is defined by user: if self.config["rois"][region]["fit_order"] is not None: fit_results, _ = exponential_fit_lmfit( @@ -635,8 +619,8 @@ def calculate_bed(self, kinetic: str) -> None: RADIOBIOLOGY_DATA_FILE = Path(this_dir, "data", "radiobiology.json") with open(RADIOBIOLOGY_DATA_FILE) as f: self.radiobiology_dic = json.load(f) - bed_df = self.df_ad[ - self.df_ad.index.isin(list(self.radiobiology_dic.keys())) + bed_df = self.results_dosimetry_organs[ + self.results_dosimetry_organs.index.isin(list(self.radiobiology_dic.keys())) ] # only organs that we know the radiobiology parameters organs = numpy.array(bed_df.index.unique()) bed = {} @@ -645,7 +629,7 @@ def calculate_bed(self, kinetic: str) -> None: t_repair = self.radiobiology_dic[organ]["t_repair"] alpha_beta = self.radiobiology_dic[organ]["alpha_beta"] AD = ( - float(self.df_ad.loc[bed_df.index == organ]["AD[Gy/GBq]"].values[0]) + float(bed_df.loc[bed_df.index == organ]["AD_total[Gy/GBq]"].values[0]) * float(self.config["InjectedActivity"]) / 1000 ) # Gy @@ -701,7 +685,9 @@ def calculate_bed(self, kinetic: str) -> None: ) print(f"{organ}", bed[organ]) - self.df_ad["BED[Gy]"] = self.df_ad.index.map(bed) + self.results_dosimetry_organs["BED[Gy]"] = ( + self.results_dosimetry_organs.index.map(bed) + ) def save_images_and_masks_at(self, time_id: int) -> None: """Save CT, NM and masks for a specific time point. @@ -773,11 +759,13 @@ def write_json_data( cycle["InjectedActivity"] = self.config["InjectedActivity"] cycle["Weight_g"] = self.config["PatientWeight_g"] cycle["Level"] = self.config["Level"] - cycle["Method"] = self.config["Method"] - cycle["OutputFormat"] = self.config["OutputFormat"] - cycle["ScaleDoseByDensity"] = self.config.get( - "ScaleDoseByDensity", cycle.get("ScaleDoseByDensity", "NA") - ) + if cycle["Level"] == "Organ": + cycle["Method"] = self.config["OrganLevel"] + elif cycle["Level"] == "Voxel": + cycle["Method"] = self.config["VoxelLevel"] + cycle["ScaleDoseByDensity"] = self.config.get( + "ScaleDoseByDensity", cycle.get("ScaleDoseByDensity", "NA") + ) cycle["ReferenceTimePoint"] = self.config["ReferenceTimePoint"] cycle["TimePoints_h"] = self.results["Time_hr"][0] @@ -884,24 +872,24 @@ def write_json_data( if "Lesion" in organ or "TTB" in organ: cycle["rois"][organ]["density_gml"]["different_tps"] = "NA" cycle["rois"][organ]["density_gml"]["uncertainty"] = "NA" - cycle["rois"][organ]["density_gml"]["mean"] = self.results_lesions.loc[ - organ, "Density_g_per_mL" - ] + cycle["rois"][organ]["density_gml"]["mean"] = ( + self.results_dosimetry_lesions.loc[organ, "Density_g_per_mL"] + ) cycle["rois"][organ]["density_gml"]["mean_uncertainty"] = "NA" cycle["rois"][organ]["mass_g"]["different_tps"] = "NA" cycle["rois"][organ]["mass_g"]["uncertainty"] = "NA" - cycle["rois"][organ]["mass_g"]["mean"] = self.results_lesions.loc[ - organ, "Mass_g" - ] + cycle["rois"][organ]["mass_g"]["mean"] = ( + self.results_dosimetry_lesions.loc[organ, "Mass_g"] + ) cycle["rois"][organ]["mass_g"]["mean_uncertainty"] = "NA" - cycle["rois"][organ]["composition"] = self.results_lesions.loc[ - organ, "Composition" - ] - cycle["rois"][organ]["total_s_value"] = self.results_lesions.loc[ - organ, "Total_S_Value" - ] + cycle["rois"][organ]["composition"] = ( + self.results_dosimetry_lesions.loc[organ, "Composition"] + ) + cycle["rois"][organ]["total_s_value"] = ( + self.results_dosimetry_lesions.loc[organ, "Total_S_Value"] + ) cycle["rois"][organ]["total_s_value_uncertainty"] = "NA" - cycle["rois"][organ]["mean_AD_Gy"] = self.results_lesions.loc[ + cycle["rois"][organ]["mean_AD_Gy"] = self.results_dosimetry_lesions.loc[ organ, "AD_Gy" ] cycle["rois"][organ]["mean_AD_Gy_uncertainty"] = "NA" @@ -915,30 +903,30 @@ def write_json_data( cycle["rois"][organ]["density_gml"]["different_tps"] = "NA" cycle["rois"][organ]["density_gml"]["uncertainty"] = "NA" cycle["rois"][organ]["density_gml"]["mean"] = ( - self.results_salivaryglands.loc[organ, "Density_g_per_mL"] + self.results_dosimetry_salivaryglands.loc[organ, "Density_g_per_mL"] ) cycle["rois"][organ]["density_gml"]["mean_uncertainty"] = "NA" cycle["rois"][organ]["mass_g"]["different_tps"] = "NA" cycle["rois"][organ]["mass_g"]["uncertainty"] = "NA" cycle["rois"][organ]["mass_g"]["mean"] = ( - self.results_salivaryglands.loc[organ, "Mass_g"] + self.results_dosimetry_salivaryglands.loc[organ, "Mass_g"] ) cycle["rois"][organ]["mass_g"]["mean_uncertainty"] = "NA" - cycle["rois"][organ]["composition"] = self.results_salivaryglands.loc[ - organ, "Composition" - ] - cycle["rois"][organ]["total_s_value"] = self.results_salivaryglands.loc[ - organ, "Total_S_Value" - ] + cycle["rois"][organ]["composition"] = ( + self.results_dosimetry_salivaryglands.loc[organ, "Composition"] + ) + cycle["rois"][organ]["total_s_value"] = ( + self.results_dosimetry_salivaryglands.loc[organ, "Total_S_Value"] + ) cycle["rois"][organ]["total_s_value_uncertainty"] = "NA" - cycle["rois"][organ]["mean_AD_Gy"] = self.results_salivaryglands.loc[ - organ, "AD_Gy" - ] + cycle["rois"][organ]["mean_AD_Gy"] = ( + self.results_dosimetry_salivaryglands.loc[organ, "AD_Gy"] + ) cycle["rois"][organ]["mean_AD_Gy_uncertainty"] = "NA" if self.config["Level"] == "Organ": - for organ in self.df_ad.index: - if organ in self.df_ad.index: + for organ in self.results_dosimetry_organs.index: + if organ in self.results_dosimetry_organs.index: cycle["Organ-level_AD"][organ] = { "AD[Gy/GBq]": {}, "AD[Gy/GBq]_uncertianty": {}, @@ -947,19 +935,21 @@ def write_json_data( "BED[Gy]": {}, "BED[Gy]_uncertianty": {}, } - cycle["Organ-level_AD"][organ]["AD[Gy/GBq]"] = self.df_ad.loc[ - organ, "AD[Gy/GBq]" - ] + cycle["Organ-level_AD"][organ]["AD[Gy/GBq]"] = ( + self.results_dosimetry_organs.loc[organ, "AD_total[Gy/GBq]"] + ) cycle["Organ-level_AD"][organ]["AD[Gy/GBq]_uncertainty"] = "NA" - cycle["Organ-level_AD"][organ]["AD[Gy]"] = self.df_ad.loc[ - organ, "AD[Gy]" - ] + cycle["Organ-level_AD"][organ]["AD[Gy]"] = ( + self.results_dosimetry_organs.loc[organ, "AD_total[Gy]"] + ) cycle["Organ-level_AD"][organ]["AD[Gy]_uncertainty"] = "NA" - if "BED[Gy]" in self.df_ad.columns: + if "BED[Gy]" in self.results_dosimetry_organs.columns: cycle["Organ-level_AD"][organ]["BED[Gy]"] = ( - self.df_ad.loc[organ, "BED[Gy]"] - if pandas.notna(self.df_ad.loc[organ, "BED[Gy]"]) + self.results_dosimetry_organs.loc[organ, "BED[Gy]"] + if pandas.notna( + self.results_dosimetry_organs.loc[organ, "BED[Gy]"] + ) else "NA" ) else: @@ -967,7 +957,9 @@ def write_json_data( cycle["Organ-level_AD"][organ]["BED[Gy]_uncertianty"] = "NA" - if "Yes" in self.config["LesionDosimetry"]: + if "Yes" in self.config["OrganLevel"]["AdditionalOptions"].get( + "LesionDosimetry" + ): cycle["Organ-level_AD"]["TTB"] = { "mass_g": {}, "volumes_mL": {}, @@ -977,22 +969,23 @@ def write_json_data( "AD[Gy/GBq]": {}, "AD[Gy/GBq]_uncertianty": {}, } - cycle["Organ-level_AD"]["TTB"]["mass_g"] = self.results_lesions.loc[ - "TTB", "Mass_g" - ] - cycle["Organ-level_AD"]["TTB"]["volumes_mL"] = self.results_lesions.loc[ - "TTB", "Volume_CT_mL" - ] - cycle["Organ-level_AD"]["TTB"]["TIA_h"] = self.results_lesions.loc[ - "TTB", "TIA_h" - ] - cycle["Organ-level_AD"]["TTB"]["AD[Gy]"] = self.results_lesions.loc[ - "TTB", "AD_Gy" - ] + cycle["Organ-level_AD"]["TTB"]["mass_g"] = ( + self.results_dosimetry_lesions.loc["TTB", "Mass_g"] + ) + cycle["Organ-level_AD"]["TTB"]["volumes_mL"] = ( + self.results_dosimetry_lesions.loc["TTB", "Volume_CT_mL"] + ) + cycle["Organ-level_AD"]["TTB"]["TIA_h"] = ( + self.results_dosimetry_lesions.loc["TTB", "TIA_h"] + ) + cycle["Organ-level_AD"]["TTB"]["AD[Gy]"] = ( + self.results_dosimetry_lesions.loc["TTB", "AD_Gy"] + ) cycle["Organ-level_AD"]["TTB"]["AD[Gy]_uncertainty"] = "NA" - cycle["Organ-level_AD"]["TTB"]["AD[Gy/GBq]"] = self.results_lesions.loc[ - "TTB", "AD_Gy" - ] / (float(self.config["InjectedActivity"]) / 1000) + cycle["Organ-level_AD"]["TTB"]["AD[Gy/GBq]"] = ( + self.results_dosimetry_lesions.loc["TTB", "AD_Gy"] + / (float(self.config["InjectedActivity"]) / 1000) + ) cycle["Organ-level_AD"]["TTB"]["AD[Gy/GBq]_uncertianty"] = "NA" with open(file_path, "w") as file: From 6604d1b3fa1f3c2414c11271190a59eb1dab4a2e Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 5 Nov 2025 12:41:34 -0800 Subject: [PATCH 04/15] add condition for when there are no lesions --- pytheranostics/dosimetry/organ_s_dosimetry.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/pytheranostics/dosimetry/organ_s_dosimetry.py b/pytheranostics/dosimetry/organ_s_dosimetry.py index 4088f7c..fc2bc37 100644 --- a/pytheranostics/dosimetry/organ_s_dosimetry.py +++ b/pytheranostics/dosimetry/organ_s_dosimetry.py @@ -250,6 +250,10 @@ def prepare_data(self) -> None: self.results_fitting_lesions = self.results[ self.results.index.str.contains("Lesion") ] + elif "No" in self.config["OrganLevel"]["AdditionalOptions"].get( + "LesionDosimetry" + ): + self.results_fitting_organs = self.results_fitting.copy() if "TotalTumorBurden" in self.results_fitting.index: self.results_fitting.drop("TotalTumorBurden", axis=0, inplace=True) From aac740356cbad1ad4fcdad3e43850beaa3b4c742 Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Thu, 13 Nov 2025 13:25:39 -0800 Subject: [PATCH 05/15] improve mass scaling for total body --- pytheranostics/dosimetry/organ_s_dosimetry.py | 64 ++++++++++--------- 1 file changed, 35 insertions(+), 29 deletions(-) diff --git a/pytheranostics/dosimetry/organ_s_dosimetry.py b/pytheranostics/dosimetry/organ_s_dosimetry.py index fc2bc37..b3515e1 100644 --- a/pytheranostics/dosimetry/organ_s_dosimetry.py +++ b/pytheranostics/dosimetry/organ_s_dosimetry.py @@ -561,28 +561,30 @@ def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: dose_df = s_values_subset.multiply(tia_series, axis=1) # ROB - print("Calculating Remainder of Body contributions...") dose_df["RemainderOfBody"] = 0 total_body_mass = self.organ_masses.loc["Total Body", "Mass_g"] - ROB_mass = tia_df.loc["RemainderOfBody", "Volume_CT_mL"] - + source_organs = tia_series.index if radiation_type == "beta": - organs = s_values_subset.index.difference( + target_organs = s_values_subset.index.difference( tia_series.index.difference(["Red Marrow", "Osteogenic Cells"]) ) + ROB_mass = tia_df.loc["RemainderOfBody", "Volume_CT_mL"] if radiation_type == "gamma": - organs = s_values_subset.index - self.contribution_from_sources_dict = {radiation_type: {}} - for target_organ in organs: - contribution_from_sources = 0 + target_organs = s_values_subset.index + total_mass_source_organs = self.organ_masses.loc[ + self.organ_masses.index.intersection(tia_series.index), "Mass_g" + ].sum() - self.contribution_from_sources_dict[radiation_type][ - target_organ - ] = {"sources": {}} + ROB_mass = total_body_mass - total_mass_source_organs - for source_organ in s_values.columns.difference(tia_series.index): + for target_organ in target_organs: + contribution_from_sources = 0 + + for source_organ in source_organs: if radiation_type == "beta": + if target_organ.split()[0] == source_organ.split()[0]: + continue if target_organ == "Red Marrow" and source_organ in [ "Trabecular Bone" ]: @@ -596,12 +598,10 @@ def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: if source_organ == "Total Body": continue - if target_organ.split()[0] == source_organ.split()[0]: - continue + if target_organ == "Total Body": continue - print(f" From source organ: {source_organ}") source_organ_mass = self.organ_masses.loc[ source_organ, "Mass_g" ] @@ -612,16 +612,18 @@ def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: contribution_from_sources += s_value_source_to_target * ( source_organ_mass / ROB_mass ) - self.contribution_from_sources_dict[radiation_type][ - target_organ - ]["sources"][ - f"contribution_from_{source_organ}" - ] = s_value_source_to_target - - s_value_ROB_to_target = ( - s_values.loc[target_organ, "Total Body"] - * (total_body_mass / ROB_mass) - ) - contribution_from_sources + + if target_organ == "Total Body": + s_value_ROB_to_target = ( + s_values.loc[target_organ, "Total Body"] + # * (total_body_mass / ROB_mass) since it is a target organ, the scaling will happen at mass scaling method, so not to have it twice its not performed here + ) - contribution_from_sources + else: + s_value_ROB_to_target = ( + s_values.loc[target_organ, "Total Body"] + * (total_body_mass / ROB_mass) + ) - contribution_from_sources + dose_df.at[target_organ, "RemainderOfBody"] = ( s_value_ROB_to_target * tia_df.loc["RemainderOfBody", "TIA_s"] ) @@ -645,10 +647,14 @@ def perform_mass_scaling( print("Performing mass scaling...") for organ in df["Target organ"]: - - if organ in model_masses_df.index and organ in self.results_fitting.index: - model_mass = model_masses_df.loc[organ, "Mass_g"] - patient_mass = self.results_fitting.loc[organ, "Volume_CT_mL"] + # Temporary mapping just for internal calculations + organ_internal = "RemainderOfBody" if organ == "Total Body" else organ + if ( + organ_internal in model_masses_df.index + and organ_internal in self.results_fitting.index + ): + model_mass = model_masses_df.loc[organ_internal, "Mass_g"] + patient_mass = self.results_fitting.loc[organ_internal, "Volume_CT_mL"] if ( pandas.notna(patient_mass) From 5d2e80dc2365406f6cf340d9f7d7e5eb7f62e8b1 Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 19 Nov 2025 12:24:51 -0800 Subject: [PATCH 06/15] get patient height --- pytheranostics/imaging_ds/cycle_loader.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/pytheranostics/imaging_ds/cycle_loader.py b/pytheranostics/imaging_ds/cycle_loader.py index b3ebee8..90d8a36 100644 --- a/pytheranostics/imaging_ds/cycle_loader.py +++ b/pytheranostics/imaging_ds/cycle_loader.py @@ -171,12 +171,19 @@ def extract_injection_from_first_tp_spect( weight_g = int(round(float(pw) * 1000.0)) except Exception: pass - + height_cm: Optional[int] = None + pw = getattr(ds, "PatientSize", None) # cm + if pw is not None: + try: + height_cm = int(round(float(pw) * 100.0)) + except Exception: + pass return { "InjectionDate": inj_date or "", "InjectionTime": inj_time or "", "InjectedActivity": injected_activity, # Bq (int) or None "PatientWeight_g": weight_g, + "PatientHeight_cm": height_cm, } From c92c374fb2e1cfc791561f9018ba01693a3f72d2 Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 19 Nov 2025 12:26:19 -0800 Subject: [PATCH 07/15] rename from rois to VOIs --- pytheranostics/dosimetry/base_dosimetry.py | 199 ++++++++++----------- 1 file changed, 99 insertions(+), 100 deletions(-) diff --git a/pytheranostics/dosimetry/base_dosimetry.py b/pytheranostics/dosimetry/base_dosimetry.py index 4358836..c7e19d5 100644 --- a/pytheranostics/dosimetry/base_dosimetry.py +++ b/pytheranostics/dosimetry/base_dosimetry.py @@ -155,8 +155,8 @@ def default_config(self) -> None: for key, value in defaults.items(): for region, _ in self.results.iterrows(): - if key not in self.config["rois"][region]: - self.config["rois"][region][key] = value + if key not in self.config["VOIs"][region]: + self.config["VOIs"][region][key] = value print( f"For {region}, the parameter '{key}' was not defined by the user, set to {value}." ) @@ -165,7 +165,7 @@ def extract_masks_and_correct_overlaps(self) -> None: """Extract masks and correct overlaps between regions.""" # Inform the user if some masks are unused and therefore excluded. for roi_name in self.nm_data.masks[0]: - if roi_name not in self.config["rois"] and roi_name != "BoneMarrow": + if roi_name not in self.config["VOIs"] and roi_name != "BoneMarrow": print( f"Although mask for {roi_name} is present, we are ignoring it because this region was not included in the" " configuration input file.\n" @@ -176,7 +176,7 @@ def extract_masks_and_correct_overlaps(self) -> None: time_id: extract_masks( time_id=time_id, mask_dataset=self.nm_data.masks, - requested_rois=list(self.config["rois"].keys()), + requested_rois=list(self.config["VOIs"].keys()), ) for time_id in self.nm_data.masks.keys() } @@ -185,13 +185,13 @@ def extract_masks_and_correct_overlaps(self) -> None: time_id: extract_masks( time_id=time_id, mask_dataset=self.ct_data.masks, - requested_rois=list(self.config["rois"].keys()), + requested_rois=list(self.config["VOIs"].keys()), ) for time_id in self.ct_data.masks.keys() } # Check availability of requested rois in existing masks - for roi_name in self.config["rois"]: + for roi_name in self.config["VOIs"]: if roi_name not in self.nm_data.masks[0] and roi_name != "BoneMarrow": raise AssertionError(f"The following mask was NOT found: {roi_name}\n") @@ -231,10 +231,10 @@ def check_mandatory_fields(self) -> None: self.config["ReferenceTimePoint"] = 0 if "Organ" in self.config["Level"]: - if "WholeBody" not in self.config["rois"]: + if "WholeBody" not in self.config["VOIs"]: raise ValueError("Missing 'WholeBody' region parameters.") - if "RemainderOfBody" not in self.config["rois"]: + if "RemainderOfBody" not in self.config["VOIs"]: raise ValueError("Missing 'RemainderOfBody' region parameters.") return None @@ -244,7 +244,7 @@ def initialize(self) -> pandas.DataFrame: tmp_results: Dict[str, List[float]] = { roi_name: [] for roi_name in self.nm_data.masks[0].keys() - if roi_name in self.config["rois"] + if roi_name in self.config["VOIs"] } cols: List[str] = ["Time_hr", "Volume_CT_mL", "Activity_MBq", "Density_HU"] @@ -294,7 +294,7 @@ def initialize_bone_marrow( ) -> Dict[str, List[float]]: """Initialize activity and times for Bone-Marrow blood-based measurements.""" if ( - "BoneMarrow" in self.config["rois"] + "BoneMarrow" in self.config["VOIs"] and self.clinical_data is not None and "BoneMarrow" not in self.nm_data.masks[0] ): @@ -387,8 +387,6 @@ def sanity_checks(self, metric: str) -> None: def compute_tia(self) -> None: """Compute Time-Integrated Activity over each source-organ.""" - # decay_constant = math.log(2) / (self.radionuclide["half_life"]) # 1/h # TODO: Check how to incorporate into bounds? (flake8) - if self.radionuclide["half_life_units"] != "hours": raise AssertionError( "Radionuclide Half-Life in Database should be in hours." @@ -403,7 +401,6 @@ def compute_tia(self) -> None: } for region, region_data in self.results.iterrows(): - fit_results = self.smart_fit_selection( region_data=region_data, region=region ) @@ -443,7 +440,7 @@ def compute_tia(self) -> None: tmp_tia_data["Lambda_eff"].append( [ fit_params[exp_params[i]] - for i in range(self.config["rois"][region]["fit_order"]) + for i in range(self.config["VOIs"][region]["fit_order"]) ] ) @@ -462,15 +459,16 @@ def smart_fit_selection( ) -> lmfit.model.ModelResult: """Select the best fit based on Akaike Information Criterion.""" # If fit_order is defined by user: - if self.config["rois"][region]["fit_order"] is not None: + if self.config["VOIs"][region]["fit_order"] is not None: + print(region) fit_results, _ = exponential_fit_lmfit( x_data=numpy.array(region_data["Time_hr"]), y_data=numpy.array(region_data["Activity_MBq"]), - fixed_params=self.config["rois"][region]["fixed_parameters"], - num_exponentials=self.config["rois"][region]["fit_order"], - bounds=self.config["rois"][region]["bounds"], - params_init=self.config["rois"][region]["param_init"], - with_uptake=self.config["rois"][region]["with_uptake"], + fixed_params=self.config["VOIs"][region]["fixed_parameters"], + num_exponentials=self.config["VOIs"][region]["fit_order"], + bounds=self.config["VOIs"][region]["bounds"], + params_init=self.config["VOIs"][region]["param_init"], + with_uptake=self.config["VOIs"][region]["with_uptake"], ) return fit_results @@ -499,7 +497,7 @@ def smart_fit_selection( y_data=numpy.array(region_data["Activity_MBq"]), fixed_params=None, num_exponentials=order, - bounds=self.config["rois"][region]["bounds"], + bounds=self.config["VOIs"][region]["bounds"], params_init={"A1": activity_init}, with_uptake=with_uptake, ) @@ -513,8 +511,8 @@ def smart_fit_selection( # If only one model fit, that is the winner. if len(aic_results) == 1: - self.config["rois"][region]["with_uptake"] = fit_config[0][0] - self.config["rois"][region]["fit_order"] = fit_config[0][1] + self.config["VOIs"][region]["with_uptake"] = fit_config[0][0] + self.config["VOIs"][region]["fit_order"] = fit_config[0][1] return all_fits[0] # If there are two more models, we check the top two models and compare their AIC. If the difference @@ -528,8 +526,8 @@ def smart_fit_selection( ): best_model_idx = aic_results[1][0] - self.config["rois"][region]["with_uptake"] = fit_config[best_model_idx][0] - self.config["rois"][region]["fit_order"] = fit_config[best_model_idx][1] + self.config["VOIs"][region]["with_uptake"] = fit_config[best_model_idx][0] + self.config["VOIs"][region]["fit_order"] = fit_config[best_model_idx][1] return all_fits[best_model_idx] @@ -704,7 +702,7 @@ def save_images_and_masks_at(self, time_id: int) -> None: time_id=time_id, out_path=self.db_dir, name="SPECT" ) self.nm_data.save_masks_to_nii_at( - time_id=time_id, out_path=self.db_dir, regions=self.config["rois"] + time_id=time_id, out_path=self.db_dir, regions=self.config["VOIs"] ) return None @@ -760,6 +758,7 @@ def write_json_data( cycle["InjectionTime"] = self.config["InjectionTime"] cycle["InjectedActivity"] = self.config["InjectedActivity"] cycle["Weight_g"] = self.config["PatientWeight_g"] + cycle["Height_cm"] = self.config["PatientHeight_cm"] cycle["Level"] = self.config["Level"] if cycle["Level"] == "Organ": cycle["Method"] = self.config["OrganLevel"] @@ -771,9 +770,9 @@ def write_json_data( cycle["ReferenceTimePoint"] = self.config["ReferenceTimePoint"] cycle["TimePoints_h"] = self.results["Time_hr"][0] - for organ in self.config["rois"].keys(): - if organ not in cycle["rois"]: - cycle["rois"][organ] = { + for organ in self.config["VOIs"].keys(): + if organ not in cycle["VOIs"]: + cycle["VOIs"][organ] = { "volumes_mL": {}, "activity_MBq": {}, "timepoints_h": {}, @@ -805,32 +804,32 @@ def write_json_data( "BED_Gy_uncertainty": {}, } - cycle["rois"][organ]["volumes_mL"]["different_tps"] = self.results.loc[ + cycle["VOIs"][organ]["volumes_mL"]["different_tps"] = self.results.loc[ organ, "Volume_CT_mL" ] - cycle["rois"][organ]["volumes_mL"]["uncertainty"] = "NA" - cycle["rois"][organ]["volumes_mL"]["mean"] = numpy.mean( + cycle["VOIs"][organ]["volumes_mL"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["volumes_mL"]["mean"] = numpy.mean( self.results.loc[organ, "Volume_CT_mL"] ) - cycle["rois"][organ]["volumes_mL"]["mean_uncertainty"] = "NA" - cycle["rois"][organ]["activity_MBq"]["values"] = self.results.loc[ - organ, "Activity_MBq" + cycle["VOIs"][organ]["volumes_mL"]["mean_uncertainty"] = "NA" + cycle["VOIs"][organ]["activity_MBq"]["values"] = [ + float(x) for x in self.results.loc[organ, "Activity_MBq"] ] - cycle["rois"][organ]["activity_MBq"]["uncertainty"] = "NA" - cycle["rois"][organ]["timepoints_h"]["values"] = self.results.loc[ + cycle["VOIs"][organ]["activity_MBq"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["timepoints_h"]["values"] = self.results.loc[ organ, "Time_hr" ] - cycle["rois"][organ]["doserate_MBq_per_h"]["values"] = "NA" - cycle["rois"][organ]["doserate_MBq_per_h"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["doserate_MBq_per_h"]["values"] = "NA" + cycle["VOIs"][organ]["doserate_MBq_per_h"]["uncertainty"] = "NA" try: - cycle["rois"][organ]["density_HU"]["different_tps"] = self.results.loc[ + cycle["VOIs"][organ]["density_HU"]["different_tps"] = self.results.loc[ organ, "Density_HU" ] except (KeyError, AttributeError): # TODO: Handle errors explicitly pass - cycle["rois"][organ]["density_HU"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["density_HU"]["uncertainty"] = "NA" try: - cycle["rois"][organ]["density_HU"]["mean"] = numpy.mean( + cycle["VOIs"][organ]["density_HU"]["mean"] = numpy.mean( self.results.loc[organ, "Density_HU"] ) except ( @@ -839,92 +838,92 @@ def write_json_data( TypeError, ): # TODO: Handle errors explicitly pass - cycle["rois"][organ]["density_HU"]["mean_uncertainty"] = "NA" - cycle["rois"][organ]["density_gml"]["different_tps"] = "NA" - cycle["rois"][organ]["density_gml"]["uncertainty"] = "NA" - cycle["rois"][organ]["density_gml"]["mean"] = "NA" - cycle["rois"][organ]["density_gml"]["mean_uncertainty"] = "NA" - cycle["rois"][organ]["mass_g"]["different_tps"] = "NA" - cycle["rois"][organ]["mass_g"]["uncertainty"] = "NA" - cycle["rois"][organ]["mass_g"]["mean"] = "NA" - cycle["rois"][organ]["mass_g"]["mean_uncertainty"] = "NA" - cycle["rois"][organ]["fitting_eq"] = self.config["rois"][organ]["fit_order"] - cycle["rois"][organ]["no_of_fit_params"] = "NA" - cycle["rois"][organ]["fit_params"] = list( + cycle["VOIs"][organ]["density_HU"]["mean_uncertainty"] = "NA" + cycle["VOIs"][organ]["density_gml"]["different_tps"] = "NA" + cycle["VOIs"][organ]["density_gml"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["density_gml"]["mean"] = "NA" + cycle["VOIs"][organ]["density_gml"]["mean_uncertainty"] = "NA" + cycle["VOIs"][organ]["mass_g"]["different_tps"] = "NA" + cycle["VOIs"][organ]["mass_g"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["mass_g"]["mean"] = "NA" + cycle["VOIs"][organ]["mass_g"]["mean_uncertainty"] = "NA" + cycle["VOIs"][organ]["fitting_eq"] = self.config["VOIs"][organ]["fit_order"] + cycle["VOIs"][organ]["no_of_fit_params"] = "NA" + cycle["VOIs"][organ]["fit_params"] = list( self.results.loc[organ, "Fit_params"] ) - cycle["rois"][organ]["fit_params_uncertainty"] = "NA" - cycle["rois"][organ]["R_2"] = self.results.loc[organ, "R_squared_AIC"][0] - cycle["rois"][organ]["AIC"] = self.results.loc[organ, "R_squared_AIC"][1] - cycle["rois"][organ]["TIA_MBqh"] = self.results.loc[organ, "TIA_MBq_h"] - cycle["rois"][organ]["TIA_MBqh_uncertainty"] = "NA" - cycle["rois"][organ]["TIA_h"] = self.results.loc[organ, "TIA_h"] - cycle["rois"][organ]["TIA_h_uncertainty"] = "NA" - cycle["rois"][organ]["mean_AD_Gy"] = "NA" - cycle["rois"][organ]["mean_AD_Gy_uncertainty"] = "NA" - cycle["rois"][organ]["min_AD_Gy"] = "NA" - cycle["rois"][organ]["max_AD_Gy"] = "NA" - cycle["rois"][organ]["peak_AD_Gy"] = "NA" - cycle["rois"][organ]["repair_halflife"] = "NA" - cycle["rois"][organ]["alpha_beta"] = "NA" - cycle["rois"][organ]["composition"] = "NA" - cycle["rois"][organ]["total_s_value"] = "NA" - cycle["rois"][organ]["total_s_value_uncertainty"] = "NA" + cycle["VOIs"][organ]["fit_params_uncertainty"] = "NA" + cycle["VOIs"][organ]["R_2"] = self.results.loc[organ, "R_squared_AIC"][0] + cycle["VOIs"][organ]["AIC"] = self.results.loc[organ, "R_squared_AIC"][1] + cycle["VOIs"][organ]["TIA_MBqh"] = self.results.loc[organ, "TIA_MBq_h"] + cycle["VOIs"][organ]["TIA_MBqh_uncertainty"] = "NA" + cycle["VOIs"][organ]["TIA_h"] = self.results.loc[organ, "TIA_h"] + cycle["VOIs"][organ]["TIA_h_uncertainty"] = "NA" + cycle["VOIs"][organ]["mean_AD_Gy"] = "NA" + cycle["VOIs"][organ]["mean_AD_Gy_uncertainty"] = "NA" + cycle["VOIs"][organ]["min_AD_Gy"] = "NA" + cycle["VOIs"][organ]["max_AD_Gy"] = "NA" + cycle["VOIs"][organ]["peak_AD_Gy"] = "NA" + cycle["VOIs"][organ]["repair_halflife"] = "NA" + cycle["VOIs"][organ]["alpha_beta"] = "NA" + cycle["VOIs"][organ]["composition"] = "NA" + cycle["VOIs"][organ]["total_s_value"] = "NA" + cycle["VOIs"][organ]["total_s_value_uncertainty"] = "NA" if "Lesion" in organ or "TTB" in organ: - cycle["rois"][organ]["density_gml"]["different_tps"] = "NA" - cycle["rois"][organ]["density_gml"]["uncertainty"] = "NA" - cycle["rois"][organ]["density_gml"]["mean"] = ( + cycle["VOIs"][organ]["density_gml"]["different_tps"] = "NA" + cycle["VOIs"][organ]["density_gml"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["density_gml"]["mean"] = ( self.results_dosimetry_lesions.loc[organ, "Density_g_per_mL"] ) - cycle["rois"][organ]["density_gml"]["mean_uncertainty"] = "NA" - cycle["rois"][organ]["mass_g"]["different_tps"] = "NA" - cycle["rois"][organ]["mass_g"]["uncertainty"] = "NA" - cycle["rois"][organ]["mass_g"]["mean"] = ( + cycle["VOIs"][organ]["density_gml"]["mean_uncertainty"] = "NA" + cycle["VOIs"][organ]["mass_g"]["different_tps"] = "NA" + cycle["VOIs"][organ]["mass_g"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["mass_g"]["mean"] = ( self.results_dosimetry_lesions.loc[organ, "Mass_g"] ) - cycle["rois"][organ]["mass_g"]["mean_uncertainty"] = "NA" - cycle["rois"][organ]["composition"] = ( + cycle["VOIs"][organ]["mass_g"]["mean_uncertainty"] = "NA" + cycle["VOIs"][organ]["composition"] = ( self.results_dosimetry_lesions.loc[organ, "Composition"] ) - cycle["rois"][organ]["total_s_value"] = ( + cycle["VOIs"][organ]["total_s_value"] = ( self.results_dosimetry_lesions.loc[organ, "Total_S_Value"] ) - cycle["rois"][organ]["total_s_value_uncertainty"] = "NA" - cycle["rois"][organ]["mean_AD_Gy"] = self.results_dosimetry_lesions.loc[ + cycle["VOIs"][organ]["total_s_value_uncertainty"] = "NA" + cycle["VOIs"][organ]["mean_AD_Gy"] = self.results_dosimetry_lesions.loc[ organ, "AD_Gy" ] - cycle["rois"][organ]["mean_AD_Gy_uncertainty"] = "NA" + cycle["VOIs"][organ]["mean_AD_Gy_uncertainty"] = "NA" if "BoneMarrow" in organ: - cycle["rois"][organ]["volumes_mL"]["different_tps"] = 1170 - cycle["rois"][organ]["volumes_mL"]["uncertainty"] = "NA" - cycle["rois"][organ]["volumes_mL"]["mean"] = 1170 + cycle["VOIs"][organ]["volumes_mL"]["different_tps"] = 1170 + cycle["VOIs"][organ]["volumes_mL"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["volumes_mL"]["mean"] = 1170 if "Gland" in organ: - cycle["rois"][organ]["density_gml"]["different_tps"] = "NA" - cycle["rois"][organ]["density_gml"]["uncertainty"] = "NA" - cycle["rois"][organ]["density_gml"]["mean"] = ( + cycle["VOIs"][organ]["density_gml"]["different_tps"] = "NA" + cycle["VOIs"][organ]["density_gml"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["density_gml"]["mean"] = ( self.results_dosimetry_salivaryglands.loc[organ, "Density_g_per_mL"] ) - cycle["rois"][organ]["density_gml"]["mean_uncertainty"] = "NA" - cycle["rois"][organ]["mass_g"]["different_tps"] = "NA" - cycle["rois"][organ]["mass_g"]["uncertainty"] = "NA" - cycle["rois"][organ]["mass_g"]["mean"] = ( + cycle["VOIs"][organ]["density_gml"]["mean_uncertainty"] = "NA" + cycle["VOIs"][organ]["mass_g"]["different_tps"] = "NA" + cycle["VOIs"][organ]["mass_g"]["uncertainty"] = "NA" + cycle["VOIs"][organ]["mass_g"]["mean"] = ( self.results_dosimetry_salivaryglands.loc[organ, "Mass_g"] ) - cycle["rois"][organ]["mass_g"]["mean_uncertainty"] = "NA" - cycle["rois"][organ]["composition"] = ( + cycle["VOIs"][organ]["mass_g"]["mean_uncertainty"] = "NA" + cycle["VOIs"][organ]["composition"] = ( self.results_dosimetry_salivaryglands.loc[organ, "Composition"] ) - cycle["rois"][organ]["total_s_value"] = ( + cycle["VOIs"][organ]["total_s_value"] = ( self.results_dosimetry_salivaryglands.loc[organ, "Total_S_Value"] ) - cycle["rois"][organ]["total_s_value_uncertainty"] = "NA" - cycle["rois"][organ]["mean_AD_Gy"] = ( + cycle["VOIs"][organ]["total_s_value_uncertainty"] = "NA" + cycle["VOIs"][organ]["mean_AD_Gy"] = ( self.results_dosimetry_salivaryglands.loc[organ, "AD_Gy"] ) - cycle["rois"][organ]["mean_AD_Gy_uncertainty"] = "NA" + cycle["VOIs"][organ]["mean_AD_Gy_uncertainty"] = "NA" if self.config["Level"] == "Organ": for organ in self.results_dosimetry_organs.index: From 890e55bec384d9093f2d7f16245fd3f1e203d1dd Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 26 Nov 2025 12:16:48 -0800 Subject: [PATCH 08/15] adjust x axis based on the highest value --- pytheranostics/plots/plots.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/pytheranostics/plots/plots.py b/pytheranostics/plots/plots.py index 1d9f1ec..eaa6a24 100644 --- a/pytheranostics/plots/plots.py +++ b/pytheranostics/plots/plots.py @@ -86,7 +86,7 @@ def plot_tac_residuals( else: weights = None # Generate x-values for plotting the fitted model starting from x=0 - x_fit = numpy.linspace(0, x_data[-1] * 3, 500) + x_fit = numpy.linspace(0, x_data.max() * 3, 500) y_fit = result.eval(x=x_fit) # First subplot: Linear scale plot @@ -96,7 +96,7 @@ def plot_tac_residuals( # Plot fitted model ax1.plot(x_fit, y_fit, color="red") ax1.set_xlim(left=0) # Start x-axis from zero - ax1.set_xlim(right=x_data[-1] * 2) # Start y-axis from zero + ax1.set_xlim(right=x_data.max() * 2) # Start y-axis from zero ax1.set_ylim(bottom=0) # Start y-axis from zero ax1.set_title(region) ax1.set_xlabel(x_label) @@ -119,7 +119,7 @@ def plot_tac_residuals( # Plot fitted model ax2.plot(x_fit, y_fit, color="red") ax2.set_xlim(left=0) # Start x-axis from zero - ax2.set_xlim(right=x_data[-1] * 2) # Start y-axis from zero + ax2.set_xlim(right=x_data.max() * 2) # Start y-axis from zero ax2.set_yscale("log") ax2.set_title(title_text) ax2.set_xlabel(x_label) From 2efa4cdcee43772d63f3c7bfe23848c103877e1d Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 26 Nov 2025 12:17:13 -0800 Subject: [PATCH 09/15] add height and change from rois to VOIs --- pytheranostics/data/output.json | 22 ++++++++++++++-------- 1 file changed, 14 insertions(+), 8 deletions(-) diff --git a/pytheranostics/data/output.json b/pytheranostics/data/output.json index db61af9..ee54f14 100644 --- a/pytheranostics/data/output.json +++ b/pytheranostics/data/output.json @@ -4,8 +4,8 @@ "ClinicalTrial": "MyClinicalTrial", "Radionuclide": "NA", "PatientID": "NA", - "Gender": "NA", - + "Gender": "NA", + "No_of_completed_cycles": "NA", "Cycle_01": [ { @@ -16,13 +16,14 @@ "InjectionTime": "NA", "InjectedActivity": "NA", "Weight_g": "NA", + "Height_cm": "NA", "Level": "NA", "Method": "NA", "OutputFormat": "NA", "ScaleDoseByDensity": "NA", "ReferenceTimePoint": "NA", "TimePoints_h": "NA", - "rois": {}, + "VOIs": {}, "Organ-level_AD": {} } ], @@ -35,13 +36,14 @@ "InjectionTime": "NA", "InjectedActivity": "NA", "Weight_g": "NA", + "Height_cm": "NA", "Level": "NA", "Method": "NA", "OutputFormat": "NA", "ScaleDoseByDensity": "NA", "ReferenceTimePoint": "NA", "TimePoints_h": "NA", - "rois": {}, + "VOIs": {}, "Organ-level_AD": {} } ], @@ -54,13 +56,14 @@ "InjectionTime": "NA", "InjectedActivity": "NA", "Weight_g": "NA", + "Height_cm": "NA", "Level": "NA", "Method": "NA", "OutputFormat": "NA", "ScaleDoseByDensity": "NA", "ReferenceTimePoint": "NA", "TimePoints_h": "NA", - "rois": {}, + "VOIs": {}, "Organ-level_AD": {} } ], @@ -73,13 +76,14 @@ "InjectionTime": "NA", "InjectedActivity": "NA", "Weight_g": "NA", + "Height_cm": "NA", "Level": "NA", "Method": "NA", "OutputFormat": "NA", "ScaleDoseByDensity": "NA", "ReferenceTimePoint": "NA", "TimePoints_h": "NA", - "rois": {}, + "VOIs": {}, "Organ-level_AD": {} } ], @@ -92,6 +96,7 @@ "InjectionTime": "NA", "InjectedActivity": "NA", "Weight_g": "NA", + "Height_cm": "NA", "ApplyBiokineticsFromPreviousCycle": "NA", "Level": "NA", "Method": "NA", @@ -99,7 +104,7 @@ "ScaleDoseByDensity": "NA", "ReferenceTimePoint": "NA", "TimePoints_h": "NA", - "rois": {}, + "VOIs": {}, "Organ-level_AD": {} } ], @@ -112,6 +117,7 @@ "InjectionTime": "NA", "InjectedActivity": "NA", "Weight_g": "NA", + "Height_cm": "NA", "ApplyBiokineticsFromPreviousCycle": "NA", "Level": "NA", "Method": "NA", @@ -119,7 +125,7 @@ "ScaleDoseByDensity": "NA", "ReferenceTimePoint": "NA", "TimePoints_h": "NA", - "rois": {}, + "VOIs": {}, "Organ-level_AD": {} } ] From 9f12424afdfb82c77533b3c9d5923bc63dfbbd47 Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 26 Nov 2025 12:18:17 -0800 Subject: [PATCH 10/15] add ICRP male phantom masses --- .../phantom/human/ICRP_mass_male_source.csv | 28 +++++++++++++++++++ .../phantom/human/ICRP_mass_male_target.csv | 27 ++++++++++++++++++ 2 files changed, 55 insertions(+) create mode 100644 pytheranostics/data/phantom/human/ICRP_mass_male_source.csv create mode 100644 pytheranostics/data/phantom/human/ICRP_mass_male_target.csv diff --git a/pytheranostics/data/phantom/human/ICRP_mass_male_source.csv b/pytheranostics/data/phantom/human/ICRP_mass_male_source.csv new file mode 100644 index 0000000..3c71886 --- /dev/null +++ b/pytheranostics/data/phantom/human/ICRP_mass_male_source.csv @@ -0,0 +1,28 @@ +Organ,Mass_g +Adrenals,14 +Brain,1450 +Esophagus,40 +Eyes,15 +Gallbladder Contents,58 +LLI Cont,75 +Small Intestine,350 +Stomach Contents,250 +ULI Cont,150 +Rectum,75 +Heart Contents,510 +Heart Wall,330 +Kidneys,310 +Liver,1800 +Lungs,1200 +Pancreas,140 +Prostate,17 +Salivary Glands,85 +Red Marrow,1170 +Cortical Bone,4400 +Trabecular Bone,1100 +Spleen,150 +Testes,35 +Thymus,25 +Thyroid,20 +Urinary Bladder Contents,211 +Total Body,73000 diff --git a/pytheranostics/data/phantom/human/ICRP_mass_male_target.csv b/pytheranostics/data/phantom/human/ICRP_mass_male_target.csv new file mode 100644 index 0000000..498bc57 --- /dev/null +++ b/pytheranostics/data/phantom/human/ICRP_mass_male_target.csv @@ -0,0 +1,27 @@ +Organ,Mass_g +Adrenals,14 +Brain,1450 +Esophagus,40 +Eyes,15 +Gallbladder Wall,10 +LLI Wall,150 +Small Intestine,650 +Stomach Wall,150 +ULI Wall,150 +Rectum,70 +Heart Wall,330 +Kidneys,310 +Liver,1800 +Lungs,1200 +Pancreas,140 +Prostate,17 +Salivary Glands,85 +Red Marrow,1170 +Cortical Bone,60 +Trabecular Bone,60 +Spleen,150 +Testes,35 +Thymus,25 +Thyroid,20 +Urinary Bladder Wall,50 +Total Body,73000 From 2a159e5b716a7ade551fba18effc5f775cee74f5 Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 26 Nov 2025 12:19:38 -0800 Subject: [PATCH 11/15] add washout ratio for biexponential washout --- pytheranostics/dosimetry/base_dosimetry.py | 30 +++++++++++++++++++--- 1 file changed, 26 insertions(+), 4 deletions(-) diff --git a/pytheranostics/dosimetry/base_dosimetry.py b/pytheranostics/dosimetry/base_dosimetry.py index c7e19d5..8af3cd5 100644 --- a/pytheranostics/dosimetry/base_dosimetry.py +++ b/pytheranostics/dosimetry/base_dosimetry.py @@ -151,6 +151,7 @@ def default_config(self) -> None: "with_uptake": False, "fit_order": 1, "bounds": None, + "washout_ratio": None, } for key, value in defaults.items(): @@ -232,10 +233,19 @@ def check_mandatory_fields(self) -> None: if "Organ" in self.config["Level"]: if "WholeBody" not in self.config["VOIs"]: - raise ValueError("Missing 'WholeBody' region parameters.") + if "No" in self.config["OrganLevel"]["AdditionalOptions"]["WholeBody"]: + pass + else: + raise ValueError("Missing 'WholeBody' region parameters.") if "RemainderOfBody" not in self.config["VOIs"]: - raise ValueError("Missing 'RemainderOfBody' region parameters.") + if ( + "No" + in self.config["OrganLevel"]["AdditionalOptions"]["RemainderOfBody"] + ): + pass + else: + raise ValueError("Missing 'RemainderOfBody' region parameters.") return None @@ -469,6 +479,7 @@ def smart_fit_selection( bounds=self.config["VOIs"][region]["bounds"], params_init=self.config["VOIs"][region]["param_init"], with_uptake=self.config["VOIs"][region]["with_uptake"], + washout_ratio=self.config["VOIs"][region]["washout_ratio"], ) return fit_results @@ -852,9 +863,20 @@ def write_json_data( cycle["VOIs"][organ]["fit_params"] = list( self.results.loc[organ, "Fit_params"] ) + cycle["VOIs"][organ]["washout_ratio"] = self.config["VOIs"][organ][ + "washout_ratio" + ] cycle["VOIs"][organ]["fit_params_uncertainty"] = "NA" - cycle["VOIs"][organ]["R_2"] = self.results.loc[organ, "R_squared_AIC"][0] - cycle["VOIs"][organ]["AIC"] = self.results.loc[organ, "R_squared_AIC"][1] + cycle["VOIs"][organ]["R_2"] = ( + "NA" + if pandas.isna(self.results.loc[organ, "R_squared_AIC"][0]) + else self.results.loc[organ, "R_squared_AIC"][0] + ) + cycle["VOIs"][organ]["AIC"] = ( + "NA" + if pandas.isna(self.results.loc[organ, "R_squared_AIC"][1]) + else self.results.loc[organ, "R_squared_AIC"][1] + ) cycle["VOIs"][organ]["TIA_MBqh"] = self.results.loc[organ, "TIA_MBq_h"] cycle["VOIs"][organ]["TIA_MBqh_uncertainty"] = "NA" cycle["VOIs"][organ]["TIA_h"] = self.results.loc[organ, "TIA_h"] From cc77931960d77c169ba71a717c38a4d28e177e15 Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 26 Nov 2025 12:20:38 -0800 Subject: [PATCH 12/15] add washout_ratio --- pytheranostics/misc_tools/tools.py | 68 +++++++++++++++++++++--------- 1 file changed, 49 insertions(+), 19 deletions(-) diff --git a/pytheranostics/misc_tools/tools.py b/pytheranostics/misc_tools/tools.py index c355428..00ad17e 100644 --- a/pytheranostics/misc_tools/tools.py +++ b/pytheranostics/misc_tools/tools.py @@ -199,9 +199,9 @@ def extract_exponential_params_from_json( and all the parameters of the fit. """ with_uptake = False - # Determine order: - parameters = json_data[cycle][0]["rois"][region]["fit_params"] + parameters = json_data[cycle][0]["VOIs"][region]["fit_params"] + washout_ratio = json_data[cycle][0]["VOIs"][region]["washout_ratio"] # Handle Legacy: if len(parameters) in [3, 5]: @@ -214,7 +214,6 @@ def extract_exponential_params_from_json( fixed_parameters: Dict[str, float] = {} all_parameters: Dict[str, float] = {} - for order in range(fit_order): fixed_parameters[f"{param_name_base[order]}2"] = parameters[ exponential_idxs[order] @@ -226,14 +225,16 @@ def extract_exponential_params_from_json( all_parameters[f"{param_name_base[order]}2"] = parameters[ exponential_idxs[order] ] - + if washout_ratio is not None: + # Fix A1 ONLY + if "B1" in all_parameters: + fixed_parameters["B1"] = all_parameters["B1"] if fit_order == 2 and parameters[0] == -parameters[2]: with_uptake = True if fit_order == 3 and parameters[4] == -(parameters[0] + parameters[2]): with_uptake = True - - return fixed_parameters, with_uptake, all_parameters + return fixed_parameters, with_uptake, all_parameters, washout_ratio def extract_exponential_params_from_json_legacy( @@ -266,7 +267,7 @@ def extract_exponential_params_from_json_legacy( When the parameter configuration is incompatible with new format. """ # Read Parameters: - parameters = json_data[cycle][0]["rois"][region]["fit_params"] + parameters = json_data[cycle][0]["VOIs"][region]["fit_params"] if len(parameters) not in [3, 5]: raise AssertionError( @@ -321,7 +322,7 @@ def initialize_biokinetics_from_prior_cycle( dict Updated configuration dictionary. """ - for roi, roi_info in config["rois"].items(): + for roi, roi_info in config["VOIs"].items(): if ( "biokinectics_from_previous_cycle" in roi_info @@ -329,15 +330,18 @@ def initialize_biokinetics_from_prior_cycle( ): # Get previous cycle parameters: - fixed_param, with_uptake, all_params = extract_exponential_params_from_json( - json_data=prior_treatment_data, cycle=cycle, region=roi + fixed_param, with_uptake, all_params, washout_ratio = ( + extract_exponential_params_from_json( + json_data=prior_treatment_data, cycle=cycle, region=roi + ) ) - config["rois"][roi] = { + config["VOIs"][roi] = { "fixed_parameters": fixed_param, "fit_order": len(all_params) // 2, "param_init": all_params, "with_uptake": with_uptake, + "washout_ratio": washout_ratio, } print( @@ -391,7 +395,7 @@ def plot_MIP(SPECT=None, vmax=300000, figsize=(10, 5), ax=None): return fig, ax -def plot_MIP_with_mask_outlines(ax, SPECT, masks=None, vmax=300000): +def plot_MIP_with_mask_outlines(ax, SPECT, masks=None, vmax=300000, label=None): """Plot Maximum Intensity Projection (MIP) of SPECT data with masks outlines. Parameters @@ -412,14 +416,20 @@ def plot_MIP_with_mask_outlines(ax, SPECT, masks=None, vmax=300000): if masks is not None: for organ, mask in masks.items(): organ_lower = organ.lower() - if "kidney" in organ_lower: - color = "lime" - elif "gland" in organ_lower: - color = "red" - elif "lesion" in organ_lower: - color = "m" - else: + print(organ_lower) + if "peak" in organ_lower: continue + else: + if "kidney" in organ_lower: + color = "lime" + elif "parotid" in organ_lower: + color = "red" + elif "submandibular" in organ_lower: + color = "red" + elif "lesion" in organ_lower: + color = "m" + else: + continue mip_mask = mask.max(axis=0) if mip_mask.shape != spect_mip.shape: @@ -433,6 +443,26 @@ def plot_MIP_with_mask_outlines(ax, SPECT, masks=None, vmax=300000): alpha=0.5, ) + # --- Add label at mask center --- + if label is True: + ys, xs = numpy.where(mip_mask > 0) + + if len(xs) > 0: + # Corrected for transpose in contour + x_center = ys.mean() - 0.2 * ys.mean() # was xs + y_center = xs.mean() # was ys + + plt.text( + x_center, + y_center, + organ, + color=color, + fontsize=8, + ha="center", + va="center", + alpha=0.7, + ) + plt.xlim(30, 100) plt.ylim(0, 234) plt.axis("off") From ce5b42c9837fde624bed26dbe0332fe572ea4aca Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 26 Nov 2025 14:05:20 -0800 Subject: [PATCH 13/15] improve the report generator --- pytheranostics/misc_tools/report_generator.py | 245 +++++++++++++----- 1 file changed, 179 insertions(+), 66 deletions(-) diff --git a/pytheranostics/misc_tools/report_generator.py b/pytheranostics/misc_tools/report_generator.py index 759540f..4c8bcf6 100644 --- a/pytheranostics/misc_tools/report_generator.py +++ b/pytheranostics/misc_tools/report_generator.py @@ -11,6 +11,7 @@ from reportlab.lib.units import inch from reportlab.platypus import ( Image, + PageBreak, Paragraph, SimpleDocTemplate, Spacer, @@ -25,21 +26,33 @@ def signature_block(person, styles, width=2.5 * inch, height=0.6 * inch): Placeholder for signature, line, and text. Returns a small table that can be inserted side-by-side with others. """ - # Empty row for signature space - sig_space = Spacer(1, height) + elements = [] + + # If a signature path is provided, insert the image + if person.get("signature"): + sig_img = Image(person["signature"]) + # Scale the image to fit width and maintain aspect ratio + sig_img.drawHeight = height + sig_img.drawWidth = width + elements.append(sig_img) + else: + # Empty space if no signature image + elements.append(Spacer(1, height)) # Line row line = Table([[""]], colWidths=[width]) line.setStyle(TableStyle([("LINEABOVE", (0, 0), (-1, -1), 1, colors.black)])) + elements.append(line) # Text row text = Paragraph( f"{person['name']}
{person['title']}
{person['affiliation']}", styles["Normal"], ) + elements.append(text) - # Wrap everything in a column table (1 col, 3 rows) - block = Table([[sig_space], [line], [text]], colWidths=[width]) + # Wrap everything in a column table (1 col, stacked) + block = Table([[e] for e in elements], colWidths=[width]) block.setStyle(TableStyle([("ALIGN", (0, 0), (-1, -1), "CENTER")])) return block @@ -69,23 +82,26 @@ def create_dosimetry_pdf(json_file, output_file, calculated_by=None, approved_by # Title title = Paragraph( - "DOSIMETRY REPORT", styles["Title"] + "DOSIMETRY ASSESSMENT", styles["Title"] ) elements.append(title) elements.append(Spacer(1, 0.5 * inch)) # Subject Information Section elements.append(Paragraph("Subject Information", styles["Heading2"])) - # Subject Information Table subject_data = [ - ["Clinical Trial", "PR.21"], + ["Clinical Trial", data.get("ClinicalTrial")], + ["Radiopharmaceutical", "177Lu-PSMA-617"], + ["Mode of administration", "I.V."], ["ID", data.get("PatientID")], ["Sex", data.get("Gender")], + ["Weight kg", data.get("Cycle_01", {})[0].get("Weight_g", "") / 1000], + ["Height cm", data.get("Cycle_01", {})[0].get("Height_cm", "")], ["Number of cycles ", data.get("No_of_completed_cycles")], ] - subject_table = Table(subject_data, colWidths=[1.5 * inch, 4 * inch]) + subject_table = Table(subject_data, colWidths=[2 * inch, 3.5 * inch]) subject_table.setStyle( TableStyle( [ @@ -146,11 +162,89 @@ def create_dosimetry_pdf(json_file, output_file, calculated_by=None, approved_by elements.append(caption) elements.append(Spacer(1, 0.2 * inch)) + elements.append(PageBreak()) + for i in range(1, data.get("No_of_completed_cycles") + 1): cycle_info(i, elements, styles, data) + elements.append(PageBreak()) + elements.append(Paragraph("Patient Summary", styles["Heading2"])) + # =============================== + # CUMULATIVE ORGAN TABLE + # =============================== + elements.append( + Paragraph("Cumulative Absorbed Dose Summary", styles["Heading3"]) + ) + + total_tia_kidneys = 0 + total_tia_salivary = 0 + total_tia_marrow = 0 + + total_ad_kidneys = 0 + total_ad_salivary = 0 + total_ad_marrow = 0 + + total_bed_kidneys = 0 + + for i in range(1, data.get("No_of_completed_cycles") + 1): + therapy_info = data.get(f"Cycle_0{i}", {})[0] + + # Kidneys + total_tia_kidneys += ( + therapy_info["VOIs"]["Kidney_Left"]["TIA_h"] + + therapy_info["VOIs"]["Kidney_Right"]["TIA_h"] + ) + total_ad_kidneys += therapy_info["Organ-level_AD"]["Kidneys"]["AD[Gy]"] + total_bed_kidneys += therapy_info["Organ-level_AD"]["Kidneys"]["BED[Gy]"] + + # Red marrow + total_tia_marrow += therapy_info["VOIs"]["BoneMarrow"]["TIA_h"] + total_ad_marrow += therapy_info["Organ-level_AD"]["Red Marrow"]["AD[Gy]"] + + # Salivary glands + total_tia_salivary += ( + therapy_info["VOIs"]["ParotidGland_Left"]["TIA_h"] + + therapy_info["VOIs"]["ParotidGland_Right"]["TIA_h"] + + therapy_info["VOIs"]["SubmandibularGland_Left"]["TIA_h"] + + therapy_info["VOIs"]["SubmandibularGland_Right"]["TIA_h"] + ) + total_ad_salivary += therapy_info["Organ-level_AD"]["Salivary Glands"]["AD[Gy]"] + + # Build the cumulative table + cumulative_data = [ + ["Organ", "Cumulative TIA (h)", "Cumulative AD (Gy)", "Cumulative BED (Gy)"], + [ + "Kidneys", + round(total_tia_kidneys, 2), + round(total_ad_kidneys, 2), + round(total_bed_kidneys, 2), + ], + ["Red Marrow", round(total_tia_marrow, 2), round(total_ad_marrow, 2), "-"], + [ + "Salivary glands", + round(total_tia_salivary, 2), + round(total_ad_salivary, 2), + "-", + ], + ] + + cumulative_table = Table(cumulative_data, colWidths=[1.5 * inch, 1.7 * inch]) + cumulative_table.setStyle( + TableStyle( + [ + ("ALIGN", (0, 0), (-1, -1), "LEFT"), + ("VALIGN", (0, 0), (-1, -1), "MIDDLE"), + ("FONTNAME", (0, 0), (-1, -1), "Helvetica"), + ("FONTSIZE", (0, 0), (-1, -1), 12), + ("GRID", (0, 0), (-1, -1), 1, colors.black), + ("BACKGROUND", (0, 0), (0, -1), colors.lightblue), + ] + ) + ) + + elements.append(cumulative_table) # Paths to your three images image_paths = [ calling_folder / "TestDoseDB/Gy_cummulative.png", @@ -162,7 +256,7 @@ def create_dosimetry_pdf(json_file, output_file, calculated_by=None, approved_by imgs = [] for path in image_paths: img = Image(str(path)) - scale = min(max_width / img.imageWidth, max_height / img.imageHeight) / 2.8 + scale = min(max_width / img.imageWidth, max_height / img.imageHeight) / 3 img.drawWidth = img.imageWidth * scale img.drawHeight = img.imageHeight * scale imgs.append(img) @@ -186,7 +280,7 @@ def create_dosimetry_pdf(json_file, output_file, calculated_by=None, approved_by # Add caption elements.append(Spacer(1, 0.2 * inch)) caption = Paragraph( - "Figure 2: Cumulative absorbed dose, absorbed dose per cycle, and absorbed dose per GBq per cycle for target organs.", + "Figure 2: Cumulative AD, AD per cycle, and AD per GBq per cycle for target organs.", styles["Normal"], ) elements.append(caption) @@ -205,7 +299,7 @@ def create_dosimetry_pdf(json_file, output_file, calculated_by=None, approved_by img = Image(str(path)) # Scale to fit 2×2 layout scale = min( - (max_width / 2.2) / img.imageWidth, (max_height / 2.2) / img.imageHeight + (max_width / 2.5) / img.imageWidth, (max_height / 2.5) / img.imageHeight ) img.drawWidth = img.imageWidth * scale img.drawHeight = img.imageHeight * scale @@ -214,7 +308,7 @@ def create_dosimetry_pdf(json_file, output_file, calculated_by=None, approved_by # Arrange in 2×2 structure trend_table_data = [[trend_imgs[0], trend_imgs[1]], [trend_imgs[2], trend_imgs[3]]] - trend_table = Table(trend_table_data, colWidths=[max_width / 2.2, max_width / 2.2]) + trend_table = Table(trend_table_data, colWidths=[max_width / 2.5, max_width / 2.5]) trend_table.setStyle( TableStyle( @@ -226,18 +320,14 @@ def create_dosimetry_pdf(json_file, output_file, calculated_by=None, approved_by ] ) ) - - elements.append(Spacer(1, 0.3 * inch)) elements.append(trend_table) # Add caption - elements.append(Spacer(1, 0.2 * inch)) caption = Paragraph( "Figure 3: Trends of hematological and renal function, and PSA.", styles["Normal"], ) elements.append(caption) - elements.append(Spacer(1, 0.2 * inch)) # =============================== # Signatures Section @@ -272,11 +362,27 @@ def create_dosimetry_pdf(json_file, output_file, calculated_by=None, approved_by ) elements.append(app_table) + # =============================== + # Appendix + # =============================== + + elements.append(PageBreak()) + title = Paragraph("Appendix", styles["Heading2"]) + elements.append(title) + fig_title = Paragraph( + "Biodistribution and kinetic analysis", + styles["Heading3"], + ) + elements.append(fig_title) + + for i in range(1, data.get("No_of_completed_cycles") + 1): + biodistribution_per_cycle(i, elements, styles, data) + # Build PDF doc.build(elements) -def cycle_info(cycle_n, elements, styles, data): +def biodistribution_per_cycle(cycle_n, elements, styles, data): """Add cycle information to the PDF report elements. Parameters @@ -290,54 +396,18 @@ def cycle_info(cycle_n, elements, styles, data): data : dict Patient data dictionary. """ - # Therapy Information Section + if cycle_n == 1: + pass + else: + elements.append(PageBreak()) + therapy_title = Paragraph(f"Cycle {cycle_n}", styles["Heading2"]) elements.append(therapy_title) - # Therapy Information Table - therapy_info = data.get(f"Cycle_0{cycle_n}", {}) - therapy_data = [ - ["Radiopharmaceutical", "177Lu-PSMA-617"], - ["Mode of administration", "I.V."], - ["Administered Activity (MBq)", therapy_info[0].get("InjectedActivity", "")], - [ - "Date of injection", - ( - datetime.strptime( - therapy_info[0].get("InjectionDate", ""), "%Y%m%d" - ).strftime("%Y-%m-%d") - if therapy_info[0].get("InjectionDate", "") - else "" - ), - ], - ] - - therapy_table = Table(therapy_data, colWidths=[2.5 * inch, 3 * inch]) - therapy_table.setStyle( - TableStyle( - [ - ("ALIGN", (0, 0), (-1, -1), "LEFT"), - ("VALIGN", (0, 0), (-1, -1), "MIDDLE"), - ("FONTNAME", (0, 0), (-1, -1), "Helvetica"), - ("FONTSIZE", (0, 0), (-1, -1), 12), - ("GRID", (0, 0), (-1, -1), 1, colors.black), - ("BACKGROUND", (0, 0), (0, -1), colors.lightgrey), - ] - ) - ) - - elements.append(therapy_table) - page_width, page_height = letter max_width = page_width - 2 * 72 # 1-inch margins max_height = page_height - 2 * 72 - fig_title = Paragraph( - "Time-activity curves, fit functions, and fit parameters", - styles["Heading3"], - ) - elements.append(fig_title) - calling_folder = Path().absolute() # notebook folder pattern = calling_folder / f"TestDoseDB/*_fit_Cycle_0{cycle_n}.png" @@ -357,19 +427,62 @@ def cycle_info(cycle_n, elements, styles, data): elements.append(Spacer(1, 0.3 * inch)) elements.append(img) + +def cycle_info(cycle_n, elements, styles, data): + """Add cycle information to the PDF report elements. + + Parameters + ---------- + cycle_n : int + The cycle number. + elements : list + List of reportlab elements to append to. + styles : dict + ReportLab styles dictionary. + data : dict + Patient data dictionary. + """ + # Therapy Information Section + therapy_title = Paragraph(f"Cycle {cycle_n}", styles["Heading2"]) + elements.append(therapy_title) + + # Therapy Information Table + therapy_info = data.get(f"Cycle_0{cycle_n}", {}) + + # Format injection date nicely + inj_date_raw = therapy_info[0].get("InjectionDate", "") + inj_date = ( + datetime.strptime(inj_date_raw, "%Y%m%d").strftime("%Y-%m-%d") + if inj_date_raw + else "" + ) + + # Build a clean paragraph instead of a table + therapy_info_para = Paragraph( + f"" + f"Administered activity: {therapy_info[0].get('InjectedActivity', '')} MBq
" + f"Date of injection: {inj_date}" + f"", + styles["Normal"], + ) + + # Add to document + elements.append(therapy_info_para) + elements.append(Spacer(1, 0.089 * inch)) + fig_title = Paragraph( "Absorbed dose results for the organs at risk", styles["Heading3"] ) elements.append(fig_title) organ_data_Gy_GBq = [ - ["Organ", "TIA (h)", "AD (Gy/GBq)", "AD(Gy)", "BED (Gy)"], + ["Organ", "TIA (h)", "AD (Gy/GBq)", "AD (Gy)", "BED (Gy)"], [ "Kidneys", round( ( - therapy_info[0]["rois"]["Kidney_Left"]["TIA_h"] - + therapy_info[0]["rois"]["Kidney_Right"]["TIA_h"] + therapy_info[0]["VOIs"]["Kidney_Left"]["TIA_h"] + + therapy_info[0]["VOIs"]["Kidney_Right"]["TIA_h"] ), 2, ), @@ -379,7 +492,7 @@ def cycle_info(cycle_n, elements, styles, data): ], [ "Red Marrow", - round((therapy_info[0]["rois"]["BoneMarrow"]["TIA_h"]), 2), + round((therapy_info[0]["VOIs"]["BoneMarrow"]["TIA_h"]), 2), round(therapy_info[0]["Organ-level_AD"]["Red Marrow"]["AD[Gy/GBq]"], 2), round(therapy_info[0]["Organ-level_AD"]["Red Marrow"]["AD[Gy]"], 2), "-", @@ -388,10 +501,10 @@ def cycle_info(cycle_n, elements, styles, data): "Salivary glands", round( ( - therapy_info[0]["rois"]["ParotidGland_Left"]["TIA_h"] - + therapy_info[0]["rois"]["ParotidGland_Right"]["TIA_h"] - + therapy_info[0]["rois"]["SubmandibularGland_Left"]["TIA_h"] - + therapy_info[0]["rois"]["SubmandibularGland_Right"]["TIA_h"] + therapy_info[0]["VOIs"]["ParotidGland_Left"]["TIA_h"] + + therapy_info[0]["VOIs"]["ParotidGland_Right"]["TIA_h"] + + therapy_info[0]["VOIs"]["SubmandibularGland_Left"]["TIA_h"] + + therapy_info[0]["VOIs"]["SubmandibularGland_Right"]["TIA_h"] ), 2, ), From 9df2dc08ce7ba8e9e36c02e24963c0bc50853893 Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 26 Nov 2025 14:06:14 -0800 Subject: [PATCH 14/15] add washout ratio --- pytheranostics/fits/fits.py | 15 +++++++++++++++ 1 file changed, 15 insertions(+) diff --git a/pytheranostics/fits/fits.py b/pytheranostics/fits/fits.py index a5cc94c..e8874dc 100644 --- a/pytheranostics/fits/fits.py +++ b/pytheranostics/fits/fits.py @@ -23,6 +23,7 @@ def exponential_fit_lmfit( bounds: Optional[Dict[str, Tuple[float, float]]] = None, params_init: Optional[Dict[str, float]] = None, with_uptake: bool = False, + washout_ratio: Optional[int] = None, ) -> Tuple[lmfit.model.ModelResult, Callable]: """Fit data to a sum of exponentials with flexible parameter fixing using lmfit. @@ -105,6 +106,8 @@ def model_func(x, **params): will_have_expr = True if num_exponentials == 3 and with_uptake and amp_name == "C1": will_have_expr = True + if num_exponentials == 2 and washout_ratio and amp_name == "B2": + will_have_expr = True if amp_name in fixed_params: params.add(amp_name, value=fixed_params[amp_name], vary=False) @@ -158,6 +161,18 @@ def model_func(x, **params): "Parameters 'A1', 'B1', and 'C1' must be present to apply the constraint 'C1 = -(A1 + B1)'." ) + if num_exponentials == 2 and washout_ratio is not None: + if "B2" in params and "A2" in params: + print( + f"Applying constant ratio between decay constants of washout phases: B2 = {washout_ratio} * A2" + ) + params["B2"].set(expr=f"{washout_ratio} * A2") + else: + raise ValueError( + "Parameters 'A2' and 'B2' must be present to apply the constraint 'B2 = washout_ratio * A2'." + ) + + print(y_data, params) # Perform the fit result = model.fit( y_data, params, x=x_data, weights=None if sigma is None else 1.0 / sigma From fd78b26ceb215f157cd304ac1459f5569e0e3c77 Mon Sep 17 00:00:00 2001 From: sarakurkowska Date: Wed, 26 Nov 2025 14:09:50 -0800 Subject: [PATCH 15/15] update organ dosimetry --- pytheranostics/dosimetry/organ_s_dosimetry.py | 83 ++++++++++--------- 1 file changed, 43 insertions(+), 40 deletions(-) diff --git a/pytheranostics/dosimetry/organ_s_dosimetry.py b/pytheranostics/dosimetry/organ_s_dosimetry.py index 647709c..0f460d5 100644 --- a/pytheranostics/dosimetry/organ_s_dosimetry.py +++ b/pytheranostics/dosimetry/organ_s_dosimetry.py @@ -55,13 +55,21 @@ def check_mandatory_fields_organ(self) -> None: return None @staticmethod - def _load_human_mass_table() -> pandas.DataFrame: + def _load_human_mass_target_organs_table(self) -> pandas.DataFrame: """Load the reference human phantom masses.""" with resource_path( - "pytheranostics.data", "phantom/human/human_phantom_masses.csv" + "pytheranostics.data", + f"phantom/human/ICRP_mass_{self.config['Gender'].lower()}_target.csv", + ) as masses_path: + self.target_organ_masses = pandas.read_csv(masses_path, index_col=0) + + def _load_human_mass_source_organs_table(self) -> pandas.DataFrame: + """Load the reference human phantom masses.""" + with resource_path( + "pytheranostics.data", + f"phantom/human/ICRP_mass_{self.config['Gender'].lower()}_source.csv", ) as masses_path: - masses = pandas.read_csv(masses_path, index_col=0) - return masses + self.source_organ_masses = pandas.read_csv(masses_path, index_col=0) def composition_and_density_from_HU(self, density: float) -> Tuple[str, float]: """Determine composition and density for a given CT HU value.""" @@ -250,10 +258,8 @@ def prepare_data(self) -> None: "LesionDosimetry" ): self.results_fitting_organs = self.results_fitting.copy() - if "TotalTumorBurden" in self.results_fitting.index: self.results_fitting.drop("TotalTumorBurden", axis=0, inplace=True) - if output_type == "Export": fmt = organ_conf["Output"]["ExportFormat"] if fmt.lower() == "olinda": @@ -273,6 +279,10 @@ def prepare_data(self) -> None: print("Not Implemented yet.") else: print(f"Export format {fmt} not recognized.") + elif "Lesion" in self.config["Level"]: + self.results_fitting_lesions = self.results[ + self.results.index.str.contains("Lesion") + ] return None @@ -309,7 +319,6 @@ def process_dosimetry(self) -> None: output_type = organ_conf["Output"]["Type"] # 'Export' or 'Calculate' self.prepare_data() - print("Processing dosimetry at the organ level of OARs") if output_type == "Export": fmt = organ_conf["Output"]["ExportFormat"] @@ -324,7 +333,8 @@ def process_dosimetry(self) -> None: print(f"Export format {fmt} not recognized.") elif output_type == "Calculate": - self.calculate_absorbed_dose() + if "Organ" in self.config["Level"]: + self.calculate_absorbed_dose() else: raise ValueError(f"Unknown output type: {output_type}") @@ -371,6 +381,9 @@ def calculate_absorbed_dose(self) -> pandas.DataFrame: print("Source organs available in the model:", svalues_beta.columns.tolist()) print("Source organs present :", self.results_fitting_organs.index.tolist()) + self._load_human_mass_target_organs_table() + self._load_human_mass_source_organs_table() + self.source_organs_missing = set(svalues_beta.columns) - set( self.results_fitting.index ) @@ -452,14 +465,14 @@ def redistribute_ROB_into_source_organs_missing( if "RemainderOfBody" not in tia_series.index: return tia_series - missing_organs_df = self.organ_masses.loc[ - self.organ_masses.index.isin(self.source_organs_missing) + missing_organs_df = self.source_organ_masses.loc[ + self.source_organ_masses.index.isin(self.source_organs_missing) ] print(f"Missing organs DataFrame:\n{missing_organs_df.index.tolist()}") missing_organs_df = missing_organs_df.drop(index="Total Body") - mass_source_organs = self.organ_masses.loc[ + mass_source_organs = self.source_organ_masses.loc[ [org for org in tia_series.index if org != "RemainderOfBody"], "Mass_g" ].sum() @@ -488,13 +501,6 @@ def redistribute_ROB_into_source_organs_missing( def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: """Multiply S-values by TIA to compute dose matrix for radiation type.""" - # Path to organ masses - masses = self._load_human_mass_table() - gender = self.config.get("Gender", "Male") - if gender not in masses.columns: - raise ValueError(f"Unknown gender '{gender}' for mass table.") - self.organ_masses = masses[[gender]].rename(columns={gender: "Mass_g"}) - print(f"Applying S-values for {radiation_type} radiation...") # Handle remainder of the body # Redistribute ROB TIA into missing source organs if needed - approach consistent with MIRDcalc software @@ -564,7 +570,7 @@ def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: # ROB print("Calculating Remainder of Body contributions...") dose_df["RemainderOfBody"] = 0 - total_body_mass = self.organ_masses.loc["Total Body", "Mass_g"] + total_body_mass = self.source_organ_masses.loc["Total Body", "Mass_g"] source_organs = tia_series.index if radiation_type == "beta": target_organs = s_values_subset.index.difference( @@ -573,8 +579,9 @@ def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: ROB_mass = tia_df.loc["RemainderOfBody", "Volume_CT_mL"] if radiation_type == "gamma": target_organs = s_values_subset.index - total_mass_source_organs = self.organ_masses.loc[ - self.organ_masses.index.intersection(tia_series.index), "Mass_g" + total_mass_source_organs = self.source_organ_masses.loc[ + self.source_organ_masses.index.intersection(tia_series.index), + "Mass_g", ].sum() ROB_mass = total_body_mass - total_mass_source_organs @@ -603,7 +610,7 @@ def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: if target_organ == "Total Body": continue - source_organ_mass = self.organ_masses.loc[ + source_organ_mass = self.source_organ_masses.loc[ source_organ, "Mass_g" ] @@ -625,10 +632,9 @@ def apply_s_value(self, tia_df, s_values, radiation_type) -> pandas.DataFrame: * (total_body_mass / ROB_mass) ) - contribution_from_sources - dose_df.at[target_organ, "RemainderOfBody"] = ( + dose_df.at[target_organ, "RemainderOfBody"] = float( s_value_ROB_to_target * tia_df.loc["RemainderOfBody", "TIA_s"] ) - else: # TODO print("No ROB option selected. Proceeding without ROB handling.") @@ -640,26 +646,24 @@ def perform_mass_scaling( df: pandas.DataFrame, ) -> pandas.DataFrame: """Apply mass scaling to absorbed dose calculations based on patient-specific organ masses.""" - masses = self._load_human_mass_table() - if self.config["Gender"] not in masses.columns: - raise ValueError( - f"Unknown gender '{self.config["Gender"]}' for mass table." - ) - model_masses_df = masses[[self.config["Gender"]]].rename( - columns={self.config["Gender"]: "Mass_g"} - ) - print("Performing mass scaling...") for organ in df["Target organ"]: - # Temporary mapping just for internal calculations - organ_internal = "RemainderOfBody" if organ == "Total Body" else organ + mass_key = ( + "Total Body" if organ in ["Total Body", "RemainderOfBody"] else organ + ) + fit_key = ( + "RemainderOfBody" + if organ in ["Total Body", "RemainderOfBody"] + else organ + ) + if ( - organ_internal in model_masses_df.index - and organ_internal in self.results_fitting.index + mass_key in self.target_organ_masses.index + and fit_key in self.results_fitting.index ): - model_mass = model_masses_df.loc[organ_internal, "Mass_g"] - patient_mass = self.results_fitting.loc[organ_internal, "Volume_CT_mL"] + model_mass = self.target_organ_masses.loc[mass_key, "Mass_g"] + patient_mass = self.results_fitting.loc[fit_key, "Volume_CT_mL"] if ( pandas.notna(patient_mass) @@ -669,7 +673,6 @@ def perform_mass_scaling( if "AD_beta[Gy/GBq]" in df.columns: scaling_factor_beta = model_mass / patient_mass - print(f"Scaling factor for {organ} [β]: {scaling_factor_beta}") df.loc[ df["Target organ"] == organ, "AD_beta[Gy/GBq]" ] *= scaling_factor_beta