library( tidyverse )
library( survival )
library( openxlsx )

## Load the data and ABSOLUTE estimates (syn3582761)
X <- read.xlsx("TCGA_MEASUREMENTS.xlsx")
P <- read.xlsx("Purity_Ploidy_All_Samples_4-17-15.xlsx") %>%
    rename( ID = individual_id ) %>% mutate_at( "purity", as.numeric )
XP <- inner_join( X, P, by="ID" )

## Compute correlations to establish internal consistency
with( X, cor.test( CAF_SCORE, percent_stromal_cells ) )    # cor = 0.20, p-value < 0.001
with( XP, cor.test( CAF_SCORE, purity ) )                  # cor = -0.54, p-value < 0.001
with( XP, cor.test( percent_stromal_cells, purity ) )      # cor = -0.12, p-value = 0.007

## Compute correlations with VGG19 predictions
with( X, cor.test( CAF_SCORE, STR ) )                      # cor = 0.26, p-value < 0.001
with( XP, cor.test( TUM, purity ) )                        # cor = 0.069, p-value = 0.14

## Compute an unbiased estimate of non-tumor components
X2 <- X %>% mutate( DeepScore = ADI + DEB + LYM + MUS + STR,
                   OS_event = vital_status, days_to_event = days_to_event/365.25 )

## Determine if there's any relationship between the unbiased score and survival
## Compare it to a similar analysis for the stromal component STR only
coxph( Surv(days_to_event, OS_event) ~ DeepScore, data = X2 )  # p = 0.01
coxph( Surv(days_to_event, OS_event) ~ STR, data = X2 )        # p = 0.85