Relatório.rmd

---
title: "Adenocarcinoma prostático"
author: "Afonso Sá (pg54434), Armindo Machado (pg52170)"
output:
  html_document:
    theme: cerulean
    toc: true
    toc_float: true
  pdf_document: default
encoding: "UTF-8"
editor_options:
  markdown:
    wrap: auto

```
```{r setup, include=FALSE, warning = FALSE, message=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

# Exploração inicial dos dados e análise de expressão diferencial e de enriquecimento

```{r,echo=FALSE,results='hide', warning=FALSE, message=FALSE}
#install.packages("jpeg")
#install.packages("TCGAbiolinks")
#install.packages("Biobase")
#install.packages("DESeq2")
#install.packages("ggbeeswarm")
#install.packages("genefilter")
#install.packages("pheatmap")
#install.packages("org.Hs.eg.db")
#install.packages("fgsea")
#install.packages("ggplot2")
#install.packages("xfun")
#install.packages("SummarizedExperiment") 
#install.packages("factoextra")
#install.packages("Rtsne")
#install.packages("EBImage")
#install.packages("BiocManager")
```

```{r,echo=FALSE,results='hide', warning = FALSE, message=FALSE,warning = FALSE, message=FALSE}
#library("jpeg")
#library("DESeq2")
#library("Biobase")
#library("DESeq2")
#library("ggbeeswarm")
#library("pheatmap")
#library("org.Hs.eg.db")
#library("fgsea")
#library("ggplot2")
#library("factoextra")
#library("limma")
#library("genefilter")
#library("SummarizedExperiment")
#library("Rtsne")
#library("party")
#library("rpart")
#library("caret") 
#library("rsample")
#library(BiocManager)

```

## Introdução
O cancro da próstata é um tipo comum de cancro que afeta pessoas do sexo masculino, incidindo sobre a glândula da próstata e provocando dificuldades urinárias. O diagnóstico envolve exames de PSA e biópsia da próstata. O tratamento pode incluir cirurgia, radioterapia, terapia hormonal ou quimioterapia. O adenocarcinoma é o tipo mais comum de cancro da próstata, originando-se nas células glandulares. Quando diagnosticado numa fase inicial, apresenta uma taxa de sobrevivência de 98%, contrastando com apenas 30% numa fase avançada.

## Origem dos dados
Os dados analisados neste trabalho provêm do projeto Pan-Cancer Atlas (TCGA), uma colaboração em larga escala apoiada pelo National Cancer Institute (NCI) e pelo National Human Genome Research Institute (NHGRI). Serão analisadas 494 amostras de mais de 11.000 tumores de 33 tipos diferentes de cancro, obtidas através de diversas técnicas que examinam alterações moleculares a vários níveis, incluindo DNA, RNA, proteínas e epigenética. A análise destes perfis moleculares visa identificar padrões comuns, diferenças e temas emergentes entre os diferentes tipos de cancro, contribuindo para o desenvolvimento de tratamentos mais eficazes e personalizados.

## Extração dos dados de transcriptómica e respetivos metadados da base de dados

Para analisar dados de transcriptómica do adenocarcinoma da próstata, foi extraída informação do cBioPortal usando o pacote TCGAbiolinks. Este método permite aceder aos dados do The Genomic Data Commons (GDC), suportando a análise de dados genómicos, clínicos e de biospecimens de programas de investigação sobre cancro.

```{r,results= 'markup'}
# Criar uma consulta ao Genomic Data Commons (GDC) e obter dados transcriptómicos do projeto TCGA sobre câncer da próstata (PRAD)
library(TCGAbiolinks)
query_PRAD <- GDCquery(
  project = "TCGA-PRAD",
  data.category = "Transcriptome Profiling", 
  data.type = "Gene Expression Quantification",
  workflow.type = "STAR - Counts"
)

# Baixar os dados do Genomic Data Commons (GDC) com base nas especificações de uma consulta criada anteriormente
GDCdownload(query = query_PRAD)

# Preparar os dados do Genomic Data Commons (GDC) para análise
rna_seq_PRAD <- GDCprepare(query = query_PRAD, save = TRUE, save.filename = "mRNA_TCGA-PRAD.rda") 
```
```{r, data_verify}
# Verificação dos dados
class(rna_seq_PRAD)
dim(rna_seq_PRAD)

# Exibir os primeiros registros de metadados
head(colData(rna_seq_PRAD))
```

A função GDCquery() é utilizada para realizar a query à base dados. Nesta função são definidos vários parâmetros que permitem filtrar o data set que posteriormente vai ser descarregado. O parâmetro "project" corresponde ao id atribuído ao tipo de cancro de modo a apenas considerar a informação disponível referente ao carcinoma endométrio do corpo uterino. A "data.category" é usado para especificar a categoria de dados que se pretende selecionar, neste caso, como pretendemos realizar a análise da transcriptómica selecionou-se a categoria "Transcriptome Profiling" para selecionar os dados relacionados com o perfil de expressão de genes. De entre os dados do perfil de expressão de genes restringiu-se apenas aos dados referentes à quantificação da expressão de genes definindo "data.type" como Gene expression quantification. Por fim, indicou-se o tipo de análise da expressão genética desejada, optando-se por recorrer ao START-Counts. Este workflow indica que o algoritmo utilizado para a contagem da expressão genética a partir do alinhamento das sequências foi realizada com a ferramenta STAR.

A função GDCdownload() é utilizada para realizar o download dos dados consoante a query realizada e estabelecida na função anterior

```{r, GDC}
#Retorna a classe do objeto (tipo de dados ou a estrutura de dados que o objeto representa, o que por sua vez determina quais funções podem ser aplicadas a ele)
class(rna_seq_PRAD)
#Retorna as dimensões do objeto(matriz ou um dataframe;objeto mais complexo)
dim(rna_seq_PRAD)
#Visualizar informações básicas do conjunto de dados
metadata(rna_seq_PRAD)

rownames(rna_seq_PRAD)
colnames(rna_seq_PRAD)

#Extrair componentes de um objeto por nome (através de colunas)
rna_seq_PRAD$paper_BCR
rna_seq_PRAD$paper_Gender
rna_seq_PRAD$paper_Grade

#Atribui a um novo objeto chamado meta_PRAD, os metadados associados ao conjunto de dados rna_seq_PRAD
meta_PRAD = colData(rna_seq_PRAD)
#Retorna as dimensões dos metadados 
dim(meta_PRAD)
#Extrair componentes de um objeto por nome (através de colunas)
meta_PRAD$patient
meta_PRAD$paper_IDH.status
```

Extração dos Metadados

```{r}
#Extrair metadados relacionados com as linhas
row_metadados <- SummarizedExperiment::rowData(rna_seq_PRAD)[1:10, c('type', 'gene_type', 'gene_name')]

#Atribuir a um novo objeto chamado meta_PRAD, os metadados associados aos dados de expressão
amostras_metadados <- colData(rna_seq_PRAD)
class(amostras_metadados)
dim(amostras_metadados)
colnames(amostras_metadados)

#Código utilizado para extrair a informação relacionada à contagem da expressão dos genes do objeto rna_seq_PRAD
geneExp <- SummarizedExperiment::assay(rna_seq_PRAD)
#Para ver os outliers - no nosso caso é zero
pre <- TCGAanalyze_Preprocessing(rna_seq_PRAD) # Faz correlação e tenta identificar outliers
dim(geneExp)
dim(pre)
```

## Pré-processamento e Filtragem 

```{r}
# Retirar as colunas dos metadados onde havia mais de 50 elementos como: “not reported/Not Reported/Reported” e/ou “NA”

cols_with_not_reported <- which(sapply(meta_PRAD, function(x) sum(x == "not reported", na.rm = TRUE)) > 50)
cols_with_Not_Reported <- which(sapply(meta_PRAD, function(x) sum(x == "Not Reported", na.rm = TRUE)) > 50)
cols_with_NA <- which(sapply(meta_PRAD, function(x) sum(is.na(x))) > 60)

# Remover as colunas baseadas nos critérios específicos de cima
metadata_matriz_clean <- meta_PRAD[, -c(cols_with_not_reported, cols_with_Not_Reported, cols_with_NA)]
dim(metadata_matriz_clean)

# Código utilizado para extrair a informação relacionada à contagem da expressão dos genes do objeto rna_seq_PRAD
geneExp <- SummarizedExperiment::assay(rna_seq_PRAD)

#Carregamento dos metadados
meta_PRAD = colData(rna_seq_PRAD)

#Os metadados são convertidos em um dataframe para facilitar a manipulação
meta_PRAD = as.data.frame(meta_PRAD)

#Retorna as dimensões dos metadados 
dim(meta_PRAD)
#Ver nomes das colunas
colnames(meta_PRAD)
#Nomes das linhas
names(meta_PRAD)
#Extrair componentes de um objeto por nome (através de colunas)
meta_PRAD$patient

#Criar vetor lógico em que diz que linhas possuem grade diferente de NA
coluna <- !is.na(meta_PRAD$definition)

#Selecionar apenas as linhas que têm grade
meta <- meta_PRAD[coluna, c("definition", "age_at_diagnosis", "disease_type")]
meta$disease_type

#Confirmar que os pacientes selecionados têm um valor de grade atribuído na tabela de metadados original
table(meta_PRAD[rownames(meta), "paper_tumor_grade"]) 
head(meta) 

#Verificação da dimensão do dataframe meta após seleção de apenas as linhas com informação relativa ao grade
sum(table(meta_PRAD$definition))
dim(meta_PRAD)

#Criação do dataframe de expressão apenas para os pacientes com informação relativa ao grade
exp_grade <- geneExp[, rownames(meta)]
dim(exp_grade)

#Nomes das linhas
names(meta_PRAD)
#Extrair componentes de um objeto por nome (através de colunas)
meta_PRAD$patient
meta_PRAD$initial_weight
```


## Análise de Expressão Diferencial com DESeq2
A análise de expressão diferencial é crucial para entender as mudanças na expressão gênica entre diferentes condições ou tratamentos. O pacote DESeq2 é amplamente utilizado para esta análise em dados de contagem de sequenciamento de RNA.

```{r, deseq2}
library(DESeq2) # DESeq2, uma ferramenta para análise de expressão diferencial de dados de contagem de sequenciamento de RNA (RNA-Seq)

#Filtra os dados de RNA para incluir apenas amostras não nulas
data_de <- rna_seq_PRAD[, !is.na(rna_seq_PRAD$colData$vital_status)]
table(as.data.frame(metadata_matriz_clean$vital_status))

countData <- assays(data_de)
#Cria um objeto DESeqDataSet para análise, especificando um design experimental que compara o status vital
ddsSE <- DESeqDataSetFromMatrix(countData = countData, 
                                colData = colData(data_de), 
                                design = ~ vital_status)
#Filtragem de Genes: Remove genes com contagens baixas (menos de 10) para melhorar a confiabilidade da análise de expressão diferencial.
keep <- rowSums(counts(ddsSE)) >= 10
ddsSE <- ddsSE[keep,]

#Executa a análise de expressão diferencial com a função DESeq
ddsSE <- DESeq(ddsSE)

#Visualizar nomes dos resultados disponíveis
resultsNames(ddsSE)

#Obter e visualizar resultados da comparação de interesse
res <- results(ddsSE)
dea <- as.data.frame(res)

#Resume os resultados para obter uma visão geral dos achados estatísticos, como o número de genes significativamente diferencialmente expressos
summary(res)
```


## Processamento de Dados Clínicos com TCGAbiolinks

Dados clínicos são essenciais para complementar as análises genômicas, permitindo uma compreensão mais aprofundada das características dos pacientes e das amostras analisadas.

```{r, TCGAbiolinks}

#Buscar dados clínicos do projeto TCGA-PRAD (adenocarcinoma prostático), especificamente suplementos clínicos, através de GDCquery
query_clin <- GDCquery(project = "TCGA-PRAD", 
                       data.category = "Clinical",
                       data.type = "Clinical Supplement",
                       data.format = "BCR Biotab")

#Baixa e prepara os dados clínicos
GDCdownload(query_clin)
clinical.PRAD <- GDCprepare(query_clin)

#Listar os nomes dos componentes do objeto
names(clinical.PRAD)

#Verificar os dados clínicos disponíveis
head(clinical.PRAD)


```

## Análise Adicional com Dados Clínicos e de Amostra 

```{r, include=FALSE}
#Verifica se o ggplot2 está instalado, senão, instala
if (!requireNamespace("ggplot2", quietly = TRUE)) {
  install.packages("ggplot2")
}
library(ggplot2)
```

Carregamento e limpeza dos dados clínicos e de amostra:

```{r}
#Carrega dados de um arquivo de texto em um formato de dataframe.
dados_sample <- read.delim("C:/Users/armin/OneDrive/Ambiente de Trabalho/dataset/data_clinical_sample.txt", header = TRUE, stringsAsFactors = FALSE)
dados_patient <- read.delim("C:/Users/armin/OneDrive/Ambiente de Trabalho/dataset/data_clinical_patient.txt", header = TRUE, stringsAsFactors = FALSE)

#Criação de novos conjuntos de dados removendo as primeiras 4 linhas dos conjuntos originais.
dados_sample_limpo <- dados_sample[-c(1:4), ]
dados_patient_limpo <- dados_patient[-c(1:4), ]

#Conversão do conjunto de dados para o tipo de dados numéricos.
dados_patient_limpo$Overall.Survival..Months. <- as.numeric(as.character(dados_patient_limpo$Overall.Survival..Months.))
dados_patient_limpo$Diagnosis.Age <- as.numeric(as.character(dados_patient_limpo$Diagnosis.Age))
dados_sample_limpo$MSI.MANTIS.Score <- as.numeric(dados_sample_limpo$MSI.MANTIS.Score)
dados_sample_limpo$TMB..nonsynonymous. <- as.numeric(dados_sample_limpo$TMB..nonsynonymous.)
```


Visualização inicial dos dados:

```{r}
#Visualizar os primeiros registros do conjunto de dados.
head(dados_sample_limpo)
head(dados_patient_limpo)

#Obter informações sobre a estrutura do conjunto de dados.
str(dados_sample_limpo)
str(dados_patient_limpo)

#Resumo estatístico do conjunto de dados.
summary(dados_sample_limpo)
summary(dados_patient_limpo)
```


Análise exploratória com gráficos:

```{r}
#Realiza uma ANOVA para comparar a idade no diagnóstico entre diferentes categorias
anova <- aov(meta_PRAD$age_at_diagnosis ~ meta_PRAD$definition)
summary(anova)

#Realiza o teste de Tukey para comparações múltiplas
TukeyHSD(anova)
```

```{r}

boxplot(meta_PRAD$age_at_diagnosis ~ meta_PRAD$definition, horizontal = TRUE)
```


```{r}
#Cria um gráfico de dispersão relacionando o MSI e o TMB.
ggplot(dados_sample_limpo, aes(x = MSI.MANTIS.Score, y = TMB..nonsynonymous.)) +
  geom_point() +
  labs(title = "Relação entre MSI e TMB", x = "Score MSI", y = "Score TMB")

#Calcula o coeficiente de correlação de Pearson.
correlation1 <- cor(dados_sample_limpo$MSI.MANTIS.Score, dados_sample_limpo$TMB..nonsynonymous.)

#Exibe o valor da correlação.
print(correlation1)
```

PCA plot

```{r}
library(ggplot2)
library(dplyr)
library(readr)

# carregar os dados de expressão genética
dados_expr_genetica <- read_csv("bla/bla/bla/nome do ficheiro")

# remover colunas não numéricas se houver
dados_expr_genetica_numeric <- select_if(dados_expr_genetica, is.numeric)

pca_result <- prcomp(dados_expr_genetica_numeric, scale. = TRUE)

#extrair os componentes principais
pca_df <- as.data.frame(pca_result$x)

#colorir os pontos no plot, substituir "coluna" pelo nome da coluna
pca_df$Coluna <- dados_expr_genetica$Coluna

#plotar o gráfico de PCA
pca_plot <- ggplot(pca_df, aes(x = PC1, y = PC2, color = Coluna)) +
  geom_point(size = 3) +
  labs(x = "PC1", y = "PC2", title = "Plot de PCA com Dados de Expressão Genética")

#gráfico
print(pca_plot)

```


Análise de sobrevivência e correlação:

```{r}
#Verifica valores únicos na coluna Overall.Survival.Status.
unique(dados_patient_limpo$Overall.Survival.Status)

#Converte os estados de sobrevivência para numéricos.
dados_patient_limpo$Overall.Survival.Status <- ifelse(dados_patient_limpo$Overall.Survival.Status == "LIVING", 0, 1)

#Ajusta a curva de sobrevivência de Kaplan-Meier.
library(survival)
library(survminer)
surv_obj <- Surv(dados_patient_limpo$Overall.Survival..Months., dados_patient_limpo$Overall.Survival.Status)
surv_fit <- survfit(surv_obj ~ 1)

#Traça a curva de sobrevivência.
ggsurvplot(surv_fit, data = dados_patient_limpo, risk.table = TRUE, xlab = "Tempo (Meses)", ylab = "Probabilidade de Sobrevivência", title = "Curva de Sobrevivência")
```


Análise de regressão linear:

```{r}
#Ajusta o modelo de regressão linear.
modelo <- lm(Diagnosis.Age ~ Overall.Survival..Months., data = dados_patient_limpo)

#Cria o gráfico de dispersão com a linha de regressão.
ggplot(dados_patient_limpo, aes(x = Overall.Survival..Months., y = Diagnosis.Age)) +
  geom_point() +
  geom_smooth(method = "lm", col = "red") +
  labs(title = "Gráfico de Regressão", x = "Sobrevivência Global (Meses)", y = "Idade do Diagnóstico")
```

Clustering dos genes

```{r}
# Ler o arquivo de texto com read.delim
dados_genes <- read.delim("C:/Users/afons/Downloads/Structural_Variant_Genes.txt", header = TRUE)

# Visualizar as primeiras linhas do conjunto de dados
head(dados_genes)

# Obter informações sobre a estrutura do conjunto de dados
str(dados_genes)

# Resumo estatístico do conjunto de dados
summary(dados_genes)

# Carregar o pacote ggplot2
library(ggplot2)

# Selecionar apenas os top 20 genes
dados_top20 <- dados_genes[1:20, c("Gene", "X..Structural.Variant")]

# Remover linhas com valores NA se houver
dados_top20 <- na.omit(dados_top20)

# Aplicar o algoritmo de clustering (k-means)
k <- 4  # Número de clusters desejado
set.seed(123) 
kmeans_result <- kmeans(dados_top20$X..Structural.Variant, centers = k)

# Adicionar os resultados do k-means aos dados originais
dados_top20$cluster <- as.factor(kmeans_result$cluster)

# Criar o gráfico de dispersão com os clusters destacados por cores
grafico_clusters <- ggplot(dados_top20, aes(x = Gene, y = X..Structural.Variant, color = cluster)) +
  geom_point() +
  labs(title = "Clusters dos Top 20 Genes",
       x = "Gene",
       y = "X..Structural.Variant") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) +
  scale_color_discrete(name = "Cluster")  # Definir legenda para as cores dos clusters

# Exibir o gráfico
print(grafico_clusters)
```

Clustering dos genes mutados

```{r}
# Ler o arquivo de texto com read.delim
dados_mutados <- read.delim("C:/Users/afons/Downloads/Mutated_Genes.txt", header = TRUE)

# Visualizar as primeiras linhas do conjunto de dados
head(dados_mutados)

# Obter informações sobre a estrutura do conjunto de dados
str(dados_mutados)

# Resumo estatístico do conjunto de dados
summary(dados_mutados)

# Selecionar apenas os top 20 genes mutados
dados_mutados_top20 <- dados_mutados[1:20, ]

# Aplicar o algoritmo de clustering (k-means)
k <- 4  # Número de clusters desejado
set.seed(123)  
kmeans_result <- kmeans(dados_mutados_top20$X..Mut, centers = k)

# Adicionar os resultados do k-means aos dados originais
dados_mutados_top20$cluster <- as.factor(kmeans_result$cluster)

# Criar o gráfico de dispersão com os clusters destacados por cores
grafico_clusters_mutados_top20 <- ggplot(dados_mutados_top20, aes(x = Gene, y = X..Mut, color = cluster)) +
  geom_point() +
  labs(title = "Clusters dos Top 20 Genes Mutados",
       x = "Gene",
       y = "X..Mut") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) +
  scale_color_discrete(name = "Cluster")  # Definir legenda para as cores dos clusters

# Exibir o gráfico
print(grafico_clusters_mutados_top20)
```