Confidence interval

[Javad]

Hi Sam 

Is there a way to calculate 95% confidence intervals for the out R2s?

[Sam Choi]

I haven't implemented it, but I guess you can do something like boostrapping in R to do so.

something like (not tested)

library(boot)
library(data.table) # For its speed
library(magrittr) # For piping
library(fmsb) # For nagelkerkeR2
prs <-
    fread("PRSice.best")[In_Regression == "Yes"]  # Change it to the best score file, filtering by In_Regression
pheno <- fread("Phenotype") # change this to your phenotype
covariates <- fread("Covariates") # change this to your covariate

data <- merge(pheno, prs, by = c("FID", "IID") %>% # Assume all structures has FID and IID

        merge(covariates, by = c, "IID")
    result <- boot(
            data = data,
            statistic = function(data, indices) {
                d <-data[indices, ]
                full.r2 <- summary(lm(Pheno ~PRS + Covariate, d))$r.squared
                null.r2 <- summary(lm(Pheno ~ Covariate, d))$r.squared
                return(full.r2 - null.r2)
            },
            R = 2000
        )
    boot.ci(result, type = "bca")

[Javad]

Hi Sam 

Thank you for the code, it worked ( with some minor changes), however, the intervals are too broad. Do you know why?

BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS

Based on 2000 bootstrap replicates

CALL : 

boot.ci(boot.out = result, type = "bca")

Intervals : 

Level       BCa          

95%   ( 0.0028,  0.0345 )  

Calculations and Intervals on Original Scale

[Sam Choi]

That depends on your sample size. And given that your estimated R2 is likely less than 0.03 (judging from your confidence interval), then it isn't too surprising, as the PRS seems to be rather weak.

[Javad]

Yes, my sample size is relatively small, as well as the R2!

Thanks Sam.