make_data.Rmd

---
title: "TCGA Breast"
output: html_notebook
---

Load the data previously-extracted from TCGA

```{r eval=FALSE}
library(TCGAbiolinks)
# get clinical data
information.clinical <- GDCquery_clinic(project = "TCGA-BRCA",type = "clinical") 

query <- GDCquery(project = "TCGA-BRCA",
                  data.category = "Clinical",
                  data.type = "Clinical Supplement",
                  data.format = "BCR Biotab")
GDCdownload(query)
clinical.all <- GDCprepare(query)
tcga_brca.clin <- clinical.all$clinical_patient_brca

readr::write_csv(tcga_brca.clin, "tcga_brca.clin.csv")

```


```{r}
library(dplyr)
library(SummarizedExperiment)
load("brcaExp.rda")
subtype <- read.csv('information_subtype.csv')
clin <- read.csv("tcga_brca.clin.csv") %>% 
  select(patient = bcr_patient_barcode,
         er_status_by_ihc,
         her2_status_by_ihc,
         pr_status_by_ihc)


table(data$shortLetterCode)
```

```{r}
set.seed(270624)
rand_t <- sample(which(data$shortLetterCode == "TP"),20)

rand_normal <- sample(which(data$shortLetterCode == "NT"),20)

data <- data[,c(rand_t,rand_normal)]

meta <- colData(data) %>% data.frame %>% 
  select(patient, shortLetterCode) %>% 
  left_join(subtype) %>% 
  select(patient,shortLetterCode,vital_status,age_at_initial_pathologic_diagnosis,BRCA_Subtype_PAM50,pathologic_stage,BRCA_Pathology,days_to_death,days_to_last_followup) %>% left_join(clin)



## find which IDs can be mapped to SYMBOLs

library(org.Hs.eg.db)

df <- data.frame(id = rownames(data)) %>% 
  tidyr::separate(id, into = c("ENSEMBL","Version"),remove = FALSE)

anno <- AnnotationDbi::select(org.Hs.eg.db, 
                              keys = df$ENSEMBL,
                              columns = "SYMBOL",
                              keytype = "ENSEMBL"
                              )

id_to_use <- left_join(df, anno) %>% 
  filter(!is.na(SYMBOL)) %>% 
  pull(id) %>% 
  unique

data <- data[id_to_use,]

counts <- assay(data,"unstranded") %>% data.frame %>% 
  tibble::rownames_to_column("ens") %>% 
  tidyr::separate(ens, into=c("ENSEMBL","V")) %>% 
  dplyr::select(-V) %>% 
  filter(!duplicated(ENSEMBL))

rownames(meta) <- colnames(counts)[-1]
write.table(counts, "brca_example.tsv",sep = "\t")
write.table(meta, "brca_example_meta.tsv", sep = "\t")
```

## Make ESR1 example for ML


```{r}
library(dplyr)
library(SummarizedExperiment)
load("brcaExp.rda")
subtype <- read.csv('information_subtype.csv')
clin <- read.csv("tcga_brca.clin.csv") %>% 
  dplyr::select(patient = bcr_patient_barcode,
         er_status_by_ihc,
         her2_status_by_ihc,
         pr_status_by_ihc,
         ethnicity)

```

```{r}

meta <- colData(data) %>% data.frame %>% 
  dplyr::select(patient, shortLetterCode) %>% 
  left_join(subtype) %>% 
  dplyr::select(patient,shortLetterCode,vital_status,age_at_initial_pathologic_diagnosis,BRCA_Subtype_PAM50,pathologic_stage,BRCA_Pathology,days_to_death,days_to_last_followup,pathologic_stage) %>% left_join(clin)



## find which IDs can be mapped to SYMBOLs

library(org.Hs.eg.db)

df <- data.frame(id = rownames(data)) %>% 
  tidyr::separate(id, into = c("ENSEMBL","Version"),remove = FALSE)

anno <- AnnotationDbi::select(org.Hs.eg.db, 
                              keys = df$ENSEMBL,
                              columns = "SYMBOL",
                              keytype = "ENSEMBL"
                              )

id_to_use <- left_join(df, anno) %>% 
  filter(!is.na(SYMBOL)) %>% 
  pull(id) %>% 
  unique

data <- data[id_to_use,]

library(DESeq2)
library(genefu)
data("pam50")
## read the data but don't need a proper design

dds <- DESeqDataSetFromMatrix(assay(data, "unstranded"),
                              meta,
                              ~1)

vsd <- vst(dds)

genes_to_use <- rownames(pam50$centroids.map)
anno <- AnnotationDbi::select(org.Hs.eg.db, 
                      keys = genes_to_use,
                      columns = "ENSEMBL",
                      keytype =  "SYMBOL") 
ens_ids <- pull(anno,"ENSEMBL") %>% na.omit()

tmp <- data.frame(ens = rownames(vsd)) %>% 
  tidyr::separate(ens,into=c("ENSEMBL","V"),remove = FALSE)
rownames(vsd) <- tmp$ENSEMBL  

ens_ids <- ens_ids[ens_ids %in% rownames(vsd)]
norm_data <- assay(vsd)[ens_ids,] %>% 
  data.frame() %>% 
  tibble::rownames_to_column("ENSEMBL") %>% 
  tidyr::pivot_longer(-ENSEMBL, names_to = "SAMPLE",values_to = "COUNT") %>% 
  left_join(anno) %>% 
  dplyr::select(-ENSEMBL) %>% 
  tidyr::pivot_wider(everything(), names_from = SYMBOL,values_from = COUNT)

export_data <- bind_cols(norm_data, meta) %>% 
  filter(shortLetterCode == "TP")

saveRDS(export_data, "brca_ML_example.rds")

```