Lab4 - Bases de Datos en Bioinformática
| Profesores | Autores |
|---|---|
| Mario Inostroza Ponta | Claudia Rivera |
| Manuel Villalobos Cid | Juan Giglio |
| Adolfo Gonzalez |
2023-12-09
Introducción
Este laboratorio se basa en el análisis de la Serie de datos que se centra en el estudio de microarrays utilizando muestras de piel de pacientes con vitíligo, abordando tanto las áreas lesionadas como las no lesionadas, nombrada serie GSE75819, estos datos fueron publicados el 1 de noviembre de 2016, en el sitio web “www.ncbi.nlm.nih.gov” que corresponde a la Biblioteca Nacional de Medicina de los Estados Unidos (National Library of Medicine, NLM), que forma parte de los Institutos Nacionales de Salud (National Institutes of Health, NIH). El organismo de interés es Homo sapiens, y el tipo de experimento realizado se clasifica como “Perfil de expresión por matriz”. En resumen, el laboratorio se enfoca en examinar y comparar los perfiles de expresión génica en la piel afectada y no afectada de individuos con vitíligo, utilizando datos provenientes de microarrays para obtener una comprensión más profunda de las diferencias genéticas asociadas con esta condición dermatológica.
Objetivo
Aplicar un algoritmo de análisis de datos de expresión génica sobre un conjunto de datos, como clustering, expresión diferencial (Selección Características) y caracterización biológica de resultados (Enriquecimiento Biológico).
Para el desarrollo de este laboratorio el primer paso es instalar las librerias necesarias.
# Instalar y cargar paquetes necesarios
if (!requireNamespace("factoextra", quietly = TRUE)) {
install.packages("factoextra")
}
if (!requireNamespace("cluster", quietly = TRUE)) {
install.packages("cluster")
}
if (!requireNamespace("Amelia", quietly = TRUE)) {
install.packages("Amelia")
}
if (!requireNamespace("dplyr", quietly = TRUE)) {
install.packages("dplyr")
}
# Instalar y cargar la librería "gplots" para usar "heatmap"
if (!requireNamespace("gplots", quietly = TRUE)) {
install.packages("gplots")
}
# Verificar si la librería DESeq2 está instalada
if (!requireNamespace("DESeq2", quietly = TRUE)) {
# Si no está instalada, instalarla usando BiocManager
if (!requireNamespace("BiocManager", quietly = TRUE)) {
install.packages("BiocManager")
}
BiocManager::install("DESeq2")
}
BiocManager::install("GEOquery")## Bioconductor version 3.18 (BiocManager 1.30.22), R 4.3.0 (2023-04-21 ucrt)## Warning: package(s) not installed when version(s) same as or greater than current; use
## `force = TRUE` to re-install: 'GEOquery'## Installation paths not writeable, unable to update packages
## path: C:/Program Files/R/R-4.3.0/library
## packages:
## class, cluster, foreign, KernSmooth, lattice, MASS, Matrix, mgcv, nlme,
## nnet, spatial, survival## Old packages: 'BiocGenerics', 'Biostrings', 'curl', 'data.table', 'dplyr',
## 'DT', 'e1071', 'fansi', 'GenomeInfoDb', 'httpuv', 'IRanges', 'jsonlite',
## 'later', 'msa', 'Rcpp', 'RCurl', 'rlang', 'rpart', 'S4Vectors', 'sass',
## 'seqinr', 'utf8', 'vctrs', 'vroom', 'XVector', 'zlibbioc'## Warning: package 'gplots' was built under R version 4.3.2##
## Attaching package: 'gplots'## The following object is masked from 'package:stats':
##
## lowess## Warning: package 'dplyr' was built under R version 4.3.2##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union## Warning: package 'factoextra' was built under R version 4.3.2## Loading required package: ggplot2## Warning: package 'ggplot2' was built under R version 4.3.2## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBa## Warning: package 'cluster' was built under R version 4.3.2## Warning: package 'Amelia' was built under R version 4.3.2## Loading required package: Rcpp## Warning: package 'Rcpp' was built under R version 4.3.1## ##
## ## Amelia II: Multiple Imputation
## ## (Version 1.8.1, built: 2022-11-18)
## ## Copyright (C) 2005-2023 James Honaker, Gary King and Matthew Blackwell
## ## Refer to http://gking.harvard.edu/amelia/ for more information
## #### Warning: package 'NbClust' was built under R version 4.3.1## Warning: package 'FactoMineR' was built under R version 4.3.2## Warning: package 'pheatmap' was built under R version 4.3.2## Warning: package 'DESeq2' was built under R version 4.3.1## Loading required package: S4Vectors## Warning: package 'S4Vectors' was built under R version 4.3.1## Loading required package: stats4## Loading required package: BiocGenerics##
## Attaching package: 'BiocGenerics'## The following objects are masked from 'package:dplyr':
##
## combine, intersect, setdiff, union## The following objects are masked from 'package:stats':
##
## IQR, mad, sd, var, xtabs## The following objects are masked from 'package:base':
##
## anyDuplicated, aperm, append, as.data.frame, basename, cbind,
## colnames, dirname, do.call, duplicated, eval, evalq, Filter, Find,
## get, grep, grepl, intersect, is.unsorted, lapply, Map, mapply,
## match, mget, order, paste, pmax, pmax.int, pmin, pmin.int,
## Position, rank, rbind, Reduce, rownames, sapply, setdiff, sort,
## table, tapply, union, unique, unsplit, which.max, which.min##
## Attaching package: 'S4Vectors'## The following objects are masked from 'package:dplyr':
##
## first, rename## The following object is masked from 'package:gplots':
##
## space## The following object is masked from 'package:utils':
##
## findMatches## The following objects are masked from 'package:base':
##
## expand.grid, I, unname## Loading required package: IRanges## Warning: package 'IRanges' was built under R version 4.3.1##
## Attaching package: 'IRanges'## The following objects are masked from 'package:dplyr':
##
## collapse, desc, slice## The following object is masked from 'package:grDevices':
##
## windows## Loading required package: GenomicRanges## Warning: package 'GenomicRanges' was built under R version 4.3.1## Loading required package: GenomeInfoDb## Warning: package 'GenomeInfoDb' was built under R version 4.3.1## Loading required package: SummarizedExperiment## Warning: package 'SummarizedExperiment' was built under R version 4.3.1## Loading required package: MatrixGenerics## Warning: package 'MatrixGenerics' was built under R version 4.3.1## Loading required package: matrixStats## Warning: package 'matrixStats' was built under R version 4.3.2##
## Attaching package: 'matrixStats'## The following object is masked from 'package:dplyr':
##
## count##
## Attaching package: 'MatrixGenerics'## The following objects are masked from 'package:matrixStats':
##
## colAlls, colAnyNAs, colAnys, colAvgsPerRowSet, colCollapse,
## colCounts, colCummaxs, colCummins, colCumprods, colCumsums,
## colDiffs, colIQRDiffs, colIQRs, colLogSumExps, colMadDiffs,
## colMads, colMaxs, colMeans2, colMedians, colMins, colOrderStats,
## colProds, colQuantiles, colRanges, colRanks, colSdDiffs, colSds,
## colSums2, colTabulates, colVarDiffs, colVars, colWeightedMads,
## colWeightedMeans, colWeightedMedians, colWeightedSds,
## colWeightedVars, rowAlls, rowAnyNAs, rowAnys, rowAvgsPerColSet,
## rowCollapse, rowCounts, rowCummaxs, rowCummins, rowCumprods,
## rowCumsums, rowDiffs, rowIQRDiffs, rowIQRs, rowLogSumExps,
## rowMadDiffs, rowMads, rowMaxs, rowMeans2, rowMedians, rowMins,
## rowOrderStats, rowProds, rowQuantiles, rowRanges, rowRanks,
## rowSdDiffs, rowSds, rowSums2, rowTabulates, rowVarDiffs, rowVars,
## rowWeightedMads, rowWeightedMeans, rowWeightedMedians,
## rowWeightedSds, rowWeightedVars## Loading required package: Biobase## Warning: package 'Biobase' was built under R version 4.3.1## Welcome to Bioconductor
##
## Vignettes contain introductory material; view with
## 'browseVignettes()'. To cite Bioconductor, see
## 'citation("Biobase")', and for packages 'citation("pkgname")'.##
## Attaching package: 'Biobase'## The following object is masked from 'package:MatrixGenerics':
##
## rowMedians## The following objects are masked from 'package:matrixStats':
##
## anyMissing, rowMedians## Warning: package 'GEOquery' was built under R version 4.3.1## Setting options('download.file.method.GEOquery'='auto')## Setting options('GEOquery.inmemory.gpl'=FALSE)Carga de datos
Una vez completada la instalación de las librerias, se realiza la carga de datos seleccionada y se lleva a cabo un paso crucial para la preparación de los datos: la eliminación de filas que contienen valores faltantes. Esta tarea es esencial para garantizar la integridad y la consistencia de los conjuntos de datos utilizados en el análisis.
metadata <- read.csv("C:/Users/clace/Documents/Bioinformatica/Tarea4/datosGeneticos.csv", stringsAsFactors = TRUE)
head(metadata )## PROBE_ID TargetID ProbeID VD13N.AVG_Signal VD13V.AVG_Signal
## 1 ILMN_1343291 EEF1A1 3450719 24091.13000 27246.04
## 2 ILMN_1343295 GAPDH 4490161 8316.76500 4792.797
## 3 ILMN_1651199 LOC643334 1770632 26.06999 28.558759
## 4 ILMN_1651209 SLC35E2 2940221 40.07432 34.587432
## 5 ILMN_1651210 DUSP22 380154 29.34726 29.91686829
## 6 ILMN_1651221 LOC642820 2490538 31.81759 36.876435
## VD18N.AVG_Signal VD18V.AVG_Signal VD19N.AVG_Signal VD19V.AVG_Signal
## 1 34437.45000 27645.24000 31980.22000 27412.39000
## 2 9885.07400 5513.59500 10703.85000 4421.30100
## 3 31.92645 28.66349 25.30298 32.70277
## 4 48.41365 50.86734 55.19860 37.60007
## 5 31.07811 24.16805 23.83542 32.11145
## 6 37.03900 32.77100 30.07098 30.08274
## VD20N.AVG_Signal VD20V.AVG_Signal VD23N.AVG_Signal VD23V.AVG_Signal
## 1 32349.11000 33625.73000 29232.30000 29365.15000
## 2 15490.89000 6156.72700 6883.87500 4733.13500
## 3 22.41737 40.01635 30.09873 25.16234
## 4 25.39535 34.68272 28.42441 35.45035
## 5 22.15163 25.55137 24.46381 28.07598
## 6 35.19772 32.41304 31.52364 28.07241
## VD24N.AVG_Signal VD24V.AVG_Signal VD22N.AVG_Signal VD22V.AVG_Signal
## 1 24838.10000 23772.65000 32266.45000 27873.5900
## 2 8330.29600 6551.99600 15001.35000 5697.1550
## 3 30.67451 28.87391 113.17528 101.0734
## 4 33.55914 32.82220 111.26250 102.0010
## 5 24.81333 25.39571 95.01883 97.9926
## 6 34.68340 32.63612 113.98072 109.2442
## VD32N.AVG_Signal VD32V.AVG_Signal VD33N.AVG_Signal VD33V.AVG_Signal
## 1 27320.54000 27523.92000 38981.8600 37865.2500
## 2 9571.12600 10533.13000 9700.5470 7674.3410
## 3 108.93726 104.48459 132.4976 130.4750
## 4 115.11550 110.52784 164.3924 149.6159
## 5 94.79396 92.96444 117.4028 105.7548
## 6 106.76630 106.80429 141.5198 132.9117
## VD37N.AVG_Signal VD37V.AVG_Signal VD39N.AVG_Signal VD39V.AVG_Signal
## 1 34312.4100 30110.08000 22330.37000 22101.17000
## 2 13501.8700 9815.60500 11060.76000 5271.77700
## 3 118.4134 115.08817 108.55663 100.60485
## 4 138.4799 118.96092 131.08160 115.50391
## 5 100.4321 94.29914 93.40465 91.84454
## 6 128.4393 126.19493 134.54600 111.66599
## VD38N.AVG_Signal VD38V.AVG_Signal VD40N.AVG_Signal VD40V.AVG_Signal
## 1 33733.51000 31024.4500 32635.85000 30112.86000
## 2 13999.62000 11655.4900 9219.03800 6223.14500
## 3 117.53503 112.6142 109.77959 104.88056
## 4 134.94400 118.4146 123.17994 116.95463
## 5 94.39503 101.1662 96.83482 92.40464
## 6 119.00259 116.9174 115.09498 114.16994
## VD44N.AVG_Signal VD44V.AVG_Signal VD49N.AVG_Signal VD49V.AVG_Signal
## 1 34140.46000 31710.83000 28029.32000 26762.01000
## 2 13696.56000 5346.32100 8550.78000 6740.67100
## 3 106.43611 100.99194 97.21587 93.06997
## 4 118.02637 113.96832 115.02520 104.28199
## 5 93.48535 90.90065 85.91356 88.01222
## 6 109.81473 108.24347 100.09927 103.81960# Manejar valores faltantes eliminando filas con NA
expr_data <- na.omit(metadata )
expr_data <- data.frame(expr_data)
num_registros <- nrow(expr_data)
cat("El número de registros en expr_data es:", num_registros, "\n")## El número de registros en expr_data es: 48803set.seed(11567)
expr_data$PROBE_ID <- as.numeric(as.factor(expr_data$PROBE_ID))
expr_data$TargetID <- as.numeric(as.factor(expr_data$TargetID))
expr_data$ProbeID <- as.numeric(as.integer(expr_data$ProbeID))
gse <- getGEO("GSE75819")## Found 1 file(s)## GSE75819_series_matrix.txt.gzPara llevar a cabo esta operación, se recurre a la función na.omit en el lenguaje de programación R, que actúa de manera efectiva al eliminar observaciones que presentan valores ausentes, identificados como NA en inglés. Esta práctica es fundamental en el ámbito del análisis de datos, ya que permite trabajar con conjuntos de datos completos y preparados para análisis estadísticos y otras operaciones, asegurando la fiabilidad y la coherencia de los resultados obtenidos.
Algoritmo de Clustering
La aplicación de técnicas de clustering en bioinformática es crucial para analizar conjuntos de datos biológicos complejos. Estas herramientas identifican patrones y similitudes entre entidades biológicas, como secuencias genéticas o perfiles de expresión. El análisis de clusters revela grupos con funciones genéticas o características similares, siendo esencial para comprender la función de genes, relaciones genéticas y clasificación de muestras biológicas. En la bioinformática, estas técnicas se utilizan para descubrir estructuras complejas, como subgrupos de pacientes en estudios genómicos o la clasificación de especies según secuencias genéticas. En resumen, el clustering en bioinformática desempeña un papel crucial en el descubrimiento de patrones biológicos, avanzando en la comprensión de procesos biológicos y el diagnóstico de enfermedades.
Método de K-means
El método de K-means, un algoritmo de clasificación no supervisada en minería de datos, agrupa objetos en k clusters con eficacia. La selección del número de clusters, k, es crucial y puede adaptarse según las necesidades del análisis. Este algoritmo minimiza la suma de distancias entre cada objeto y el centroide de su cluster asignado. En su proceso iterativo, se seleccionan aleatoriamente k centroides, los objetos se asignan al cluster más cercano inicialmente, y luego los centroides se recalculan como promedio de observaciones en cada cluster. Este ciclo se repite hasta que los centroides convergen. Ampliamente usado en diversas disciplinas como marketing y bioinformática, el método de K-means es valioso para analizar grandes conjuntos de datos, revelando patrones sin necesidad de etiquetas de clases previas debido a su simplicidad y eficiencia.
Para determinar el número óptimo de clusters en un proceso de clustering, se utiliza la herramienta fviz_nbclust, especialmente en R. Esta función emplea métodos como el de codo y la silueta para evaluar la variabilidad intra-cluster y la cohesión entre clusters. Identificar el “codo” en la curva de variabilidad ayuda a determinar el número adecuado de clusters, mientras que la silueta ofrece una evaluación detallada de la cohesión y separación de los clusters. Integrar fviz_nbclust en el proceso de selección de clusters mejora la toma de decisiones, asegurando que el número elegido sea óptimo para la interpretación y utilidad de los resultados del clustering.
Resultados
a. Generalidades del cluster
Se genera un subconjunto de datos, tomando las primeras 15000 filas de los datos originales. Dada la cantidad de registros de 48803, no es posible realizar el laboratorio con el total por limitaciones de hardware.
b. Índice silueta general
## Cluster 1 - Tamaño del grupo: 13391
## Cluster 2 - Tamaño del grupo: 124
## Cluster 3 - Tamaño del grupo: 54
## Cluster 4 - Tamaño del grupo: 318
## Cluster 5 - Tamaño del grupo: 1059
## Cluster 6 - Tamaño del grupo: 54## Índice silueta general: 0.82779## Grupo con mejor índice silueta: 14304c. Grupo con mejor índice silueta
Del proceso anterior se obtienen los 3 grupos con mejor índice silueta:
## Los tres mejores grupos son: 14304 5963 7444## Las etiquetas de los tres mejores clusters son: 1 1 1
Selección de genes diferencialmente expresados
## [1] "title" "geo_accession"
## [3] "status" "submission_date"
## [5] "last_update_date" "type"
## [7] "channel_count" "source_name_ch1"
## [9] "organism_ch1" "characteristics_ch1"
## [11] "treatment_protocol_ch1" "molecule_ch1"
## [13] "extract_protocol_ch1" "label_ch1"
## [15] "label_protocol_ch1" "taxid_ch1"
## [17] "hyb_protocol" "scan_protocol"
## [19] "description" "description.1"
## [21] "data_processing" "platform_id"
## [23] "contact_name" "contact_email"
## [25] "contact_laboratory" "contact_institute"
## [27] "contact_address" "contact_city"
## [29] "contact_state" "contact_zip/postal_code"
## [31] "contact_country" "supplementary_file"
## [33] "data_row_count" "tissue:ch1"## converting counts to integer mode## Warning in DESeq(dds): the design is ~ 1 (just an intercept). is this intended?## estimating size factors## estimating dispersions## gene-wise dispersion estimates## mean-dispersion relationship## final dispersion estimates## fitting model and testing## -- replacing outliers and refitting for 77 genes
## -- DESeq argument 'minReplicatesForReplace' = 7
## -- original counts are preserved in counts(dds)## estimating dispersions## fitting model and testing## log2 fold change (MLE): Intercept
## Wald test p-value: Intercept
## DataFrame with 48803 rows and 6 columns
## baseMean log2FoldChange lfcSE stat pvalue
## <numeric> <numeric> <numeric> <numeric> <numeric>
## ILMN_1343291 18103.2637 14.14396 0.147479 95.90503 0.00000e+00
## ILMN_1343295 5217.4749 12.34914 0.158579 77.87361 0.00000e+00
## ILMN_1651199 22.8675 4.51523 0.399152 11.31206 1.14377e-29
## ILMN_1651209 28.0975 4.81237 0.243826 19.73688 1.04008e-86
## ILMN_1651210 19.5927 4.29224 0.441011 9.73274 2.18633e-22
## ... ... ... ... ... ...
## ILMN_2415911 90.3068 6.49676 0.1122460 57.87967 0.00000e+00
## ILMN_2415926 391.2154 8.61182 0.1265364 68.05805 0.00000e+00
## ILMN_2415949 158.8202 7.31125 0.0958825 76.25218 0.00000e+00
## ILMN_2415979 381.1401 8.57418 0.2164718 39.60875 0.00000e+00
## ILMN_2416019 20.3323 4.34570 0.4809462 9.03573 1.62913e-19
## padj
## <numeric>
## ILMN_1343291 0.00000e+00
## ILMN_1343295 0.00000e+00
## ILMN_1651199 1.53228e-29
## ILMN_1651209 2.31228e-86
## ILMN_1651210 2.49013e-22
## ... ...
## ILMN_2415911 0.00000e+00
## ILMN_2415926 0.00000e+00
## ILMN_2415949 0.00000e+00
## ILMN_2415979 0.00000e+00
## ILMN_2416019 1.66252e-19Justificación:
Selección de Grupos: Se seleccionaron los 3 grupos con el mejor índice de silueta para garantizar una alta cohesión y separación entre clusters, lo que indica que estos clusters son más homogéneos internamente y más distintos entre sí.
Análisis de Expresión Diferencial: Se utilizó DESeq2, una herramienta robusta para el análisis de expresión diferencial en datos de secuenciación de ARN. DESeq2 tiene en cuenta la variabilidad intrínseca en los datos de ARN-Seq y ajusta estadísticamente para mejorar la detección de genes diferencialmente expresados.
Filtrado de Genes Diferencialmente Expresados: Se seleccionaron los genes con un ajuste p-valor (padj) menor que 0.05 y un cambio en la expresión mayor a 2 veces (en valor absoluto) para resaltar los genes más significativos y biológicamente relevantes.
Heatmap: Se generó un heatmap para visualizar la expresión de los genes diferencialmente expresados en cada grupo seleccionado, proporcionando una representación gráfica de los patrones de expresión génica.