Yo tengo un PhViD

Luego de instalar y conocer los rudimentos de R y de instalar los componentes necesarios, podemos utilizar el paquete "PhViD" de R para detectar señales.
Es conveniente primero, si no lo han hecho ya, revisar las nociones básicas de tablas de 2 x 2 y el concepto de independencia en las mismas.
Bueno, voy a desarrollar el uso básico de los algoritmos para detección de señales en bases de datos de Eventos Adversos.

Lo primero que necesitamos es tener una tabla con datos para analizar. Las tablas que se utilizan habitualmente son sumamente extensas y esto hace que sea difícil manipular los datos de una forma entendible en primera instancia. Por eso, utilizaré una tabla pequeña con datos ficticios. La tabla puede descargarse del siguiente link: ListitaEventosDroga.

Descargamos y guardamos la tabla en cualquier carpeta, recordando dónde y vamos a importarla a R.

La importación y exportación de datos a y desde R es una de las tareas menos intuitivas en este software, pero con algo de práctica se le toma la mano.

Para hacerlo, usaremos la instrucción “read.csv”, la cual permite leer un archivo del tipo “csv” que consta de un archivo de texto común y que tiene los datos de las columnas separados por algún carácter (generalmente una coma, por eso “csv” quiere decir “comma separated variable”) y una línea para cada fila de la tabla.

El método read.csv recibe algunos parámetros. El más importante es el que indica la ruta del archivo a importar. Aquí, lo más importante es entrecomillar la ruta del archivo y que, a diferencia de la forma de indicar la ruta en Windows, se usa el carácter “/” en lugar de “\” para separar los nombres de las carpetas.

Los otros parámetros que usaremos serán “header = TRUE” que indica que nuestro archivo tiene los títulos de las columnas, y “sep = “,””, indicando que el separador de columnas es la coma. Finalmente, hay que saber que convertiremos la tabla en un objeto que sea manipulable por el nombre del mismo, asignando la instrucción para importar la tabla a una variable.

La estructura del comando debería seguir este patrón:

variableObjeto = read.csv(“rutaArchivo”, header=TRUE, sep = “,”)

donde “variableObjeto” será reemplazado por el nombre del objeto a crear (el mismo es totalmente arbitrario y pueden elegir el que quieran, siempre y cuando lo recuerden para luego poderlo invocar) y “rutaArchivo” por la ruta para llegar al archivo alojado en el disco. En mi caso la instrucción será la siguiente:

x1=read.csv("C:/Users/Francisco/Documents/Repositorios/00-PhVid/listitaeventodroga.csv", header = TRUE, sep = ";")

Ahora podemos visualizar la tabla con solo escribir el nombre del objeto, como en el siguiente ejemplo:

> x1

        Droga    Evento  n

1 lamotrigina   prurito 10

2   omeprazol   prurito  1

3 lamotrigina       tos  1

4   omeprazol   diarrea  1

5    aspirina sialorrea  1

6  ibuprofeno depresion  1

Vemos que R muestra una tabla con 4 columnas. La primera de ellas es solamente una etiqueta con el número de línea, colocado por R, pero no tiene ningún valor durante el análisis a realizar. De tal modo, la primer columna útil es la que indica el nombre de la droga, la siguiente la que indica el evento ocurrido, y la última indicando el número de reportes para esa combinación Droga-Evento.

Todas las tablas que tengamos que analizar con el paquete “PhViD” deberán tener esta estructura, respetando esas 3 columnas. El número de asociaciones Droga-Evento no tiene limitaciones (a excepción de la memoria del sistema que estemos usando).

Como se puede ver fácilmente, la relación entre los pares Droga-Evento se halla intencionalmente exagerada, resaltando la ocurrencia de uno de ellos (lamotrigina-prurito). Es conveniente repetir que los datos son totalmente inventados a los fines de este tutorial.

Bueno, como expliqué previamente, al iniciar R se cargan solamente los paquetes básicos. Para dotar a R de las capacidades que necesitamos hay que cargar este paquete. Escribimos la siguiente instrucción:

> library(PhViD)

Así el paquete se carga en el espacio de trabajo actual y podemos utilizar sus métodos para manipular los objetos creados.

Si les llega a dar algún error avisando que no puede instalar "LBE", tipean lo siguiente:

source("http://bioconductor.org/biocLite.R")
biocLite("LBE")

Y luego vuelven a tipear "library(PhViD)". Con eso debería solucionarse.

Si bien antes dije que la tabla debe tener una estructura determinada, PhViD necesita cambiar el formato que le entregamos para realizar el análisis que queremos. Para esto, tiene el método “as.PhViD”. Como siempre, asignamos el resultado del método a un nuevo objeto, con la siguiente instrucción:

> l1=as.PhViD(x1)

Ahora el objeto “l1” contiene mis datos en un formato que es analizable mediante los métodos que aporta PhViD.

Voy a mostrar los pasos usando la metodología usada por el Uppsala Monitoring Centre, con el algoritmo llamado “Bayes Confidence Propagation Neural Network” (BCPNN). El método a usar es justamente el “BCPNN”. Los mismos pasos pueden usarse para aplicar cualquiera de los otros algoritmos disponibles en el paquete.

Crearemos un objeto usando el método BCPNN. Podemos nombrar al objeto de cualquier manera, yo lo llamaré “bl1.1”. Para ello simplemente pasamos el objeto “l1” previamente creado, por el método BCPNN asignándolo a la variable nueva.

> bl1.1=BCPNN(l1)

El nuevo objeto creado tiene varios miembros capaces de entregarnos los resultados obtenidos. Si queremos ver los miembros de bl1.1 usamos el método “names” de la siguiente manera.

> names(bl1.1)

[1] "INPUT.PARAM" "ALLSIGNALS"  "SIGNALS"     "NB.SIGNALS"

Los miembros que nos interesan son “SIGNALS” y “ALLSIGNALS”, que nos entregan la lista de pares detectados como señales positivas ordenados por ranking y la lista ordenada por ranking de los valores obtenidos para todos los pares, respectivamente.

> bl1.1$SIGNALS

Nos muestra la tabla con las señales positivas:

par

drug code

event effect

count

expected count

post.H0

n11/E

drug margin

event margin

FDR

FNR

Se

Sp

lamotrigina prurito

lamotrigina

prurito

10

8.066666667

0.180741235

1.239669421

11

11

0.180741235

0.735913946

0.222650697

1

Que en nuestro caso consta solamente del par “lamotrigina prurito”. 

A su vez: 
> bl1.1$ALLSIGNALS

par	drug code	event effect	count	expected count	post.H0	n11/E	drug margin	event margin	FDR	FNR	Se	Sp
lamotrigina prurito	lamotrigina	prurito	10	8.066666667	0.180741235	1.239669421	11	11	0.180741235	0.735913946	0.222650697	1
aspirina sialorrea	aspirina	sialorrea	1	0.066666667	0.295550161	15	1	1	0.238145698	0.715077741	0.414099668	0.922108741
ibuprofeno depresion	ibuprofeno	depresion	1	0.066666667	0.295550161	15	1	1	0.257280519	0.718620375	0.605548639	0.794740052
omeprazol diarrea	omeprazol	diarrea	1	0.133333333	0.330254358	7.5	2	1	0.275523979	0.725705643	0.787566019	0.667371363
lamotrigina tos	lamotrigina	tos	1	0.733333333	0.541284471	1.363636364	11	1	0.328676077	0.781665645	0.912231554	0.525046743
omeprazol prurito	omeprazol	prurito	1	1.466666667	0.677049884	0.681818182	2	11	0.386738379	Inf	1	0.29177773

Nos da la tabla con todos los pares ordenados por intensidad de señal. 

Obviamente todos estos números son difíciles de leer, pero analizaremos alguna de las columnas y les mostraré, que en el corazón del método residen mecanismos que ya conocemos. 
En primer lugar, la primera columna contiene las etiquetas con los pares combinados. La segunda y tercera columnas son los nombres de las drogas y eventos, respectivamente. La tercera columna es solamente el número de veces en las que aparece esa combinación. Hasta aquí, salvo por el orden, no es más que lo que se hallaba contenido en la tabla original. 
La columna “expected count” corresponde a lo que conocemos como “valor esperado” en las tablas de 2 x 2 desarrollado previamente. Para verlo más fácil, construiremos la tabla para el par “lamotrigina prurito”. Necesitamos las marginales y el total general. Sabemos que el total de reportes es de 15. El número de reportes para ese par es de 10. Continuamos colocando el número de reportes que contienen “prurito” y el número de reportes que  contienen “lamotrigina”. La tabla quedaría así:

		lamotrigina
		si	no
prurito	si	10		11
	no
		11		15

Podemos ver que estos dos últimos valores son los que están presentes en las columnas “drug margin” y “event margin”. Es interesante que con sólo conocer éstos valores, podemos calcular el resto sin necesidad de volver a los datos y es por esto que BCPNN solamente guarda estos valores. La tabla completa sería:

		lamotrigina
		si	no
prurito	si	10	1	11
	no	1	3	4
		11	4	15

Sabemos que estos son los “valores observados” y que podemos calcular a través de ellos los “valores esperados”. Esta vez lo haremos sin redondear los números. La tabla es:

		lamotrigina
		si	no
prurito	si	8.06666667	2.93333333	11
	no	2.93333333	1.06666667	4
		11	4	15

Chan! Resulta que el valor esperado para la celda correspondiente al par “lamotrigina prurito” es el mismo que el de la columna “expected count”. Exactamente lo mismo obtendremos calculando las tablas correspondientes a cada par. 
Agregaré que la columna “n11/E” no es nada más que el valor obtenido dividiendo “count” sobre “expected.count”. El resto de las columnas serán analizadas en entregas posteriores. 

El paquete PhViD trae un set de datos propio para testear. Esta tabla se llama “PhViDdata.frame” y contiene 102.483 pares Droga-Evento, totalizando 1.196.487 reportes. ¿Suficiente como para practicar, no? Lo único que hay que hacer es 1- Cargar el dataset con la instrucción data(PhViDdata.frame) y 2- repetir las instrucciones previas, reemplazando el nombre de variable “x1” por “PhViDdata.frame” y listo.

La mayoría de los resultados en datasets reales dará como resultado un enorme número de señales, imposible de manejar en forma práctica. El truco está en que es una lista ordenada por importancia de mayor a menor. El enfoque adecuado es tomar un número manejable de las primeras señales de la lista.

Con esto ya podemos usar el paquete PhViD para generación de señales. Claro que otro tema es que hacer luego de tenerlas, ¿no? 
Bueno, los que se animen, a probar. Espero sus comentarios.