In this Notebook we use scCNR to reconstruct the same signaling network of the scMRA_Tutorial comparatively between two cell populations. Next to the wildtype populations that we analyzed in the scMRA_Tutorial we here have as well a simulated BRAF cell population. ScCNR tries to not only identify the network topology, but also identify the main differences in signaling between the wildtype and BRAF mutant population.
import scmra
import pandas as pd
import numpy as np
import seaborn as sns; sns.set_theme()
import matplotlib.pyplot as plt
from IPython.display import display, HTML
import networkx as nx
#read in data
dataWT = pd.read_csv("../data/scData_wildtype.tsv", sep='\t')
dataBRAF = pd.read_csv("../data/scData_BRAFmut.tsv", sep='\t')
Similar to scMRA we need to calculate deviations of single cells to the population mean. Since we now consider several populations of single cells, we additionally need to assign group labels to map a cellID to its population.
#add a SC-ID
dataWT.index = ["CTR."+str(i) for i in range(dataWT.shape[0])]
dataBRAF.index = ["BRAFmut."+str(i) for i in range(dataBRAF.shape[0])]
#calculate input: deviatipon of protein from cell population mean
populations = [dataWT, dataBRAF]
#later we also have to access the populations specific rloc values by these labels
labels = ["Ctr", "BRAFmut"]
dataScaled, groupAnnot = scmra.prepare_data(populations, labels)
#get active protein deviations (frmo those we infer interaction strengths)
rglob = dataScaled.loc[[n for n in list(dataScaled.index) if ("Active" in n)]]
rglob.index = [i.replace("Active", "") for i in rglob.index]
#get total protein deviations: 'natural perturbations'
rtot = dataScaled.loc[[n for n in list(dataScaled.index) if ("Tot" in n)]]
rtot.index = [i.replace("Tot", "") for i in rtot.index]
display(rglob.head)
<bound method NDFrame.head of CTR.0 CTR.1 CTR.2 CTR.3 CTR.4 CTR.5 CTR.6 \ Sos -0.051155 0.085232 0.008321 0.073530 0.037336 0.077462 0.092302 Ras -0.076395 0.054095 -0.007728 0.025145 0.028869 0.058157 0.088900 Raf1 -0.103680 0.069777 0.053183 0.041268 0.066881 0.064294 0.005852 Mek -0.033101 -0.152314 -0.050093 -0.128995 0.069491 -0.148370 0.166537 Erk -0.106257 -0.177075 0.030958 -0.113891 0.051923 -0.161771 0.092711 P90Rsk -0.129504 -0.117604 -0.007665 -0.074405 -0.126830 -0.019424 0.120107 PI3K 0.142959 0.021876 -0.143707 0.023346 0.029151 0.198616 0.008005 Akt 0.116549 -0.098978 -0.143954 -0.084988 0.046098 0.104219 0.110371 C3G 0.007717 -0.082415 -0.006071 -0.138128 0.056060 -0.083842 -0.024886 Rap1 -0.005878 -0.069144 0.016994 -0.137141 0.027891 -0.114642 0.005731 bRaf 0.072555 -0.084724 -0.030218 -0.107014 0.048830 -0.116599 0.097450 CTR.7 CTR.8 CTR.9 ... BRAFmut.240 BRAFmut.241 \ Sos -0.008603 -0.222532 0.126810 ... -0.146622 -0.093920 Ras 0.004750 -0.215548 0.117134 ... -0.129316 -0.127258 Raf1 0.066933 -0.259652 0.135834 ... -0.123249 -0.081300 Mek 0.058568 0.107128 -0.142194 ... 0.111068 -0.166231 Erk 0.024951 0.221935 -0.127559 ... 0.102615 0.071023 P90Rsk 0.015984 0.237644 -0.188018 ... 0.238576 -0.028478 PI3K 0.087303 -0.074450 -0.020944 ... 0.052805 0.006293 Akt -0.176313 0.044074 -0.082037 ... 0.036224 -0.175243 C3G 0.005292 0.115988 0.113998 ... 0.052111 0.107812 Rap1 0.009414 0.123415 0.086880 ... 0.055672 0.117495 bRaf 0.006757 0.086650 0.085719 ... 0.031120 0.183854 BRAFmut.242 BRAFmut.243 BRAFmut.244 BRAFmut.245 BRAFmut.246 \ Sos 0.041602 0.083655 0.090291 -0.179314 0.160474 Ras 0.041264 0.118002 0.074339 -0.171745 0.168231 Raf1 0.017540 0.185428 0.098197 -0.200488 0.163796 Mek -0.066659 -0.113333 -0.065334 0.011382 0.044878 Erk 0.223207 -0.092159 -0.069101 -0.063405 -0.034418 P90Rsk 0.026048 -0.060250 -0.011676 0.124937 0.005208 PI3K 0.030364 -0.048787 -0.154532 -0.080240 0.048716 Akt 0.104959 -0.075421 -0.116958 0.152433 -0.011654 C3G 0.138902 0.092160 -0.031306 -0.015758 -0.005661 Rap1 0.164348 0.098692 -0.030684 0.013867 -0.006914 bRaf -0.019847 0.222383 -0.149625 -0.095272 0.202139 BRAFmut.247 BRAFmut.248 BRAFmut.249 Sos 0.366688 0.029405 -0.023816 Ras 0.407325 0.011223 -0.040518 Raf1 0.496338 -0.056441 -0.009805 Mek -0.043485 -0.000050 0.031329 Erk -0.059881 -0.025874 0.281526 P90Rsk -0.180146 0.017860 0.090534 PI3K -0.188240 -0.041037 -0.089324 Akt 0.003887 0.010586 -0.002628 C3G 0.099245 -0.083501 -0.032298 Rap1 0.110646 -0.056336 -0.041837 bRaf 0.052489 -0.135774 0.252340 [11 rows x 500 columns]>
scCNR penalizes the number of non-equal edges in the network(additional to the number of edges in the network). Similar to the estimation of a suitable η-parameter, the θ-penalty on non-equal edges can also be estimated (see scMRA_Tutorial). Here we continue with parameters that we estimated for the method evaluation in the paper.
The result of the scCNR ananlysis return not only one weighted adjacency matrix of the inferred network, but one matrix for each cell population. Depending on the θ-penalty we force edges to be similar between networks, and would ideally retain 'rloc' matrices, that are fairly similar and differ in the main interactions that explain hte main signaling differences between the populations. Here we display the differences between the wildtype and BRAFmut population and see that most differences are detected downstream of BRAF.
%%capture
#make MRA simulation
eta = 0.02
theta = 0.1
scd = scmra.ScData(rglob=rglob, rtot = rtot, group_annot= groupAnnot)
p = scmra.ScCnrProblem(scd, eta=eta, theta=theta)
p.cpx.solve()
result = scmra.ScCnrResult(p)
#display difference between both adjacency matrices
display(result.rloc["Ctr"] - result.rloc["BRAFmut"])
differenceRloc = result.rloc["Ctr"] - result.rloc["BRAFmut"]
#diagonal rloc entries are per definition -1, replace those with 0 to plot network
A = differenceRloc.transpose()
A[A == -1] = 0
G = nx.from_pandas_adjacency(A, create_using=nx.DiGraph)
pos=nx.shell_layout(G)
nx.draw_networkx(G, with_labels=True, node_size = 750, pos = pos, node_color = "yellow")
#create edge labels
edgeLabels = {}
edgeColor = {}
for row_label, row_data in A.iterrows():
for col_label, value in row_data.items():
key = (col_label, row_label)
if(np.abs(value) > 0):
edgeLabels[key] = round(value,1)
nx.draw_networkx_edge_labels(
G, pos,
edge_labels=edgeLabels,
font_color="black"
)
plt.show()
Sos | Ras | Raf1 | Mek | Erk | P90Rsk | PI3K | Akt | C3G | Rap1 | bRaf | |
---|---|---|---|---|---|---|---|---|---|---|---|
Sos | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.000000 |
Ras | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.000000 |
Raf1 | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.000000 |
Mek | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.791585 |
Erk | 0.0 | 0.0 | 0.0 | 0.965511 | 0.000000 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.000000 |
P90Rsk | 0.0 | 0.0 | 0.0 | 0.000000 | 0.701235 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.000000 |
PI3K | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.000000 |
Akt | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.000000 |
C3G | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.000000 |
Rap1 | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.000000 |
bRaf | 0.0 | 0.0 | 0.0 | 0.000000 | 0.000000 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.000000 |