Na Prática: NetworkX + Banco de Dados STRING

pip install networkx matplotlib pandas requests python-louvain

import networkx as nx
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.cm as cm
import requests
import json
from collections import Counter
import community as community_louvain  # python-louvain

def get_string_interactions(proteins: list, species: int = 9606,
                             min_score: int = 400) -> pd.DataFrame:
    """
    Busca interações de proteínas do STRING.
    
    species: 9606 = Homo sapiens, 10090 = Mus musculus, 4932 = S. cerevisiae
    min_score: pontuação mínima de confiança (0-1000)
    """
    url = "https://string-db.org/api/json/network"
    params = {
        "identifiers": "\r".join(proteins),
        "species": species,
        "required_score": min_score,
        "caller_identity": "bio_for_devs_tutorial"
    }
    
    response = requests.post(url, data=params)
    response.raise_for_status()
    
    data = response.json()
    
    if not data:
        return pd.DataFrame()
    
    df = pd.DataFrame(data)
    return df[["preferredName_A", "preferredName_B", "score"]]

# Rede de resposta a dano de DNA: genes chave
dna_repair_genes = [
    "TP53", "BRCA1", "BRCA2", "ATM", "ATR", "CHEK1", "CHEK2",
    "RAD51", "RPA1", "PCNA", "PARP1", "MDM2", "CDKN1A",
    "H2AFX", "NBN", "MRE11", "RAD50", "FANCD2", "XRCC1"
]

interactions = get_string_interactions(dna_repair_genes, min_score=400)
print(f"Interações encontradas: {len(interactions)}")
print(interactions.head(10))

def build_protein_network(interactions_df: pd.DataFrame,
                           min_score: float = 0.4) -> nx.Graph:
    """
    Constrói um grafo de rede proteica a partir de interações STRING.
    
    A pontuação STRING (0-1000) é convertida para (0-1) como peso da borda.
    """
    G = nx.Graph()
    
    for _, row in interactions_df.iterrows():
        weight = row["score"] / 1000.0
        if weight >= min_score:
            G.add_edge(
                row["preferredName_A"],
                row["preferredName_B"],
                weight=weight
            )
    
    print(f"Nós: {G.number_of_nodes()}")
    print(f"Bordas: {G.number_of_edges()}")
    print(f"Componente conectado: {nx.is_connected(G)}")
    
    # Componente maior
    if not nx.is_connected(G):
        largest_cc = max(nx.connected_components(G), key=len)
        G = G.subgraph(largest_cc).copy()
        print(f"Componente maior: {G.number_of_nodes()} nós")
    
    return G

G = build_protein_network(interactions, min_score=0.4)

def compute_centrality_metrics(G: nx.Graph) -> pd.DataFrame:
    """
    Calcula múltiplas métricas de centralidade para identificar proteínas hub.
    
    Diferentes métricas capturam diferentes noções de "importância":
    - Grau: quantas conexões diretas
    - Intermediação: quão frequentemente na rota mais curta entre outros
    - Proximidade: quão rapidamente alcança todos os outros
    - PageRank: importância recursiva (algoritmo original do Google)
    """
    metrics = pd.DataFrame(index=G.nodes())
    
    # Centralidade de grau: fração de nós conectados
    degree_centrality = nx.degree_centrality(G)
    metrics["degree"] = pd.Series(G.degree()).values
    metrics["degree_centrality"] = pd.Series(degree_centrality)
    
    # Centralidade de intermediação: fração de caminhos mais curtos que passam por aqui
    betweenness = nx.betweenness_centrality(G, weight="weight")
    metrics["betweenness"] = pd.Series(betweenness)
    
    # Centralidade de proximidade: recíproca da soma das distâncias
    closeness = nx.closeness_centrality(G)
    metrics["closeness"] = pd.Series(closeness)
    
    # PageRank: importância com base em conexões importantes
    pagerank = nx.pagerank(G, weight="weight")
    metrics["pagerank"] = pd.Series(pagerank)
    
    metrics = metrics.sort_values("degree", ascending=False)
    return metrics

centrality = compute_centrality_metrics(G)
print("\nTop 10 proteínas por grau:")
print(centrality.head(10)[["degree", "betweenness", "closeness", "pagerank"]].round(4))

def detect_communities(G: nx.Graph) -> dict:
    """
    Detecta comunidades funcionais usando algoritmo Louvain.
    
    Comunidades = grupos de proteínas densamente conectadas entre si.
    Em redes biológicas, frequentemente correspondem a complexos proteicos
    ou módulos funcionais.
    """
    # Algoritmo Louvain: otimiza modularidade
    partition = community_louvain.best_partition(G, weight="weight")
    
    # Contar tamanho das comunidades
    community_sizes = Counter(partition.values())
    print(f"\nNúmero de comunidades: {len(community_sizes)}")
    
    # Listar proteínas em cada comunidade
    communities = {}
    for node, comm_id in partition.items():
        if comm_id not in communities:
            communities[comm_id] = []
        communities[comm_id].append(node)
    
    print("\nComunidades (ordenadas por tamanho):")
    for comm_id, members in sorted(communities.items(), 
                                    key=lambda x: -len(x[1])):
        print(f"  Comunidade {comm_id} ({len(members)} proteínas): {', '.join(sorted(members))}")
    
    return partition, communities

partition, communities = detect_communities(G)

def visualize_network(G: nx.Graph, partition: dict, 
                       centrality: pd.DataFrame,
                       title: str = "Rede de Interação Proteica") -> None:
    """
    Visualiza a rede com:
    - Cor dos nós = comunidade
    - Tamanho dos nós = grau (proteínas hub são maiores)
    - Espessura da borda = peso da interação (pontuação STRING)
    """
    fig, ax = plt.subplots(figsize=(14, 10))
    
    # Layout: spring layout (força dirigida)
    pos = nx.spring_layout(G, k=2, seed=42)
    
    # Propriedades dos nós
    n_communities = len(set(partition.values()))
    cmap = cm.get_cmap("tab20", n_communities)
    node_colors = [cmap(partition.get(node, 0)) for node in G.nodes()]
    node_sizes = [300 + centrality.loc[node, "degree"] * 50 
                  if node in centrality.index else 300 
                  for node in G.nodes()]
    
    # Propriedades das bordas
    edge_weights = [G[u][v].get("weight", 0.5) for u, v in G.edges()]
    
    # Desenhar
    nx.draw_networkx_edges(G, pos, alpha=0.3, 
                           width=[w * 3 for w in edge_weights],
                           edge_color="gray", ax=ax)
    nx.draw_networkx_nodes(G, pos, node_color=node_colors,
                           node_size=node_sizes, alpha=0.9, ax=ax)
    
    # Rotular apenas hubs (proteínas com grau >= 5)
    hub_labels = {node: node for node in G.nodes() 
                  if G.degree(node) >= 5}
    nx.draw_networkx_labels(G, pos, labels=hub_labels,
                            font_size=8, font_weight="bold", ax=ax)
    
    ax.set_title(title, fontsize=14, fontweight="bold")
    ax.axis("off")
    
    plt.tight_layout()
    plt.savefig("protein_network.png", dpi=150, bbox_inches="tight")
    plt.show()

visualize_network(G, partition, centrality, 
                  "Rede de Resposta a Dano de DNA")

def analyze_paths(G: nx.Graph, source: str, targets: list) -> None:
    """
    Analisa caminhos mais curtos de uma proteína de interesse a alvos.
    
    Útil para entender cadeias de sinalização: como uma perturbação
    em 'source' propaga para 'targets' através da rede.
    """
    print(f"\nCaminhos de {source}:")
    
    for target in targets:
        if target not in G.nodes():
            print(f"  {target}: não no grafo")
            continue
        
        try:
            path = nx.shortest_path(G, source, target)
            length = nx.shortest_path_length(G, source, target)
            print(f"  → {target}: {' → '.join(path)} (comprimento {length})")
        except nx.NetworkXNoPath:
            print(f"  → {target}: sem caminho")

if "TP53" in G.nodes():
    analyze_paths(G, "TP53", ["BRCA1", "RAD51", "CHEK1", "MDM2"])

def get_string_enrichment(proteins: list, species: int = 9606) -> pd.DataFrame:
    """
    Obtém enriquecimento funcional para uma lista de proteínas via STRING.
    
    Retorna termos GO, vias KEGG/Reactome enriquecidos.
    """
    url = "https://string-db.org/api/json/enrichment"
    params = {
        "identifiers": "\r".join(proteins),
        "species": species,
        "caller_identity": "bio_for_devs_tutorial"
    }
    
    response = requests.post(url, data=params)
    response.raise_for_status()
    
    data = response.json()
    if not data:
        return pd.DataFrame()
    
    df = pd.DataFrame(data)
    
    # Filtrar e ordenar por valor-p
    result = df[["category", "term", "description", "fdr", "preferredNames"]].copy()
    result["n_proteins"] = result["preferredNames"].apply(
        lambda x: len(x) if isinstance(x, list) else 0
    )
    result = result.sort_values("fdr").head(20)
    
    return result

# Enriquecimento para a comunidade maior
largest_community = max(communities.values(), key=len)
enrichment = get_string_enrichment(largest_community)

print("\nEnriquecimento funcional da maior comunidade:")
if not enrichment.empty:
    print(enrichment[["category", "description", "fdr", "n_proteins"]]
          .to_string(index=False))

def build_differential_network(normal_interactions: pd.DataFrame,
                                disease_interactions: pd.DataFrame,
                                min_score: float = 0.4) -> dict:
    """
    Compara redes entre dois estados (normal vs doença).
    
    Identifica:
    - Interações ganhas na doença
    - Interações perdidas na doença
    - Nós com grau alterado significativamente
    """
    G_normal = build_protein_network(normal_interactions, min_score)
    G_disease = build_protein_network(disease_interactions, min_score)
    
    normal_edges = set(G_normal.edges())
    disease_edges = set(G_disease.edges())
    
    # Normalizar (bordas não direcionadas podem ter ordem revertida)
    normal_edges_norm = {(min(u,v), max(u,v)) for u,v in normal_edges}
    disease_edges_norm = {(min(u,v), max(u,v)) for u,v in disease_edges}
    
    gained = disease_edges_norm - normal_edges_norm
    lost = normal_edges_norm - disease_edges_norm
    
    print(f"Interações ganhas na doença: {len(gained)}")
    print(f"Interações perdidas na doença: {len(lost)}")
    
    # Nós com maior mudança de grau
    degree_change = {}
    all_nodes = set(G_normal.nodes()) | set(G_disease.nodes())
    for node in all_nodes:
        normal_deg = G_normal.degree(node) if node in G_normal else 0
        disease_deg = G_disease.degree(node) if node in G_disease else 0
        degree_change[node] = disease_deg - normal_deg
    
    df = pd.DataFrame.from_dict(degree_change, orient="index", 
                                 columns=["degree_change"])
    df = df.sort_values("degree_change")
    
    print("\nProteínas com maior perda de grau na doença:")
    print(df.head(5))
    print("\nProteínas com maior ganho de grau na doença:")
    print(df.tail(5))
    
    return {"gained": gained, "lost": lost, "degree_change": df}

def full_ppi_analysis(gene_list: list, species: int = 9606,
                       min_score: int = 400,
                       output_prefix: str = "network_analysis"):
    """Pipeline completo: STRING → rede → análise → visualização"""
    print(f"Analisando {len(gene_list)} proteínas...")
    
    # 1. Buscar interações
    interactions = get_string_interactions(gene_list, species, min_score)
    
    if interactions.empty:
        print("Nenhuma interação encontrada. Tente reduzir min_score.")
        return
    
    # 2. Construir rede
    G = build_protein_network(interactions, min_score / 1000)
    
    # 3. Centralidade
    centrality = compute_centrality_metrics(G)
    centrality.to_csv(f"{output_prefix}_centrality.tsv", sep="\t")
    
    # 4. Comunidades
    partition, communities = detect_communities(G)
    
    # 5. Visualizar
    visualize_network(G, partition, centrality, 
                      f"Rede PPI: {', '.join(gene_list[:3])}...")
    
    # 6. Enriquecimento
    enrichment = get_string_enrichment(gene_list, species)
    if not enrichment.empty:
        enrichment.to_csv(f"{output_prefix}_enrichment.tsv", sep="\t", index=False)
    
    print(f"\nArquivos salvos com prefixo: {output_prefix}")
    return G, centrality, partition, communities

# Executar pipeline completo
G, centrality, partition, communities = full_ppi_analysis(
    dna_repair_genes,
    output_prefix="dna_repair_network"
)