#!/usr/bin/R
# ------------------------------------------
# Load expression data for the joint cohort.
# Updated: 08/30/22
# ------------------------------------------
library(cbrbase)
set_proj('DEVTRAJ', 'mosaicism')


# Load full expression matrices (from both cohorts):
# --------------------------------------------------
mostrecent = '20211104' # Fixed for paper
merged.mat.rds = paste0(dbdir, 'merged_ALZ_cohorts_expression_matrix_', mostrecent, '.rds')
merged.mati.rds = paste0(dbdir, 'merged_ALZ_cohorts_expression_matrix_intr_', mostrecent, '.rds')
if (!file.exists(merged.mat.rds)){
    fullmat = NULL
    fullmati = NULL
    for (project in c('ALZ', 'ALZ_repl')){
        commandArgs <- function(trailingOnly=TRUE){ c(mostrecent) }
        source(paste0(bindir, 'mosaicism/load_count_data.R'))
        cat(project, dim(allmat), '\n')
        if (is.null(fullmat)){
            fullmat = allmat 
            fullmati = allmati
        } else { 
            fullmat = cbind(fullmat, allmat)
            fullmati = cbind(fullmati, allmati)
        }
        rm(allmat, allmati)
    }
    saveRDS(fullmat, file=merged.mat.rds)
    saveRDS(fullmati, file=merged.mati.rds)
} else {
    fullmat = readRDS(merged.mat.rds)
    fullmati = readRDS(merged.mati.rds)
}


# Load in and process metadata
# -------------------------------------------
metalist = list()
for (project in c('ALZ', 'ALZ_repl')){
    commandArgs <- function(trailingOnly=TRUE){ c(project) }
    source('~/DEVTRAJ/bin/mosaicism/load_metadata.R')
    meta$cohort = project
    metalist[[project]] = meta
}
mcols = lapply(metalist, colnames)
keep.cols = intersect(mcols[[1]], mcols[[2]])
metalist = lapply(metalist, function(x){ x[, keep.cols] })
meta = do.call(rbind, metalist)