LAB 6 - Paquets de R per a la bioinformàtica
Carles Castellà Barnaus
2026-05-16
EXERCICIS
EXERCICI 1
A partir del conjunt de dades heart del paquet survival, feu una aplicació Shiny que a partir de la selecció de dues variables del conjunt de dades representi un gràfic de dispersió i un de línies, a més d’un histograma sol per a la variable horitzontal seleccionada. Per fer aquest exercici podeu basar-vos en l’Exemple 2 realitzat en l’apartat 1.1. Ecosistema de R: Entorns reproduïbles i visualització de dades. S’aconsella ampliar coneixements en l’ús de diferents paràmetres en les funcions definides, així com explorar un altre tipus de funcionalitats per aprofundir en l’estudi del conjunt de dades, a criteri de l’estudiant.
library(shiny)
library(survival)
data("heart")
# Definim la interficie de Shiny
ui <- fluidPage(
titlePanel("Explorador del conjunt de dades Heart"),
sidebarLayout(
sidebarPanel(
selectInput("x", "Selecciona la variable X:", choices = colnames(heart)),
selectInput("y", "Selecciona la variable Y:", choices = colnames(heart)),
sliderInput("bins", "Nombre de barres per l'histograma:", min = 5, max = 50, value = 20) # Afegim una personalització, en que permetem concretar el nombre de barres de l'histograma
),
mainPanel(
# Organitzem els gràfics en pestanyes Tabs
tabsetPanel(
tabPanel("Dispersió", plotOutput("dispersio")),
tabPanel("Línies", plotOutput("linies")),
tabPanel("Histograma", plotOutput("histograma"))
)
)
)
)
server <- function(input, output){
output$dispersio <- renderPlot({
plot(heart[[input$x]], heart[[input$y]],
xlab = input$x, ylab = input$y,
main = paste("Gràfic de Dispersió"))
})
output$linies <- renderPlot({
# Ordenem les dades per l'eix X per a que el gràfic de línies no es creuï de manera caòtica
dades_ordenades <- heart[order(heart[[input$x]]), ]
plot(dades_ordenades[[input$x]], dades_ordenades[[input$y]],
type = "l", xlab = input$x, ylab = input$y,
main = paste("Gràfic de Línies:", input$x, "vs", input$y),
)
})
# Gràfic 3: Histograma sol per a la variable X
output$histograma <- renderPlot({
# Ens assegurem que la variable seleccionada sigui numèrica per fer l'histograma
validate(
need(is.numeric(heart[[input$x]]), "Per favor, selecciona una variable numèrica per a l'histograma.")
)
hist(heart[[input$x]],
breaks = input$bins,
col = "#3498db",
xlab = input$x,
main = paste("Distribució de la variable", input$x))
})
}
shinyApp(ui = ui, server = server)EXERCICI 2
Mostreu un conjunt de dades a partir d’un arxiu .csv a la vostra elecció que mostri també el resum estadístic del conjunt de dades triat. D’altra banda, afegiu la selecció d’una variable del conjunt de dades a partir de la qual es representi un gràfic del tipus histograma.
library(shiny)
library(MASS)
data("Melanoma")
library(DT) ##
## Adjuntando el paquete: 'DT'## The following objects are masked from 'package:shiny':
##
## dataTableOutput, renderDataTableui <- fluidPage(
titlePanel("Explorador de Dades"),
sidebarLayout(
sidebarPanel(
fileInput("file", "Especifica el conjunt de dades a mostrar"),
hr(),
selectInput("variable", "Selecciona una variable per a l'histograma", choices = NULL),
),
mainPanel(
tabsetPanel(
tabPanel("Taula de Dades", DTOutput("table")),
tabPanel("Resum Estadístic",
h4("Resum Estadístic del Conjunt de Dades"),
verbatimTextOutput("resum")),
tabPanel("Histograma",
h4("Representació Gràfica"),
plotOutput("histograma"))
)
)
)
)
server <- function(input, output, session){
dataset <- reactive ({
req(input$file)
read.csv(input$file$datapath)
})
observeEvent(dataset(), {
# Només agafem les columnes que siguin numèriques, que són les vàlides per a un histograma
columnes_numeriques <- names(dataset())[sapply(dataset(), is.numeric)]
updateSelectInput(session, "variable", choices = columnes_numeriques)
})
output$table <- renderDT({
datatable(dataset())
})
output$resum <- renderPrint({
tryCatch({
summary(dataset())
}, error = function(e){
"Per mostrar el resultat, cal carregar previament l'arxiu de dades"
})
})
output$histograma <- renderPlot({
req(input$variable) # Demana que hi hagi una variable seleccionada
dades_columna <- dataset()[[input$variable]]
hist(dades_columna,
xlab = input$variable,
main = paste("Histograma de la variable:", input$variable))
})
}
shinyApp(ui, server)EXERCICI 3
A partir del conjunt de dades Iris, feu una aplicació Shiny per mostrar, a partir d’una de les variables per a seleccionar del conjunt de dades, els gràfics boxplot i barplot.
library(shiny)
data("iris")
ui <- fluidPage(
titlePanel("Explorador Interactiu de Iris: Boxplot i Barplot"),
sidebarLayout(
sidebarPanel(
# Selector de variables (excloem la columna 5 que és 'Species' perquè és el factor de grup)
selectInput("variable", "Selecciona la Variable:",
choices = colnames(iris)[1:4], selected = "Sepal.Length")
),
mainPanel(
# Organitzem els dos gràfics en pestanyes
tabsetPanel(
tabPanel("Gràfic de Caixes (Boxplot)",
plotOutput("box")),
tabPanel("Gràfic de Barres (Barplot)",
plotOutput("bar"))
)
)
)
)
server <- function(input, output) {
# Boxplot
output$box <- renderPlot({
boxplot(iris[[input$variable]] ~ iris$Species,
data = iris,
col = c("#7bf1a8", "#b5e2fa", "#edafb8"),
xlab = "Espècie (Species)",
ylab = input$variable,
main = paste("Distribució de les espècies a IRIS"))
})
# Barplot
output$bar <- renderPlot({
# Calculem les mitjanes de la variable triada per a cada espècie
mitjanes_grup <- tapply(iris[[input$variable]], iris$Species, mean)
barplot(mitjanes_grup,
col = c("#1b9e77", "#d95f02", "#7570b3"),
xlab = "Espècie (Species)",
ylab = paste("Mitjana de", input$variable),
main = paste("Mitjana de les espècies"))
})
}
shinyApp(ui = ui, server = server)EXERCICI 4
Utilitzant el paquet GenomicFeatures, definiu un vector amb la informació de cada gen, la posició inicial i final, la direcció del bri (strand), l’identificador del gen (id) i el tipus de gen (type).
Una orientació per definir algunes d’aquestes dades les podeu trobar a l’Exemple 4. Bioconductor utilitzant el paquet GenomicFeatures, i podem definir la direcció del bri i el tipus de gen de la manera següent:
Mostreu els detalls utilitzant un objecte de tipus GRanges.
library(GenomicFeatures)## Cargando paquete requerido: BiocGenerics## Cargando paquete requerido: generics##
## Adjuntando el paquete: 'generics'## The following objects are masked from 'package:base':
##
## as.difftime, as.factor, as.ordered, intersect, is.element, setdiff,
## setequal, union##
## Adjuntando el paquete: 'BiocGenerics'## 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, is.unsorted, lapply, Map, mapply, match, mget,
## order, paste, pmax, pmax.int, pmin, pmin.int, Position, rank,
## rbind, Reduce, rownames, sapply, saveRDS, table, tapply, unique,
## unsplit, which.max, which.min## Cargando paquete requerido: S4Vectors## Cargando paquete requerido: stats4##
## Adjuntando el paquete: 'S4Vectors'## The following object is masked from 'package:utils':
##
## findMatches## The following objects are masked from 'package:base':
##
## expand.grid, I, unname## Cargando paquete requerido: IRanges##
## Adjuntando el paquete: 'IRanges'## The following object is masked from 'package:grDevices':
##
## windows## Cargando paquete requerido: Seqinfo## Cargando paquete requerido: GenomicRanges## Cargando paquete requerido: AnnotationDbi## Cargando paquete requerido: Biobase## Welcome to Bioconductor
##
## Vignettes contain introductory material; view with
## 'browseVignettes()'. To cite Bioconductor, see
## 'citation("Biobase")', and for packages 'citation("pkgname")'.##
## Adjuntando el paquete: 'AnnotationDbi'## The following object is masked from 'package:MASS':
##
## selectlibrary(GenomicRanges)
# Definim gens inventats per fer l'exercici
cromosomes <- c("chr1", "chrX", "chr6")
inici <- c(10050, 25000, 48000)
final <- c(12100, 26500, 51200)
gene_id <- c("GN1", "HR2", "BCN3")
strand <- c("+", "-", "+")
gene_type <- c("protein_coding", "non-coding", "protein_coding")
# Construim el vector genòmic amb GRanges
vector_gens <- GRanges(
seqnames = cromosomes,
ranges = IRanges(start = inici, end = final),
strand = strand,
id = gene_id,
type = gene_type
)
print(vector_gens)## GRanges object with 3 ranges and 2 metadata columns:
## seqnames ranges strand | id type
## <Rle> <IRanges> <Rle> | <character> <character>
## [1] chr1 10050-12100 + | GN1 protein_coding
## [2] chrX 25000-26500 - | HR2 non-coding
## [3] chr6 48000-51200 + | BCN3 protein_coding
## -------
## seqinfo: 3 sequences from an unspecified genome; no seqlengthsEXERCICI 5
El paquet limma és un paquet de R que s’utilitza per analitzar dades d’expressió gènica i fer anàlisis estadístiques en experiments de bioxips i seqüenciació d’ARN (ARN-Seq). La sigla limma significa linear models for biochip data. A partir de les dades següents:
Feu una anàlisi d’expressió gènica diferencial amb el paquet limma.
library(limma)##
## Adjuntando el paquete: 'limma'## The following object is masked from 'package:BiocGenerics':
##
## plotMAdata <- read.table(header = TRUE, text = "
Gene Sample1 Sample2 Sample3 Sample4 Sample5 Sample6
Gene1 5 6 7 8 9 10
Gene2 4 5 6 7 8 9
Gene3 3 4 5 6 7 8
")
# Definim dos grups, del mateix nombre.
group <- factor(c("Control", "Control", "Control", "Tractament", "Tractament",
"Tractament"))
# Definim el nom de les línies, i ens quedem amb el dataframe excloent la primera columna
rownames(data) <- data$Gene
data_expres <- data[,-1]
data_expres <- as.matrix(sapply(data_expres, as.numeric))
rownames(data_expres) <- data$Gene
design <- model.matrix(~0 + group)
colnames(design) <- levels(group)
# Procedim a fer limma
fit <- lmFit(data_expres, design)
contrast <- makeContrasts(Tractament - Control, levels = design)
fit2 <- contrasts.fit(fit, contrast)
fit2 <- eBayes(fit2)
results <- topTable(fit2, coef = 1)
print(results)## logFC AveExpr t P.Value adj.P.Val B
## Gene1 3 7.5 3.674235 0.003182214 0.003182214 0.2754422
## Gene3 3 5.5 3.674235 0.003182214 0.003182214 0.2754422
## Gene2 3 6.5 3.674235 0.003182214 0.003182214 0.2754422En tenir els tres gens un logFC = 3, voldria dir que els tres s’han sobreexpressat en rebre el tractament, en unes 3 unitats en escala logarítmica. Pel que fa al valor de t, un valor de 3,67 és considerablement alt, pel que els indica que hi ha una separació clara entre controls i tractament.
El p-value ajustat és molt menor que el valor llindar estàndard (0.05), pel que els tres gens s’expressen de forma diferencial de forma molt significativa.
EXERCICI 6
El paquet ComplexHeatMap és útil per visualitzar dades de matriu en forma de mapes de calor complexos amb anotacions addicionals. A partir d’una matriu fictícia, creeu un objecte HeatMap, mostreu-lo i interpreteu-ne els resultats.
library(ComplexHeatmap)## Cargando paquete requerido: grid## ========================================
## ComplexHeatmap version 2.26.1
## Bioconductor page: http://bioconductor.org/packages/ComplexHeatmap/
## Github page: https://github.com/jokergoo/ComplexHeatmap
## Documentation: http://jokergoo.github.io/ComplexHeatmap-reference
##
## If you use it in published research, please cite either one:
## - Gu, Z. Complex Heatmap Visualization. iMeta 2022.
## - Gu, Z. Complex heatmaps reveal patterns and correlations in multidimensional
## genomic data. Bioinformatics 2016.
##
##
## The new InteractiveComplexHeatmap package can directly export static
## complex heatmaps into an interactive Shiny app with zero effort. Have a try!
##
## This message can be suppressed by:
## suppressPackageStartupMessages(library(ComplexHeatmap))
## ========================================library(circlize)## ========================================
## circlize version 0.4.18
## CRAN page: https://cran.r-project.org/package=circlize
## Github page: https://github.com/jokergoo/circlize
## Documentation: https://jokergoo.github.io/circlize_book/book/
##
## If you use it in published research, please cite:
## Gu, Z. circlize implements and enhances circular visualization
## in R. Bioinformatics 2014.
##
## This message can be suppressed by:
## suppressPackageStartupMessages(library(circlize))
## ========================================# Generem una matriu aleatoria amb rnorm
matriu <- matrix(rnorm(64), ncol= 8)
# I utilitzem Heatmap per generar el gràfic
heatmap <- Heatmap(matriu)
draw(heatmap)
EXERCICI 7
El paquet Biostrings de Bioconductor en R està dissenyat per a la manipulació i l’anàlisi de seqüències biològiques. Suposem que definim les seqüències següents:
Combineu les seqüències i obtingueu el complement invers de cada seqüència.
library(Biostrings)## Cargando paquete requerido: XVector##
## Adjuntando el paquete: 'Biostrings'## The following object is masked from 'package:grid':
##
## pattern## The following object is masked from 'package:base':
##
## strsplitsequencia1 <- DNAString("ATCGAATCGA")
sequencia2 <- DNAString("TAGCTAGCTA")
sequencies <- DNAStringSet(c(sequencia1, sequencia2))
inversa <- reverseComplement(sequencies)
print(inversa)## DNAStringSet object of length 1:
## width seq
## [1] 20 TAGCTAGCTATCGATTCGATEXERCICI 8
A partir del paquet caret i del conjunt de dades Iris, creeu un particionament del conjunt de dades sobre la variable species i creeu i entreneu un model de k-NN per fer prediccions. Finalment, mostreu la matriu de confusió i comenteu el resultat.
library(caret)## Cargando paquete requerido: ggplot2## Cargando paquete requerido: lattice##
## Adjuntando el paquete: 'caret'## The following object is masked from 'package:generics':
##
## train## The following object is masked from 'package:survival':
##
## clusterdata("iris")
set.seed(384)
# Dividim les dades
trainIndex <- createDataPartition(iris$Species, p = .8,
list = FALSE,
times = 1)
# Definim a quina variable van les dades dividides
dades_train <- iris[trainIndex,]
dades_test <- iris[-trainIndex,]
# Definim el controlador d'entrenament fent servir el mètode de classificació k-NN
controlador <- trainControl(method = "cv")
model_knn <- train(Species ~ .,
data = dades_train,
method = "knn",
trControl = controlador,
)
prediccions <- predict(model_knn, dades_test)
conf_matrix <- confusionMatrix(prediccions, dades_test$Species)
print(conf_matrix)## Confusion Matrix and Statistics
##
## Reference
## Prediction setosa versicolor virginica
## setosa 10 0 0
## versicolor 0 10 1
## virginica 0 0 9
##
## Overall Statistics
##
## Accuracy : 0.9667
## 95% CI : (0.8278, 0.9992)
## No Information Rate : 0.3333
## P-Value [Acc > NIR] : 2.963e-13
##
## Kappa : 0.95
##
## Mcnemar's Test P-Value : NA
##
## Statistics by Class:
##
## Class: setosa Class: versicolor Class: virginica
## Sensitivity 1.0000 1.0000 0.9000
## Specificity 1.0000 0.9500 1.0000
## Pos Pred Value 1.0000 0.9091 1.0000
## Neg Pred Value 1.0000 1.0000 0.9524
## Prevalence 0.3333 0.3333 0.3333
## Detection Rate 0.3333 0.3333 0.3000
## Detection Prevalence 0.3333 0.3667 0.3000
## Balanced Accuracy 1.0000 0.9750 0.9500El model té una precisió o encert del 96,67%, i un valor kappa de 95%, pel que molt probablement hi ha una concordança significativa entre les prediccions i les espècies originals.
El model té una sensitivitat i especificitat perfecta per setoses, sensitivitat perfecta en versicolor (i especificitat molt elevada, del 95%) i especificitat perfecta en virginica (amb sensibilitat alta del 90%).
Per tant, aquest model k-NN sembla bastant adequat per aquesta tasca.
EXERCICI 9
Suposem que tenim les següents dades de 5 pacients als quals se’ls han practicat algunes proves mèdiques i on, entre altres dades, s’indica si tenen o no la malaltia estudiada:
Es demana entrenar un model de classificació supervisada, per exemple, el random forest.
library(caret)
library(randomForest)## randomForest 4.7-1.2## Type rfNews() to see new features/changes/bug fixes.##
## Adjuntando el paquete: 'randomForest'## The following object is masked from 'package:ggplot2':
##
## margin## The following object is masked from 'package:Biobase':
##
## combine## The following object is masked from 'package:BiocGenerics':
##
## combinedades <- data.frame(
Edat = c(30, 45, 50, 22, 65),
Sexe = c("M", "F", "M", "F", "M"),
ResultatProva1 = c(0.1, 0.5, 0.8, 1.2, 0.9),
ResultatProva2 = c(12, 9, 15, 8, 10),
Malaltia = c("No", "No", "Sí", "Sí", "Sí"))
# Establim una seed per a que els resultats siguin reproduibles
set.seed(596)
dades$Sexe <- as.factor(dades$Sexe)
controlador <- trainControl(method = "cv", number = 5)
model_forest <- train(Malaltia ~ ., data = dades, method = "rf", trControl = controlador)## Warning in nominalTrainWorkflow(x = x, y = y, wts = weights, info = trainInfo,
## : There were missing values in resampled performance measures.print(model_forest)## Random Forest
##
## 5 samples
## 4 predictors
## 2 classes: 'No', 'Sí'
##
## No pre-processing
## Resampling: Cross-Validated (5 fold)
## Summary of sample sizes: 4, 4, 4, 4, 4
## Resampling results across tuning parameters:
##
## mtry Accuracy Kappa
## 2 0.4 0
## 3 0.4 0
## 4 0.6 0
##
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was mtry = 4.EXERCICI 10
Suposem que tenim la següent matriu d’adjacència que representa les interaccions entre proteïnes:
Es demana representar gràficament la xarxa d’interaccions entre aquestes proteïnes.
library(igraph)##
## Adjuntando el paquete: 'igraph'## The following object is masked from 'package:Biostrings':
##
## union## The following object is masked from 'package:XVector':
##
## path## The following object is masked from 'package:circlize':
##
## degree## The following object is masked from 'package:GenomicRanges':
##
## union## The following object is masked from 'package:IRanges':
##
## union## The following object is masked from 'package:S4Vectors':
##
## union## The following objects are masked from 'package:BiocGenerics':
##
## normalize, path, union## The following objects are masked from 'package:generics':
##
## components, union## The following objects are masked from 'package:stats':
##
## decompose, spectrum## The following object is masked from 'package:base':
##
## unionmatrix_adjacencia <- matrix(c( 0, 1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 0, 0, 1, 0, 0), nrow= 4, byrow = TRUE)
graph_matrix <- graph_from_adjacency_matrix(matrix_adjacencia, mode = "undirected")
plot(graph_matrix, main="Gràfic d'interaccions entre proteïnes")
CASOS PRÀCTICS
CAS PRÀCTIC 1
A partir de Bioconductor, i utilitzant el paquet Biostrings, es demana fer una aplicació senzilla amb Shiny que sigui un comptador d’ocurrències de patrons d’ADN. Es requereix que la interfície que s’obtingui sigui similar a la següent:
Per fer proves, es poden fer servir els valors següents: seqüència: ATCGATCGATCGATCGATCGATCGATCGATCG patró: ATCG.
library(shiny)
library(Biostrings)
ui <- fluidPage(
titlePanel("Contador d'ocurrències de patrons d'ADN"),
sidebarLayout(
sidebarPanel(
textInput("sequencia", "Introdueix una seqüència d'ADN:"),
textInput("diana", "Introdueix el patró a buscar:"),
actionButton("buscar", "Contar ocurrències")
),
mainPanel(
h4("Resultat del recompte de patrons:"),
verbatimTextOutput("resultat")
)
)
)
server <- function(input, output){
recompte <- eventReactive(input$buscar, {
req(input$sequencia, input$diana) # Requerim que els dos camps tinguin text
seq_adn <- DNAString(input$sequencia)
patro_adn <- DNAString(input$diana)
num_ocurrencies <- countPattern(patro_adn, seq_adn)
return(num_ocurrencies)
})
output$resultat <- renderPrint({
if (input$buscar == 0) {
return("Introdueix les dades i fes clic al botó per iniciar el cerca-patrons.")
}
cat("El patró '", input$diana, "' s'ha trobat un total de:", recompte(), "vegades.")
})
}
shinyApp(ui=ui, server=server)