funciones_network_medicine.py

#! /usr/bin/env python

"""
# ---------------------------------------------------------------------------
#
# funciones_network_medicine.py
#
# Archivo con todas las funciones que he utilizado para calcular
# los módulos de las enfermedades, su tamaño, y su significancia 
# estadísticas. También incluye la función que nos permite representar
# tanto la distribución de degrees del módulo de la enfermedad 
# como la distribución de LCCs de redes generadas aleatoriamente (ndp y dp).
#
# TFG: Ciencia de redes y reposicionamiento de fármacos: potencial a través 
# de la medicina de redes
#
# María Marín Tercero
# ----------------------------------------------------------------------------
"""

import pandas as pd
import numpy as np
import seaborn as sns
import matplotlib.pyplot as plt
import networkx as nx
from tabulate import tabulate
from networkx.algorithms import bipartite
import random
from scipy.stats import norm
from itertools import combinations
import re
from itertools import product

# =================================================================================

def arch_distances(PPI):
    
    """
    Esta función calcula los SPL desde cada nodo de la red PPI que demos de entrada al resto de nodos de la red. Después,crea un Dataframe que contiene una
    columna con todos los nodos, y otra columna con la lista de SPL para cada nodo. Por último, crea un archivo llamado distances.csv con ese DataFrame.
    """
    
    nodos = []
    d = []
    # Calculo el Shortest Path Length para cada nodo
    for source_node in PPI.nodes(): #selecciono un nodo de la red PPI
        nodos.append(source_node) #lo añado a la lista de nodos
        path_lengths = [] #creo una lista con todos los path lengths para ese nodo
        for target_node in PPI.nodes(): #selecciono otro nodo de la red PPI
            if source_node != target_node: #me aseguro de que sean nodos diferentes
                if nx.has_path(PPI, source_node, target_node):
                    pl = nx.shortest_path_length(PPI, source=source_node, target=target_node) #calculo el SPL entre ello
                    path_lengths.append(pl) #añado esa distancia a la lista de todos los path lengths de ese nodo
                else:
                    continue
                
        d.append(path_lengths) #lista con los SPL de cada nodo
    
    # Crear DataFrame
    data = {'Nodo': nodos, 'Lista Shortest Path Lenght': d}
    distances = pd.DataFrame(data)
    #Llevo el dataframe a un archivo
    distances.to_csv('distances.csv')

# =================================================================================

def degrees_lista(G):
    """
    Esta función nos devuelve una lista con los nodos y otra lista con sus degrees
    """
    nodos = list(G.nodes())
    degrees = list(dict(G.degree()).values())
    
    return nodos, degrees
    
# =================================================================================

def average_spl(arch):
    """
    Esta función calcula el Average SPL para cada nodo a partir de un archivo que contenga todos los SPL entre todos los nodos de la red PPI
    """
    ave_spl = []
    average_spl_por_nodo = arch.mean(axis=1)
    # Crear un nuevo DataFrame con los resultados
    df_average_spl = pd.DataFrame({'Nodo': average_spl_por_nodo.index, 'AverageSPL': average_spl_por_nodo.values})
    return df_average_spl

# =================================================================================

def genes_enf(enf, archivo):
    """
    Esta función crea una lista con los genes asociados a la enfermedad "enf" en el archivo dis_gen
    
    """
    genes=[]
    for i, dis in enumerate(archivo["dis"]): #obtengo "posición de la enfermedad en la hoja de datos dis_gen (i) y el id de la enfermedad (dis2)"
        if dis == enf:
            gen = archivo["gen"][i] #obtengo el valor de la columna "gen" con la misma posición que cada enfermdad 
            genes.append(gen) #añado los genes a la lista de genes
    return genes

# =================================================================================

def pro_gen_dict(lista_genes, archivo): 
    """
    Esta función crea un diccionario a partir de la lista de genes asociados a la enfermedad con:
    key: proteína asociada a cada gen del archivo gen_pro
    value: gen relacionado a la proteína key en el archivo gen_pro
  """
    dict1={}
    for i, gen in enumerate(archivo["gen"]): #bucle que opera sobre gen_pro, que relaciona genes y proteínas. Me guardo la posición del gen (i) y el id del gen (gen)
        if gen in lista_genes: #busco cada gen de gen_pro en la lista de genes correspondiente de cada enfermedad 
            prot = archivo["pro"][i] #si ese gen está en la lista de genes de cada enfermedad, busco la proteína asociada a la misma posición
            dict1[prot]= gen #añado al diccionario de cada enfermedad la proteína como key y el gen relacionado como value
    return dict1

# =================================================================================

def gen_pro_PPI(dict1, archivo): 
    """
    A partir de un diccionario con las relaciones entre las proteínas y genes se asocian con cada una de nuestras enfermedades,
    esta función se queda con la relación prot:gen del diccionario solo si dicha prot aparece en la red de PPI del archivo pro_pro
    key: proteínas que aparecen en la red PPI
    value: genes relacionados con la proteína key
    """
    dict2={}
    for prot in dict1.keys(): #itero sobre todas las proteínas del diccionario general prot:gen
        if prot in archivo["prA"].tolist() or prot in archivo["prB"].tolist(): #selecciono las prots que aparecen en la PPI
            dict2[prot] = dict1[prot] #añado al diccionario prot:gen de PPI solo las relaciones prot:gen para prots que estén en la PPI
    return dict2   

# =================================================================================

def SG(dic, PPI): 
    """ 
    Datos de entrada: diccionario con las proteínas de la PPI (keys) y los genes (values) asociados para una enfermedad, red     PPI
    Esta función crea una subred solo con las proteínas de la red PPI asociadas a mi enfermedad como nodos
    """

    #Creo una subred solo con las proteínas de la red PPI asociadas a mi enfermedad como nodos 
    SG = nx.subgraph(PPI, dic.keys())
    return SG

# =================================================================================

def lcc(SG): 
    """
    Esta función nos da el LCC de las proteínas de la red PPI asociadas a una enfermedad a partir de una subred
    formada solo con las proteínas asociadas a la enfermedad
    """
    lcc = max(nx.connected_components(SG), key=len) #calculo el lcc (módulo que comprende el mayor número de poteínas de una enfermedad asociadas entre sí)
    
    #Nuestro objetivo es obtener el número de genes que formen parte del LCC de la enfermedad: 
    #El número de proteínas de la enfermedad en el LCC es el mismo número que los genes en el LCC
    #(porque hemos sacado la lista de proteínas del diccionario en el que forman una tupla con sus genes asociados)
    
    
    return lcc

# =================================================================================

def nodes_by_degree(G):
    """
    Esta función nos devuelve un diccionario en el que obtendremos los los degrees como key y, en los values, todos los nodos de la red que contengan ese degree
    """
    degree_dict = {}
    for node in G.nodes():
        degree = G.degree(node)

        if degree not in degree_dict:
            degree_dict[degree] = []

        degree_dict[degree].append(node)

    return degree_dict 

# =================================================================================

def lcc_simulation(G, lcc, PPI, dp=False):
    """
    Datos de entrada: subred módulo de la enfermedad (G), PPI y LCC de la enfermedad (lcc). Dp puede ser False para
    calcular redes preservar la distribución de degrees del módulo de la enfermedad o True en el caso contrario.
    
    Esta función devuelve la media y la desviación estándar del LCC de 1000 redes aleatorias con el mismo número 
    de nodos y links que el gráfico G que demos de entrada. 
    
    Para la creación de estas redes en el caso de dp = False, se distribuyen aleatoriamente las asociaciones entre 
    proteínas de la enfermedad de forma en la red. Por lo tanto, las redes no tendrán la misma estructuctura que el
    módulo de la enfermedad. 
    
    Para la creación de redes en el caso dp = True, se seleccionan nodo en la red PPI con el mismo degree que los
    nodos del módulo de la enfermedad. De esta forma, se crran redes con la misma estructura que e módulo de la 
    enfermedad.
    
    """
    
    #Cálculos previos:
    
    #Obtengo el total de proteínas de la PPI
    nodos_ppi=PPI.nodes()
    
    #Obtengo el total de las proteínas en el módulo de la enfermedad
    nodos_enf = G.nodes()
    
    #Obtengo el número de proteínas del módulos de la enfermedad
    numero_nodos_enf=len(nodos_enf & set(nodos_ppi))
    
    #agrupo los nodos según su degree
    
    grupo_nodos_degree = nodes_by_degree(PPI)


    #1000 simulaciones aleatorias para calcular 1000 LCC
    lista_random= []
    for i in range(1000):
        if dp == False: #simulación para non degree preserving
            #Obtengo un set aleatorio de proteínas dentro del total de nodos de la red PPI.
            #Dicho set presenta un número de nodos que equivale al número de nodos del módulo de la enfermedad
            nodos_random = set(random.sample(list(nodos_ppi), numero_nodos_enf)) #mismo número de nodos que el módulo de la enfermedad, tomados de la lista de prots totales de la PPI aleatoriamente
        if dp == True: #simulación para degree preserving
            nodos_random = set()
            for nodo in nodos_enf: #para cada nodo del módulo de la enfermedad
                degree = PPI.degree(nodo) #obtengo su degree en la PPI
                nodos_disponibles = grupo_nodos_degree[degree] #obtengo un grupo de nodos de la PPI con igual degree al nodo iterado
                control = True
                while(control): #bucle para elegir los nodos sampleados solo una vez
                    nodo_elegido = random.choice(nodos_disponibles) #eligo un nodo entre los nodos de la PPI seleccionados en el paso anterior
                    control = nodo_elegido in nodos_random #compruebo si ese nodo está entre los nodos sampleados
                nodos_random.add(nodo_elegido) #añado ese nodo a la lista de nodos que usaré para crear la red aleatoria
                #solo se añaden nodos que no estén anteriormente entre los nodos sampleados gracias al bucle while  
        r = nx.subgraph(PPI,nodos_random) #subred de la PPI con los nodos aleatorios selccionados
        r = nx.Graph(r) #para remover los links paralelos
        r.remove_edges_from(nx.selfloop_edges(r)) #para eliminar los self-loops que conectan un nodo consigo mismo
        lista_random.append(len(max(nx.connected_components(r), key=len)))
    media = np.mean(lista_random) #media
    std = np.std(lista_random) #desviación estándar
    zscore = (len(lcc) - media)/std
        
    return media, std, zscore, lista_random

# =================================================================================

def rep(nombre_enf, G, PPI, lista_ndp,lista_dp, LCC, media_ndp, std_ndp, zscore_ndp, media_dp, std_dp, zscore_dp):
    """
    Datos de entrada: 
    1. Nombre de la enfermedad
    2. SG: LCC de las proteínas de la red PPI asociadas a una enfermedad
    3. G: Red PPI 
    4. lista_ndp: lista de LCCs (ndp)
    5. lista_dp: lista de LCCs (dp)
    6.LCC: LCC observado de la enfermedad
    7. media_ndp: media de la lista de LCCs (ndp) 
    8. std_ndp: desviación estándar de la lista de LCCs (ndp) 
    9.zscore_ndp: z-score de la lista de LCCs (ndp)
    10. media_dp: media de la lista de LCCs (dp) 
    11. std_dp: desviación estándar de la lista de LCCs (dp)
    12. zscore_dp: z-score de la lista de LCCs (dp)
    
    Esta función realiza una representación de: 
    1. Distribución de degrees
    2. Distribución de LCCs (ndp)
    3. Distribución de LCCs (dp)
    """
     #obtengo una lista con los nodos y los degrees de la PPI
    G_lista_ppi = degrees_lista(PPI)
    
    degree_ppi = pd.DataFrame(list(zip(G_lista_ppi[0], G_lista_ppi[1])), columns=['node','degree'])
     
    #obtengo una lista con los nodos y los degrees de la PPI
    G_lista_enf = degrees_lista(G)
    
    degree_enf = pd.DataFrame(list(zip(G_lista_enf[0], G_lista_enf[1])), columns=['node','degree'])

      
    #Agrupo la lista de nodos y degrees de la PPI en función del degree y cuento cuántos nodos de la red tienen ese degree    
    G_plot_ppi = degree_ppi.groupby('degree').count()
    
    #Agrupo la lista de nodos y degrees de la enfermedad en función del degree y cuento cuántos nodos de la red tienen ese degree
    G_plot_enf = degree_enf.groupby('degree').count()
    
    #Creo una figura en la que añadiré 3 subplots
    fig, axs = plt.subplots(1, 3, figsize=(15, 4))

    #Representación Distribución de Degrees de la PPI en color azul y degrees de los nodos de la enfermedad en color rojo
    axs[0].scatter(G_plot_ppi.index, G_plot_ppi['node'], label='Nodos interactoma', marker='o', color='#79C4FF',s=10)
    axs[0].scatter(G_plot_enf.index, G_plot_enf['node'], label='Nodos '+str(nombre_enf) + " module", marker='o', color='#FF7A7A', s=10)
    axs[0].set_yscale('log')
    axs[0].set_xscale('log')
    axs[0].set_xlabel('Grado', fontsize=12)
    axs[0].set_ylabel('Nº de nodos', fontsize=12)
    axs[0].set_title('Distribución de grados en el interactoma', fontsize=14)
    axs[0].legend(loc='upper right')
    
    #Representación de la distribución de LCC ndp
    axs[1].hist(lista_ndp, color = '#79C4FF', bins=20, edgecolor='black')
    axs[1].set_xlabel('Nº de nodos en LCC', fontsize=12)
    axs[1].set_ylabel('Nº de redes', fontsize=12)
    axs[1].set_title('Distribución del LCC en 1000\n redes aleatoria (ndp)', fontsize=14)
    axs[1].axvline(x=len(LCC), color='#FF7A7A', linestyle='--')
    axs[1].legend(["LCC "+str(nombre_enf)], loc='upper right')

    
    #Representación de la distribución de LCC dp
    axs[2].hist(lista_dp, color = '#79C4FF', bins=20, edgecolor='black')
    axs[2].set_xlabel('Nº de nodos en LCC', fontsize=12)
    axs[2].set_ylabel('Nº de redes', fontsize=12)
    axs[2].set_title('Distribución del LCC en 1000\n redes aleatoria (dp)', fontsize=14)
    axs[2].axvline(x=len(LCC), color='#FF7A7A', linestyle='--')
    axs[2].legend(["LCC "+str(nombre_enf)], loc='upper right')

    
    # Ajusto el diseño y espacio entre gráficos
    plt.tight_layout()
    
    # Leyenda 1
    texto_leyenda_1 = r'$LCC (ndp)_{\mathrm{obs}}$' + ' (' + r'$\mathrm{mean\ LCC} - \mathrm{STD}$, z-score' + ')'
    leyenda_valor_1 = f'{len(LCC)} ({round(media_ndp,2)} - {round(std_ndp,2)}, {round(zscore_ndp,2)})'
    fig.text(0.51, -0.05, texto_leyenda_1+str(" = ")+leyenda_valor_1, fontsize=12, ha='center', va='center', bbox=dict(facecolor='white', alpha=0.5, boxstyle='round'))
    
    #Leyenda 2 
    texto_leyenda_2 = r'$LCC (dp)_{\mathrm{obs}}$' + ' (' + r'$\mathrm{mean\ LCC} - \mathrm{STD}$, z-score' + ')'
    leyenda_valor_2 = f'{len(LCC)} ({round(media_dp,2)} - {round(std_dp,2)}, {round(zscore_dp,2)})'
    fig.text(0.51, -0.15, texto_leyenda_2+str(" = ")+leyenda_valor_2, fontsize=12, ha='center', va='center', bbox=dict(facecolor='white', alpha=0.5, boxstyle='round'))
    
    # Título general encima de los tres gráficos
    plt.suptitle(str(nombre_enf), fontsize=20, y=1.1, fontweight='bold', style='italic')

    # Muestro la figura 
    plt.show()

# =================================================================================   


def archivo_modulos(nombre_enf, nodes_lcc, arch_datos, SG):
    """
    Esta función permite obtener dos archivos: un dataframe con los datos sobre el subgrafo 
    correspondiente al LCC de una enfermedad, y otro con el Id, el grado en el interactoma, 
    y el símbolo del gen para cada nodo del subgrafo.
    
    Datos de entrada: 
    1. Nombre_enf: nombre enfermedad
    2. Nodes_lcc: nodos en el componente más conectado del subgrafo de la enfermedad
    3. Archivo_datos: archivo con los datos sobre el grado de cada nodo en el interactoma
    4. SG: subgrafo de la enfermedad
    """
    
    #obtengo un diccionario con los nodos del módulo de la enfermedad y sus degrees en la PPI
    dict_degrees={} #diccionario vacío
    for node1 in nodes_lcc: #itero sobre los nodos en el módulo de la enfermedad
        for i,node2 in enumerate(arch_datos["Nodo"]): #itero sobre el total de nodos de la PPI
            if node1 == node2: #si un nodo de la PPI está en el módulo de la enfermedad
                dict_degrees[node1]=arch_datos["Degree"][i] #añado al diccionario el nodo como key y su degree como value
    
    #Creo una red con las proteínas del LCC de la enfermedad:
    mod = SG.subgraph(nodes_lcc)
    mod = nx.Graph(mod) #para remover los links paralelos
    mod.remove_edges_from(nx.selfloop_edges(mod)) #para eliminar los self-loops que conectan un nodo consigo mismo
    
    #Convierto el grafo a un DataFrame
    mod_df = nx.to_pandas_edgelist(mod)
    
    #Guardo el df en un archivo
    mod_df.to_csv(str(nombre_enf) + " módulo.csv", index = False)
    
              
    #Junto los datos de id, degrees y símbolos (que rellenaré a continuación)
    data = {'Id': list(dict_degrees.keys()), 'Degrees': list(dict_degrees.values()), "Symbols":[None] * len(list(dict_degrees.keys()))}

    #Añado los símbolos
    for j,prot1 in enumerate(data["Id"]):
        for i, prot2 in enumerate(sim["protein_id"]):
            if prot1==prot2:
                data["Symbols"][j]=sim["gene_symbol"][i]
    
    # Guarda el DataFrame actualizado en el mismo archivo o en uno nuevo
    df = pd.DataFrame(data)
    df.to_csv(str(nombre_enf)+' degrees.csv', index=False)

# ================================================================================= 

def drugs_enf(enf, arch): 
    """
    Esta función crea una lista con los fármacos asociados a cada enfermedad según el archivo dis_dru_the
    """
    drugs=[]
    for x, enf2 in enumerate(arch["dis"]): #para la columna "dis" del archivo dis_dru_the, que recoge los id de las enfermedades, me quedo con el id (enf2) y su posición en la columna(x)
        if enf == enf2: #si se encuentra una coincidencia entre el id de una de nuestras cuatro enfermedades y el id de una enfermedad en el archivo dis_dru_the
            drug = arch["dru"][x] # selecciono la droga asociada al id de esa enfermedad (que se encuentre en su misma fila) 
            drugs.append(drug) #añado ese fármaco a la lista correspondiente (la que se encuentra en la misma posición en la lista drugs que el id de la enfermedad en la lista enf)
    return drugs


# ================================================================================= 

def targets(lista_drugs, arch):
    """
    Datos de entrada: lista con los fármacos que tratan una enfermedad (lista_drugs) y archivo con las relaciones fármaco-diana (arch)
    Esta función nos permite obtener un DataFrame con fármacos en la columna "Fármacos" y sus dianas en la columna "Dianas"
    """
    targets_total = [] #lista en la que añadiré los targets de cada fármaco separados por comas
    farmacos = [] #lista en la que añadiré los fármacos que estén en el archivo arch
    for drug1 in lista_drugs: #itero sobre cada fármaco de la lista de fármacos
        targets=[] #lista vacía en la que añadiré los targets de cada fármaco
        for i, drug2 in enumerate(arch["dru"]): #itero sobre los fármacos del archivo me quedo con el fármaco (drugs2) y fila (i) en la columna de fármacos del archivo dru_pro
            if drug1 == drug2: #si un fármaco de la enfermedad está en el archivo de fármacos-diana
                targets.append(arch["pro"][i]) #añado su diana a la lista de targets, que estará en la misma fila (i) en la columna de dianas                    
        
        if len(targets) > 0: #comrpuebo que la lista de targets no esté vacía
            farmacos.append(drug1) #añado el fármaco a la lista de fármacos para así solamente guardar los fármacos que aparecen en el archivo arch
            targets_total.append(targets) #añado la lista de targets de ese fármaco a la lista de targets de todos los fármacos
            
    data = {"Fármacos": farmacos, "Dianas": targets_total} #junto los datos y los clasifico en Fármacos y Dianas
    df = pd.DataFrame(data) #Creo un DataFrame con los resultados
    return df

# =================================================================================

def calcular_dc_farmaco(lista_targets, arch_dist,  prots_enf, PPI):
    """
    Datos de entrada: lista con las dianas de ese fármaco (lista_targets), proteínas del módulo de la enfermedad 
    (prots_enf), archivo con la matriz de distrancias entre todos los nodos de la PPI (arch_dist), red PPI (PPI)
    Función que devuelve el closest measure (dc) de un fármaco y una enfermedad.
    """
    targets_en_ppi = set(lista_targets) & set(PPI.nodes()) #lista con las dianas del fármaco que también están en la red PPI
    targets_en_enf = set(targets_en_ppi) & set(prots_enf) ##lista con las dianas del fármaco que también forman parte del módulo de la enfermedad
    
    distancias_enf_target = arch_dist.loc[list(targets_en_ppi), list(prots_enf)].values  #genero una matriz con los SPLs entre todas las diana del fármaco y todas las proteínas del módulo de la enfermedad según el archivo de distancias
    
    filas_no_vacias = ~np.isnan(distancias_enf_target).all(axis=1) #me quedo con las filas que no estén vacías, es decir, elimino las dianas que no tienen ningún camino a ninguna proteína de la enfermedad
    
    if np.isnan(distancias_enf_target).all(): #si la matriz anterior está vacía
        return np.nan #no hay camino entre ninguna diana del fármaco fármaco y la enfermedad
    
    elif len(targets_en_enf) == len(targets_en_ppi): #si todas las dianas del fármaco forman parte del módulo de la enfermedad
        return 0 #el valor dc será 0
    
    else: 
        return  np.nanmean(np.nanmin(distancias_enf_target[filas_no_vacias], axis=1)) #en caso contrario, calculo la media entre los SPLs mínimos de las dianas que sí tengan algún camino al módulo de la enfermedad (la media de los valores mínimos de cada fila de mi matriz)

# ================================================================================= 

def proximity(df, arch_dist, prots_enf, PPI):
    """
    Datos de entrada: DataFrame (df) con los "Fármacos" y sus "Dianas; archivo con las 
    distancias entre todos los nodos de la red (arch_dist); lista de proteínas relacionadas con una enfermedad 
    (prots_enf); y red PPI (PPI)
    Función que nos permite obtener un DataFrame con los valores de proximidad (dc) para una lista de fármacos
    y una enfermedad.
    """
    dc_farmaco_total = [] #lista en la que guardaré el valor dc de todos los fármacos
    enfermedades = [] #lista en la que añadiré si el medicamento pertenece a una enfermedad en concreto
    for i, drug in enumerate(df["Fármacos"]): #para cada fármaco del df de fármacos y dianas me quedo con fila (i)
        lista_targets = df["Dianas"][i] #obtengo la lista de dianas de ese fármaco, que estará en la misma fila que el fármaco pero en la columna de dianas
        dc_farmaco_total.append(calcular_dc_farmaco(lista_targets, arch_dist,  prots_enf, PPI))

    data = {'Fármacos': list(df["Fármacos"]), 'dc': dc_farmaco_total} #guardo la relación entre la lista de fármacos y la de dcs
    result_tabla = pd.DataFrame(data) #convierto los reultados en un dataframe

    return result_tabla

# =================================================================================

def proximity_random(df, arch_dist, prots_enf, PPI, num_iteraciones):
    """
    Datos de entrada: DataFrame (df) con los "Fármacos" y sus "Dianas", archivo con las 
    distancias entre todos los nodos de la red (arch_dist), lista de proteínas relacionadas con una enfermedad 
    (prots_enf), red PPI (PPI) y el número de iteraciones (num_iteraciones)
    que serán realizadas para calcular los dc (num_iteraciones)
    Función que nos permite obtener un DataFrame con los valores de proximidad (dc) medios calculados a partir de:
    - 1000 módulos aleatorios de proteínas con el mismo número de proteínas y la misma distribución de degrees
    que el módulo de una enfermedad
    - 1000 módulos aleatorios de dianas con el mimso número de proteínas y la misma distribución de degrees que 
    cada fármaco del DataFrame de fármacos y dianas
    """
    grupo_nodos_degree = nodes_by_degree(PPI) #agrupo los nodos según su degree
    
    # Inicializo una matriz con valores nulos que tenga el mismo número de filas que el total de fármacos de Dataframe
    #y el mismo número de columnas que de iteraciones
    matriz_proximidad = np.full((len(df["Fármacos"]), num_iteraciones), None, dtype=object)

    for i in range(num_iteraciones): #para cada iteración
        
        # Módulo enfermedad aleatorio
        prots_random = set() #creo un set de proteínas aleatorias
        for prot in prots_enf: #para cada proteína del módulo de la enfermedad
            degree_prot = PPI.degree(prot) #calculo su degree
            prots_disponibles = grupo_nodos_degree[degree_prot] #elijo las proteínas de la PPI con su mismo degree
            prots_random.add(np.random.choice(prots_disponibles))  #mismo número de nodos que el módulo de la enfermedad, tomados de la lista de prots totales de la PPI aleatoriamente

        # Módulo diana
        targets_random = set() #creo un set de proteínas dianas aleatorias
        for j, drug in enumerate(df["Fármacos"]): #para cada fármaco me quedo con su fila (j)
            lista_targets = df["Dianas"][j] #me quedo con la lista de dianas de ese fármaco
            lista_targets_PPI = set(lista_targets) & set(PPI.nodes()) #me quedo con las dianas del fármaco que estén en la PPI
            for target in lista_targets_PPI: #para cada diana de la lista de dianas
                degree_target = PPI.degree(target) #calculo su degree
                targets_disponibles = grupo_nodos_degree[degree_target] #elijo las proteínas de la PPI con su mismo degree
                targets_random.add(np.random.choice(targets_disponibles)) #mismo número de nodos que el módulo de la enfermedad, tomados de la lista de prots totales de la PPI aleatoriament
             
            #calculo el dc de cada fármaco con el módulo aleatrio de enfermedad y de dianas
            dc_farmaco_total = calcular_dc_farmaco(targets_random, arch_dist, prots_random, PPI)
            matriz_proximidad[j, i] = dc_farmaco_total #añado el dc de cada fármaco (en la fila j) en la columa de la matriz correspondiente a la iteración (i)
        print(i)
    
    medias_proximidad_farmaco = [] #lista en la que añadiré el promedio de dc tras 1000 iteraciones para cada fármaco
    desviacion = [] #lista en la que añadiré la desviación del dc tras 1000 iteraciones para cada fármaco
    
    for fila in matriz_proximidad: #para cada fila (para cada fármaco)
        if all(x is np.nan for x in fila): #si la fila completa es none
            medias_proximidad_farmaco.append(np.nan) #añado none a la lista medias_proximidad_farmaco
        else:
            medias_proximidad_farmaco.append(np.nanmean(fila)) #añado la media de esa fila
            desviacion.append(np.nanstd(fila)) #desviación estándar
                
    data = {'Fármacos': list(df["Fármacos"]), 'dc_mean': medias_proximidad_farmaco, 'dc_std' : desviacion}
    result_tabla = pd.DataFrame(data)

    return  result_tabla


# ================================================================================= 

def determinar_tratamiento(row, arch):
    """
    Esta función toma las filas de un archivo tipo DataFrame y añade una columna indicando 
    si dichos fármacos se emplean para el tratamiento de la demencia, la bipolaridad,
    la epilepsia o la esquizofrenia. Finalmente, devuelve el DataFrame original con 
    la nueva columna Tratamientos añadida.
    """
    enfermedad = row['Enfermedades'] #selecciona cada fila de la columna enfermedades
    farmaco = row['Fármacos'] #selecciona cada fila de la columna fármacos
    codigo_dis = { #código que relaciona nombre y id de la enfermedad
        'Demencia': 'C0497327',
        'Bipolaridad': 'C0005586',
        'Epilepsia': 'C0014544',
        'Esquizofrenia': 'C0036341'
    }

    # Añado "yes" a tratamiento si ese fármaco se usa para tratar una de las enfermedades del codigo_dis
    codigo_tratamiento = 'yes' if arch[(arch['dis'] == codigo_dis[enfermedad]) & (arch['dru'] == farmaco)].shape[0] > 0 else 'unknown' #si no, añado "unkwown"

    return codigo_tratamiento

# =================================================================================

def combinar_proximidades(enfermedades,farmacos, proximity_files, proximity_random_files, arch):
    """
    Esta función crea un DataFrame que contiene los resultados para la proximidad, 
    la media, la desviación típica y el zscore de la proximidad (obtenidos a partir
    de los DataFrames aportados como datos de entrada) para todas las combinaciones 
    entre la lista de fármacos y la lista de enfermedades proporcionada como datos de entrada.
    Columnas DataFrame resultado:
    1. Enfermedades (Demencia, Bipolaridad, Epilepsia o Esquizofrenia)
    2. Fármacos
    3. Tratamiento: columna que indica si el el fármaco de la columna 2 trata a la enfermedad
    de la columna 1 ("yes" en caso afirmativo, "unkwown" en el caso contrario
    4. Valor de proximidad a la Enfermedad de la columna 1
    5. Proximidad aleatoria media a la Enfermedad de la columna 1
    6. Desviación estándar de en las proximidades aleatorias a la Enfermedad de la columna 1
    7. Z-Score de la proximidad observada (indicada en la columna 3) en comparación con
    la proximidad aleatoria media (en la columna 4) y la desviación estándar de la 
    proximidad aleatoria (columna 5)
    
    Datos de entrada:
    1. Lista fármacos
    2. Lista enfermedades (Demencia, Bipolaridad, Epilepsia o Esquizofrenia)
    3. Diccionario que contiene un Dataframe para cada enfermedad. Dichos DataFrames almacenan la 
    información sobre la proximidad observada de cada fármaco a la enfermedad (proximity_files)
    4. Diccionario que contiene un Dataframe para cada enfermedad. Dichos DataFrames almacenan la 
    información sobre la la proximidad aleatoria media y su desviación estándar de cada fármaco a 
    la enfermedad (proximity_random_files)
    5. Archivo con información sobre las relaciones entre fármacos y enfermedades (arch)
    """
    
    combinaciones = list(product(enfermedades, farmacos)) #todas las combinaciones enfermedad-fármaco posibles
    df_enfermedades_farmacos = pd.DataFrame(combinaciones, columns=['Enfermedades', 'Fármacos']) #DataFrame a partir de las combinaciones 
    df_enfermedades_farmacos['Tratamiento'] = df_enfermedades_farmacos.apply(lambda row: determinar_tratamiento(row, arch), axis=1)
 #Añado la columna Tratamiento
    
    filas_combinadas = [] #fila vacía en la que añadiré cada fila combinada de los dos DataFrames de entrada

    for index, row in df_enfermedades_farmacos.iterrows(): #itero sobre las filas y columnas del DataFrame creado
        enfermedad = row['Enfermedades'] #selecciona cada fila de la columna Enfermedades
        farmaco = row['Fármacos'] #selecciona cada fila de la columna Fármacos
        tratamiento = row['Tratamiento'] #selecciona cada fila de la columna Tratamiento

        #Obtengo el archivo de proximidad correspondiente a cada enfermedad 
        if enfermedad in proximity_files:
            df_proximidad = proximity_files[enfermedad]

            # Buscar el farmaco en el archivo de proximidad y obtener el valor dc
            #Busco cada fármaco en el DataFrame con los valores de proximidad observada para cada enfermedad
            filtro = (df_proximidad['Fármacos'] == farmaco)
            if filtro.any(): #si se encuentra el fármaco en el DataFrame
                dc = df_proximidad[filtro]['dc'].values[0] #guardo su valor de proximidad a la enfermedad (dc)

                #Obtengo el archivo de proximidad correspondiente a cada enfermedad
                df_proximidad_random = proximity_random_files[enfermedad]

                #Busco cada fármaco en el DataFrame con los valores de proximidad aleatoria para cada enfermedad
                filtro_random = (df_proximidad_random['Fármacos'] == farmaco)
                if filtro_random.any(): #si se encuentra el fármaco en el DataFrame
                    dc_mean = df_proximidad_random[filtro_random]['dc_mean'].values[0] #guardo su valor de proximidad aleatoria media (dc_mean)
                    dc_std = df_proximidad_random[filtro_random]['dc_std'].values[0] #guardo su valor de desviación de la proximidad aleatoria (dc_std)

        # Añado la fila a la lista
        filas_combinadas.append({'Enfermedades': enfermedad, 'Fármacos': farmaco, 'Tratamiento': tratamiento, 'dc': dc, 'dc_mean': dc_mean, 'dc_std': dc_std})

    # Crear el DataFrame con todos los datos combinados
    df_combinado = pd.DataFrame(filas_combinadas)
    df_combinado['dc_zscore'] = (df_combinado['dc'] - df_combinado['dc_mean']) / df_combinado['dc_std'] #añado una columna con el valor del zscore

    return df_combinado

# =================================================================================

def rep_prox(df_combinado, nombre_enf):
    """
    Esta función representa en un boxplot la distribución de la cercanía y a una enfermedad y
    la distribución de su z-score para el grupo de fármacos que se emplean para tratar dicha 
    enfermedad y el grupo de fármacos que no se emplean para su tratamiento (unkwown)
    
    Datos de entrada: 
    1. df_combinado: Dataframe combinado con los fármacos, su cercanía observada, su cercanía aleatoria media,
    la devsiación estándar de la proximidad aleatoria, su z-score, una columna que indica la enfermedad
    (la demencia, epilepsia, la bipolaridad o la esquizofrenia) y otra columna que indica si el fármaco
    se usa para el tratamiento de la enfermedad en su misma fila.
    2. nombre_enf: Nombre de la enfermedad
    """
    # Filtro los datos para generar dos nuevos dataframes, uno con los datos de los fármacos que se emplean como tratamiento de la enfermedad, y otro con los datos para el resto de fármacos
    farmacos_con_enf = df_combinado[(df_combinado['Enfermedades'] == nombre_enf) & (df_combinado['Tratamiento'] == 'yes')]
    farmacos_sin_enf = df_combinado[(df_combinado['Enfermedades'] == nombre_enf) & (df_combinado['Tratamiento'] == 'unknown')]
    
    combined_data = pd.concat([farmacos_con_enf.assign(Tratamiento='Tratamiento'), farmacos_sin_enf.assign(Tratamiento='Unknown')])
    
    # Combino los dos conjuntos de datos en un mismo subplot
    fig, axes = plt.subplots(1, 2, figsize=(14, 6))  # Creo un subplot con 1 fila y 2 columnas
    
    # Represento el boxplot con ambas distribuciones de la proximidad
    sns.boxplot(x='Tratamiento', y='dc', data=combined_data, hue='Tratamiento', ax=axes[0], palette={'Tratamiento': '#FF7A7A', 'Unknown': '#79C4FF'}, dodge=False, medianprops=dict(linewidth=2), legend=False)
    axes[0].set_ylabel('Cercanía ($\mathregular{d_c}$)', fontsize=12)
    axes[0].set_xlabel('')
    for label in axes[0].get_xticklabels():
        label.set_fontsize(12)

    # Represento el boxplot con ambas distribuciones de la proximidad
    sns.boxplot(x='Tratamiento', y='dc_zscore', data=combined_data, hue='Tratamiento', ax=axes[1], palette={'Tratamiento': '#FF7A7A', 'Unknown': '#79C4FF'}, dodge=False, medianprops=dict(linewidth=2), legend=False)
    axes[1].set_ylabel('Proximidad [z-score ($\mathregular{d_c}$)]', fontsize=12)
    axes[1].set_xlabel('')
    for label in axes[1].get_xticklabels():
        label.set_fontsize(12)
    
    plt.suptitle(str(nombre_enf), fontsize=20, y=0.98, fontweight='bold', style='italic')


    plt.tight_layout()  # Ajusto la disposición del subplot para evitar superposiciones
    plt.show()
    
# =================================================================================

def targets_modulo(df_combinado, targets, lcc_files):
    """
    Esta función permite obtener un DataFrame con las siguientes columnas:
    1. Enfermedades (Demencia, Bipolaridad, Epilepsia o Esquizofrenia)
    2. Fármacos
    3. Tratamiento: columna que indica si el el fármaco de la columna 2 trata a la enfermedad
    de la columna 1 ("yes" en caso afirmativo, "unkwown" en el caso contrario
    4. Valor de proximidad a la Enfermedad de la columna 1
    5. Proximidad aleatoria media a la Enfermedad de la columna 1
    6. Desviación estándar de en las proximidades aleatorias a la Enfermedad de la columna 1
    7. Z-Score de la proximidad observada (indicada en la columna 3) en comparación con
    la proximidad aleatoria media (en la columna 4) y la desviación estándar de la 
    proximidad aleatoria (columna 5)
    8. Target en módulo: indica si el fármaco en de la columna 2 presenta alguna de sus
    proteínas diana en el módulo de la enfermedad en la columna 1.
    
    Datos de entrada: 
    1. df_combinado: DataFrame que contiene las 7 primeras columnas del archivo que se 
    obtiene como resultado (enfermedades, fármacos, tratamiento, dc, dc_mean, dc_std, dc_zscore)
    2. targets: DataFrame con la lista de proteínas diana de cada fármaco
    3. lcc_files: diccionario con el nombre de la enfermedad (key) y la lista de proteínas que
    forman parte de su módulo (values).
    """
    
    targets_modulo = []  # lista para almacenar si el fármaco de la columna 2 en cada fila tiene una proteína diana en el módulo de la enfermedad de la columna 1 situada en su misma fila.
    
    # Itero sobre cada fila del DataFrame df_combinado
    for index, row in df_combinado.iterrows():
        farmaco = row['Fármacos'] #selecciono cada fila de la columna Fármacos
        enfermedad = row['Enfermedades'] #selecciono cada fila de la columna Enfermedades

        # Obtengo los targets de cada fármaco según el DataFrame de targets
        farmaco_targets = set(targets[targets['Fármacos'] == farmaco]['Dianas'].iloc[0])

        # Obtengo el módulo (lcc) de la enfermedad desde el diccionario de lcc_files
        lcc_enf = lcc_files[enfermedad]

        # Verifico si alguno de los targets está en el módulo (lcc) de la enfermedad
        resultado_farmaco = 'Sí' if farmaco_targets.intersection(lcc_enf) else 'No'

        # Agrego el resultado a la lista de resultados
        targets_modulo.append(resultado_farmaco)

    # Añado la lista de resultados como una nueva columna del DataFrame df_combinado
    df_combinado['Target en módulo'] = targets_modulo

    return df_combinado

# =================================================================================