Skip to content

Statistics

PLS-DA

hide_deconv.statistic.run_plsda(data, sample_sheet, sample_id_col, cohort_col, out_path, datasets_to_map=[], labels_data_map=[])

Runs a PLS-DA (PLS2) model and saves the corresponding plots

Parameters:

Name Type Description Default
data DataFrame

Estimated composition to be used for PLS-DA

 required
sample_sheet DataFrame

Samples sheet holding clinical metainformation

 required
sample_id_col str

Column name linking the sample sheet with the estimated compositions

 required
cohort_col str

Column name holding the different cohorts

 required
out_path Path

Path, where the created plots will be stored

 required

Returns:

Type Description
DataFrame

Estimated Scores
Source code in src/hide_deconv/statistic/plsda.py 
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
def run_plsda(
    data: pd.DataFrame,
    sample_sheet: pd.DataFrame,
    sample_id_col: str,
    cohort_col: str,
    out_path: Path,
    datasets_to_map: list[pd.DataFrame] = [],
    labels_data_map: list[str] = [],
) -> pd.DataFrame:
    """
    Runs a PLS-DA (PLS2) model and saves the corresponding plots


    Parameters
    ----------
    data : pd.DataFrame
        Estimated composition to be used for PLS-DA
    sample_sheet : pd.DataFrame
        Samples sheet holding clinical metainformation
    sample_id_col : str
        Column name linking the sample sheet with the estimated compositions
    cohort_col : str
        Column name holding the different cohorts
    out_path : Path
        Path, where the created plots will be stored

    Returns
    -------
    pd.DataFrame
        Estimated Scores
    """

    data, labels = prepare_plsda_inputs(data, sample_sheet, sample_id_col, cohort_col)

    classes = sorted(labels.unique())
    y = pd.get_dummies(labels, dtype=float).reindex(columns=classes, fill_value=0.0)

    if len(datasets_to_map) != len(labels_data_map):
        raise ValueError("Mapped datasets and labels must have the same length.")

    x = data.T.to_numpy(dtype=float)
    scaler = StandardScaler()
    x = scaler.fit_transform(x)

    n_components = min(2, x.shape[0], x.shape[1], y.shape[1])
    if n_components < 2:
        raise ValueError("PLS-DA requires at least two samples and two features.")

    model = PLSRegression(n_components=2, scale=False)
    model.fit(x, y.to_numpy(dtype=float))

    scores = pd.DataFrame(
        model.x_scores_[:, :2],
        index=data.columns,
        columns=["PLS1", "PLS2"],
    )
    scores[cohort_col] = labels.reindex(scores.index).values

    mapped_scores = []
    for dataset, label in zip(datasets_to_map, labels_data_map):
        if set(dataset.index) != set(data.index):
            raise ValueError("Mapped datasets must have the same cell type labels.")
        mapped_data = dataset.reindex(data.index)

        mapped_x = scaler.transform(mapped_data.T.to_numpy(dtype=float))
        mapped_scores.append(
            pd.DataFrame(
                model.transform(mapped_x)[:, :2],
                index=mapped_data.columns,
                columns=["PLS1", "PLS2"],
            ).assign(**{cohort_col: label})
        )

    if mapped_scores:
        scores = pd.concat([scores, *mapped_scores])

    out_path = Path(out_path)
    out_path.parent.mkdir(parents=True, exist_ok=True)

    scores.to_csv(out_path.with_suffix(".csv"))

    vip = pd.Series(calculate_vip(model), index=data.index, name="VIP")

    loadings = pd.DataFrame(
        model.x_loadings_[:, :2], index=data.index, columns=["PLS1", "PLS2"]
    )

    plot_plsda_score(scores, out_path.with_suffix(".png"), cohort_col)
    plot_plsda_vip(vip, out_path.with_name(f"{out_path.stem}_vip.png"))
    plot_plsda_loading(loadings, out_path.with_name(f"{out_path.stem}_loading.png"))
    plot_plsda_biplot(
        scores, loadings, out_path.with_name(f"{out_path.stem}_biplot.png"), cohort_col
    )

    return scores

Cohort differences

hide_deconv.statistic.run_mann_whitney_u(bulks, sample_list, sample_id_col, cohort_col, celltypes_to_normalize_to=[])

Performs a Man Whitney U Test with FDR correction for two cohorts.

Parameters:

Name Type Description Default
bulks DataFrame

bulk Samples (either Celltype x Sample or Gene x Sample)

 required
sample_list DataFrame

sample list (samples x clinical variables)

 required
sample_id_col str

Name of the column, that links to the bulks

 required
cohort_col str

Name of the column containing the cohort identifiers.

 required
celltypes_to_normalize_to list = []

Cell types to exclude before sample wise renormalization.

 []

Returns:

Type Description
DataFrame

DataFrame containing the Mean, Standard Deviation, p-Value and adjusted p-Value
Source code in src/hide_deconv/statistic/mann_whitney_u.py 
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
def run_mann_whitney_u(
    bulks: pd.DataFrame,
    sample_list: pd.DataFrame,
    sample_id_col: str,
    cohort_col: str,
    celltypes_to_normalize_to: list = [],
) -> pd.DataFrame:
    """
    Performs a Man Whitney U Test with FDR correction for two cohorts.

    Parameters
    ----------
    bulks : pd.DataFrame
        bulk Samples (either Celltype x Sample or Gene x Sample)
    sample_list : pd.DataFrame
        sample list (samples x clinical variables)
    sample_id_col : str
        Name of the column, that links to the bulks
    cohort_col: str
        Name of the column containing the cohort identifiers.
    celltypes_to_normalize_to : list = []
        Cell types to exclude before sample wise renormalization.

    Returns
    -------
    pd.DataFrame
        DataFrame containing the Mean, Standard Deviation, p-Value and adjusted p-Value
    """

    samples = sample_list[[sample_id_col, cohort_col]].set_index(sample_id_col)

    # Subset samples, as sample list might have more samples, than in deconvolution
    samples = samples.reindex(bulks.columns)
    samples = samples.dropna(subset=[cohort_col])
    cohorts = samples[cohort_col].unique()
    bulks = bulks.loc[:, samples.index]

    if len(cohorts) > 2 or len(cohorts) < 2:
        raise Exception(
            "MWU can only be performed between two cohorts. Consider using a Kruskal Wallis Test."
        )

    group1, group2 = cohorts[0], cohorts[1]

    if len(celltypes_to_normalize_to) > 0:
        bulks = bulks.drop(index=celltypes_to_normalize_to).copy()
        sample_sums = bulks.sum(axis=0)
        bulks = bulks.div(sample_sums, axis=1)
    else:
        bulks = bulks.copy()

    pvals = []
    mean1, std1, mean2, std2 = [], [], [], []

    for ct in bulks.index:
        values = bulks.loc[ct]
        x = values[samples[cohort_col] == group1]
        y = values[samples[cohort_col] == group2]

        mean1.append(x.mean())
        mean2.append(y.mean())

        std1.append(x.std())
        std2.append(y.std())

        _, p = mannwhitneyu(x, y)
        pvals.append(p)

    pvals = np.array(pvals)

    pvals_adj = false_discovery_control(pvals)

    result = pd.DataFrame(
        {
            "celltype": bulks.index,
            f"mean[{group1}]": mean1,
            f"std[{group1}]": std1,
            f"mean[{group2}]": mean2,
            f"std[{group2}]": std2,
            "p": pvals,
            "p_adj": pvals_adj,
        }
    ).set_index("celltype")

    return result

hide_deconv.statistic.run_kruskal_wallis(bulks, sample_list, sample_id_col, cohort_col, celltypes_to_normalize_to=[])

Performs a Kruskal Wallis Test with FDR correction.

Parameters:

Name Type Description Default
bulks DataFrame

bulk Samples (either Celltype x Sample or Gene x Sample)

 required
sample_list DataFrame

sample list (samples x clinical variables)

 required
sample_id_col str

Name of the column, that links to the bulks

 required
cohort_col str

Name of the column containing the cohort identifiers.

 required
celltypes_to_normalize_to list = []

Cell types to exclude before sample wise renormalization.

 []

Returns:

Type Description
DataFrame

DataFrame containing p-Value and adjusted p-Value
Source code in src/hide_deconv/statistic/kruskal_wallis.py 
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
def run_kruskal_wallis(
    bulks: pd.DataFrame,
    sample_list: pd.DataFrame,
    sample_id_col: str,
    cohort_col: str,
    celltypes_to_normalize_to: list = [],
) -> pd.DataFrame:
    """
    Performs a Kruskal Wallis Test with FDR correction.

    Parameters
    ----------
    bulks : pd.DataFrame
        bulk Samples (either Celltype x Sample or Gene x Sample)
    sample_list : pd.DataFrame
        sample list (samples x clinical variables)
    sample_id_col : str
        Name of the column, that links to the bulks
    cohort_col: str
        Name of the column containing the cohort identifiers.
    celltypes_to_normalize_to : list = []
        Cell types to exclude before sample wise renormalization.

    Returns
    -------
    pd.DataFrame
        DataFrame containing p-Value and adjusted p-Value
    """

    samples = sample_list[[sample_id_col, cohort_col]].set_index(sample_id_col)

    # Subset samples, as sample list might have more samples, than in deconvolution
    samples = samples.reindex(bulks.columns)
    samples = samples.dropna(subset=[cohort_col])
    cohorts = samples[cohort_col].unique()
    bulks = bulks.loc[:, samples.index]

    fWarning = False
    for cohort in cohorts:
        if (samples[cohort_col] == cohort).sum() <= 1:
            console.print(f"[red]Cohort {cohort} contains only 1 sample.[/red]")
            fWarning = True
        elif (samples[cohort_col] == cohort).sum() < 5:
            console.print(
                f"[yellow]Cohort {cohort} contains less than 5 samples.[/yellow]"
            )
            fWarning = True
    if fWarning:
        console.print(
            "[dim]It is recommended, that each cohort contains at least 5 samples.[/dim]"
        )

    if len(celltypes_to_normalize_to) > 0:
        bulks = bulks.drop(index=celltypes_to_normalize_to).copy()
        sample_sums = bulks.sum(axis=0)
        bulks = bulks.div(sample_sums, axis=1)
    else:
        bulks = bulks.copy()

    pvals = []
    valid_celltypes = []

    for ct in bulks.index:
        values = bulks.loc[ct]
        data = []

        for cohort in cohorts:
            cohort_values = values[samples[cohort_col] == cohort].dropna()
            data.append(cohort_values)

        if any(len(group) == 0 for group in data):
            console.print(
                f"[yellow]Celltype {ct} cannot be tested because at least one cohort is empty.[/yellow]"
            )
            continue

        all_values = pd.concat(data, ignore_index=True)
        if all_values.nunique(dropna=True) < 2:
            console.print(
                f"[yellow]Celltype {ct} cannot be tested because all cohort values are identical.[/yellow]"
            )
            continue

        try:
            _, p = kruskal(*data)
        except ValueError as exc:
            console.print(
                f"[yellow]Celltype {ct} could not be tested by Kruskal-Wallis and was skipped: {exc}[/yellow]"
            )
            continue

        if not np.isfinite(p):
            console.print(
                f"[yellow]Celltype {ct} produced an invalid p-value and was skipped.[/yellow]"
            )
            continue

        valid_celltypes.append(ct)
        pvals.append(p)

    if len(pvals) == 0:
        console.print(
            "[yellow]No testable celltypes found for Kruskal-Wallis analysis.[/yellow]"
        )
        return pd.DataFrame(columns=["p", "p_adj"]).rename_axis("celltype")

    pvals = np.array(pvals, dtype=float)
    pvals_adj = false_discovery_control(pvals)

    result = pd.DataFrame(
        {
            "celltype": valid_celltypes,
            "p": pvals,
            "p_adj": pvals_adj,
        }
    ).set_index("celltype")

    return result

hide_deconv.statistic.run_dunn(kruskal_results, bulks, sample_list, sample_id_col, cohort_col, sign_level=0.05)

Performs a Posthoc Dunn Test on significant Kruskal Wallis results.

Parameters:

Name Type Description Default
kruskal_results DataFrame

Results of obtained by the run_kruskal_wallis() method

 required
bulks DataFrame

bulk Samples (either Celltype x Sample or Gene x Sample)

 required
sample_list DataFrame

sample list (samples x clinical variables)

 required
sample_id_col str

Name of the column, that links to the bulks

 required
cohort_col str

Name of the column containing the cohort identifiers.

 required
sign_level float = 0.05

Float, below which results are considered as significant.

 0.05

Returns:

Type Description
DataFrame

DataFrame containing the results of the Dunn Test.
Source code in src/hide_deconv/statistic/posthoc_dunn.py 
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
def run_dunn(
    kruskal_results: pd.DataFrame,
    bulks: pd.DataFrame,
    sample_list: pd.DataFrame,
    sample_id_col: str,
    cohort_col: str,
    sign_level: float = 0.05,
) -> pd.DataFrame:
    """
    Performs a Posthoc Dunn Test on significant Kruskal Wallis results.

    Parameters
    ----------
    kruskal_results : pd.DataFrame
        Results of obtained by the run_kruskal_wallis() method
    bulks : pd.DataFrame
        bulk Samples (either Celltype x Sample or Gene x Sample)
    sample_list : pd.DataFrame
        sample list (samples x clinical variables)
    sample_id_col : str
        Name of the column, that links to the bulks
    cohort_col: str
        Name of the column containing the cohort identifiers.
    sign_level : float = 0.05
        Float, below which results are considered as significant.

    Returns
    -------
    pd.DataFrame
        DataFrame containing the results of the Dunn Test.
    """

    samples = sample_list[[sample_id_col, cohort_col]].set_index(sample_id_col)
    samples = samples.reindex(bulks.columns)

    # get significant results
    sig_celltypes = kruskal_results.index[kruskal_results["p_adj"] < sign_level]

    rows = []
    for ct in sig_celltypes:
        values = bulks.loc[ct]

        df = pd.DataFrame(
            {
                "value": values.values,
                "cohort": samples.loc[values.index, cohort_col].values,
            },
            index=values.index,
        ).dropna()

        dunn = sp.posthoc_dunn(
            df, val_col="value", group_col="cohort", p_adjust="fdr_bh"
        )

        dunn.index.name = "cohort_1"
        dunn.columns.name = "cohort_2"

        cohort_means = df.groupby("cohort")["value"].mean()

        long_df = (
            dunn.where(np.triu(np.ones(dunn.shape, dtype=bool), k=1))
            .stack()
            .reset_index(name="p_adj")
        )

        long_df["mean[cohort_1]"] = long_df["cohort_1"].map(cohort_means)
        long_df["mean[cohort_2]"] = long_df["cohort_2"].map(cohort_means)

        long_df.insert(0, "celltype", ct)
        rows.append(long_df)

    if len(rows) == 0:
        return pd.DataFrame(
            columns=[
                "celltype",
                "cohort_1",
                "cohort_2",
                "p_adj",
                "mean[cohort_1]",
                "mean[cohort_2]",
            ]
        )
    return pd.concat(rows, ignore_index=True)

Clustering

hide_deconv.statistic.run_clustering(data, is_bulk=False)

Performs a clustering using greedy modular communities of the entered bulk or composition data.

Parameters:

Name Type Description Default
data DataFrame

Dataframe containing either bulks or composition (celltype/genes x samples)

 required
is_bulk bool = False

If set to true, uses correlation as distance measure for the neighborhood calculation.

 False

Returns:

Type Description
DataFrame

DataFrame containing the sample ids and the assigned clusters
Source code in src/hide_deconv/statistic/clustering.py 
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
def run_clustering(data: pd.DataFrame, is_bulk: bool = False) -> pd.DataFrame:
    """
    Performs a clustering using greedy modular communities of the entered bulk or composition data.

    Parameters
    ----------
    data : pd.DataFrame
        Dataframe containing either bulks or composition (celltype/genes x samples)
    is_bulk : bool = False
        If set to true, uses correlation as distance measure for the neighborhood calculation.

    Returns
    -------
    pd.DataFrame
        DataFrame containing the sample ids and the assigned clusters

    """

    obs = pd.DataFrame(index=data.columns)
    obs["id"] = data.columns
    var = pd.DataFrame(index=data.index)
    adata = sp.AnnData(X=data.T.to_numpy(), obs=obs, var=var)

    if is_bulk:
        sp.pp.neighbors(adata, metric="correlation")
    else:
        sp.pp.neighbors(
            adata, metric="braycurtis"
        )  # braycurtis seems more approriate here

    # Sparse arrays have distances instead of connectivities property
    if "connectivities" in adata.obsp:
        conn = adata.obsp["connectivities"]
    elif "distances" in adata.obsp:
        conn = adata.obsp["distances"]
    else:
        raise ValueError("Calculation of neighborhood graph has failed!")

    if not sps.issparse(conn):
        conn = sps.csr_matrix(conn)

    coo = conn.tocoo()
    G = nx.Graph()
    for u, v, w in zip(coo.row, coo.col, coo.data):
        if w != 0:
            G.add_edge(int(u), int(v), weight=float(w))

    if G.number_of_nodes() == 0:
        labels = []
    else:
        communities = nx.algorithms.community.greedy_modularity_communities(G)
        labels = [None] * G.number_of_nodes()
        for i, comm in enumerate(communities):
            for node in comm:
                labels[node] = i

    return pd.DataFrame({"id": list(adata.obs_names), "assigned_cluster": labels})

Survival analysis

hide_deconv.statistic.run_cox_regression(bulks, sample_sheet, sample_id_col, time_col, event_col, covariates)

Perform Cox Regression for each cell type composition.

Parameters:

Name Type Description Default
bulks DataFrame

Cell type compositions (either Celltype x Sample or Gene x Sample)

 required
sample_sheet DataFrame

sample list (samples x clinical variables)

 required
sample_id_col str

Name of the column, that links to the bulks

 required
time_col str

Name of the column with survival times

 required
event_col str

Name of column that indicates event (0: censored, 1: event)

 required
covariates list[str]

List of covariate column names

 required

Returns:

Type Description
DataFrame

DataFrame with estimated impact on survival for each cell type
Source code in src/hide_deconv/statistic/survival_analysis.py 
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
def run_cox_regression(
    bulks: pd.DataFrame,
    sample_sheet: pd.DataFrame,
    sample_id_col: str,
    time_col: str,
    event_col: str,
    covariates: list[str] | None,
) -> pd.DataFrame:
    """
    Perform Cox Regression for each cell type composition.

    Parameters
    ----------
    bulks : pd.DataFrame
        Cell type compositions (either Celltype x Sample or Gene x Sample)
    sample_sheet : pd.DataFrame
        sample list (samples x clinical variables)
    sample_id_col : str
        Name of the column, that links to the bulks
    time_col : str
        Name of the column with survival times
    event_col : str
        Name of column that indicates event (0: censored, 1: event)
    covariates : list[str], optional
        List of covariate column names

    Returns
    -------
    pd.DataFrame
        DataFrame with estimated impact on survival for each cell type
    """

    samples = sample_sheet.set_index(sample_id_col)
    samples = samples.reindex(bulks.columns)

    # Remove entries with missing time or event
    non_none_samples = samples[[time_col, event_col]].notna().all(axis=1)

    # Remove entries with missing covariates, if given
    if covariates:
        non_none_samples = non_none_samples & samples.loc[:, covariates].notna().all(
            axis=1
        )

    # Warning if samples have been removed
    n_removed = (~non_none_samples).sum()
    if n_removed > 0:
        console.print(
            f"[yellow]Removed {n_removed} samples due to missing values.[/yellow]"
        )

    samples = samples.loc[non_none_samples]
    bulks = bulks.loc[:, non_none_samples]

    results = []
    pvals = []

    with console.status(
        "[bold blue]Fitting CoxPH-models...[/bold blue]",
        spinner="dots",
    ):
        for ct in bulks.index:
            try:
                X = pd.DataFrame({"ct_comp": bulks.loc[ct]}, index=bulks.columns)

                X[time_col] = samples.loc[X.index, time_col].values
                X[event_col] = samples.loc[X.index, event_col].values

                if covariates:
                    for cov in covariates:
                        cov_data = samples.loc[X.index, cov].values

                        # Check if covariate is numeric or categorical
                        if pd.api.types.is_numeric_dtype(cov_data):
                            X[cov] = cov_data
                        else:
                            # One-hot encode categorical variables
                            cov_dummies = pd.get_dummies(
                                samples.loc[X.index, cov],
                                prefix=cov,
                                drop_first=True,  # Avoid multicollinearity
                            )
                            X = X.join(cov_dummies)

                X = X.loc[samples.index]

                # Train CoxPH model
                cph = CoxPHFitter()
                cph.fit(X, duration_col=time_col, event_col=event_col)

                coef = cph.params_.loc["ct_comp"]
                hr = np.exp(coef)
                ci = cph.confidence_intervals_.loc["ct_comp"]
                p_val = cph.summary.loc["ct_comp", "p"]

                concordance = cph.concordance_index_

                results.append(
                    {
                        "celltype": ct,
                        "coef": coef,
                        "hr": hr,
                        "ci_lower": np.exp(ci.iloc[0]),
                        "ci_upper": np.exp(ci.iloc[1]),
                        "p_value": p_val,
                        "concordance_index": concordance,
                    }
                )
                pvals.append(p_val)

            except Exception:
                console.print(f"[yellow]Failed to fit CoxPH model for {ct}.[/yellow]")
                results.append(
                    {
                        "celltype": ct,
                        "coef": np.nan,
                        "hr": np.nan,
                        "ci_lower": np.nan,
                        "ci_upper": np.nan,
                        "p_value": np.nan,
                        "concordance_index": np.nan,
                    }
                )
                pvals.append(np.nan)

        result_df = pd.DataFrame(results)
        p_vals = np.array(pvals)
        p_vals_valid = ~np.isnan(p_vals)

        pvals_adj = np.full_like(p_vals, np.nan, dtype=float)
        pvals_adj[p_vals_valid] = false_discovery_control(p_vals[p_vals_valid])

        result_df["p_value_adj"] = pvals_adj

        return result_df

Differential expression

hide_deconv.statistic.pydeseq2_preprocess(bulk, sample_sheet, sample_id_col, condition_col, covariates)

Prepare bulk and sample metainfo for PyDESeq2.

Parameters:

Name Type Description Default
bulk DataFrame

Bulk RNA-seq file (genes x samples)

 required
sample_sheet DataFrame

Sample sheet containing the metainformation on samples (samples x info)

 required
sample_id_col str

Column name of the sample sheet containing the sample ids of the bulk file

 required
condition_col str

Column name of the sample sheet containing the conditions

 required
covariates list[str] | None

List of column names containing possible covariates

 required

Returns:

Type Description
tuple[DataFrame, DataFrame]

Counts, Metainformation to be used with the run_pydeseq2 function
Source code in src/hide_deconv/statistic/pydeseq2.py 
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
def pydeseq2_preprocess(
    bulk: pd.DataFrame,
    sample_sheet: pd.DataFrame,
    sample_id_col: str,
    condition_col: str,
    covariates: list[str] | None,
) -> tuple[pd.DataFrame, pd.DataFrame]:
    """
    Prepare bulk and sample metainfo for PyDESeq2.

    Parameters
    ----------
    bulk : pd.DataFrame
        Bulk RNA-seq file (genes x samples)
    sample_sheet : pd.DataFrame
        Sample sheet containing the metainformation on samples (samples x info)
    sample_id_col : str
        Column name of the sample sheet containing the sample ids of the bulk file
    condition_col : str
        Column name of the sample sheet containing the conditions
    covariates : list[str] | None
        List of column names containing possible covariates

    Returns
    -------
    tuple[pd.DataFrame, pd.DataFrame]
        Counts, Metainformation to be used with the run_pydeseq2 function
    """

    sample_cols = [sample_id_col, condition_col]
    if covariates:
        sample_cols.extend(covariates)

    samples = sample_sheet[sample_cols].copy()
    samples = samples.dropna(subset=[sample_id_col, condition_col])

    if covariates:
        samples = samples.dropna(subset=covariates)

    samples[sample_id_col] = samples[sample_id_col].astype(str)
    samples = samples.drop_duplicates(subset=[sample_id_col], keep="first")
    bulk = bulk.copy()
    bulk.columns = bulk.columns.astype(str)
    samples = samples[samples[sample_id_col].isin(bulk.columns)]

    if len(samples) == 0:
        raise ValueError("No matching sample ids were found in the bulk file.")

    samples = samples.set_index(sample_id_col).reindex(bulk.columns)
    samples = samples.dropna(subset=[condition_col])

    if covariates:
        samples = samples.dropna(subset=covariates)

    counts = bulk.loc[:, samples.index].T.copy()
    counts.index = counts.index.astype(str)
    counts.columns = counts.columns.astype(str)

    metadata = samples.copy()
    metadata.index = metadata.index.astype(str)
    metadata[condition_col] = metadata[condition_col].astype(str)

    return counts, metadata

hide_deconv.statistic.run_pydeseq2(bulk, metadata, condition_col, tested_condition, reference_condition, covariates, out_path)

Run PyDESeq2

Parameters:

Name Type Description Default
bulk DataFrame

Raw count bulk RNA-seq data (Samples x Genes)

 required
metadata DataFrame

Sample metainformation

 required
condition_col str

Name of the column in metadata holding the condition

 required
tested_condition str

Name of the condition that will be tested

 required
reference_condition str

Name of the condition that is used as reference

 required
covariates list[str]

Column names of covariates to include

 required
out_path Folder path, where all results and plots are stored
required

Returns:

Type Description
DataFrame

pydeseq2 result dataframe
Source code in src/hide_deconv/statistic/pydeseq2.py 
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
def run_pydeseq2(
    bulk: pd.DataFrame,
    metadata: pd.DataFrame,
    condition_col: str,
    tested_condition: str,
    reference_condition: str,
    covariates: list[str] | None,
    out_path: Path,
) -> pd.DataFrame:
    """
    Run PyDESeq2

    Parameters
    ----------
    bulk : pd.DataFrame
        Raw count bulk RNA-seq data (Samples x Genes)
    metadata : pd.DataFrame
        Sample metainformation
    condition_col : str
        Name of the column in metadata holding the condition
    tested_condition : str
        Name of the condition that will be tested
    reference_condition : str
        Name of the condition that is used as reference
    covariates : list[str]
        Column names of covariates to include
    out_path : Folder path, where all results and plots are stored

    Returns
    -------
    pd.DataFrame
        pydeseq2 result dataframe
    """

    design_terms = [condition_col]
    if covariates:
        design_terms.extend(covariates)

    dds = DeseqDataSet(
        counts=bulk,
        metadata=metadata,
        design="~ " + " + ".join(design_terms),
        quiet=True,
    )
    dds.deseq2()

    deseq_stats = DeseqStats(
        dds,
        contrast=[condition_col, tested_condition, reference_condition],
        quiet=True,
    )
    deseq_stats.summary()

    results = deseq_stats.results_df.copy()

    out_path = Path(out_path)
    out_path.parent.mkdir(parents=True, exist_ok=True)

    results.to_csv(out_path.with_name(f"{out_path.name}_results.csv"))
    deseq_stats.plot_MA(save_path=str(out_path.with_name(f"{out_path.name}_ma.png")))
    plot_volcano(results, out_path.with_name(f"{out_path.name}_volcano.png"))

    return results