Blood cells regions

Last updated: 2020-11-22

Checks: 7 0

Knit directory: finemap-uk-biobank/

This reproducible R Markdown analysis was created with workflowr (version 1.6.2). The Checks tab describes the reproducibility checks that were applied when the results were created. The Past versions tab lists the development history.

R Markdown file: up-to-date

Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.

Environment: empty

Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it's best to always run the code in an empty environment.

Seed: set.seed(20191114)

The command set.seed(20191114) was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.

Session information: recorded

Great job! Recording the operating system, R version, and package versions is critical for reproducibility.

Cache: none

Nice! There were no cached chunks for this analysis, so you can be confident that you successfully produced the results during this run.

File paths: relative

Great job! Using relative paths to the files within your workflowr project makes it easier to run your code on other machines.

Repository version: e3da525

Great! You are using Git for version control. Tracking code development and connecting the code version to the results is critical for reproducibility.

The results in this page were generated with repository version e3da525. See the Past versions tab to see a history of the changes made to the R Markdown and HTML files.

Note that you need to be careful to ensure that all relevant files for the analysis have been committed to Git prior to generating the results (you can use wflow_publish or wflow_git_commit). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:


Ignored files:
    Ignored:    .DS_Store
    Ignored:    .Rhistory
    Ignored:    .Rproj.user/
    Ignored:    analysis/.Rhistory
    Ignored:    data/.DS_Store
    Ignored:    scripts/.DS_Store

Untracked files:
    Untracked:  data/bloodcells1.csv
    Untracked:  data/height.chr3.matrix
    Untracked:  data/susie_ss_input_sex.rds

Unstaged changes:
    Modified:   analysis/compare_result_more.Rmd

Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.

These are the previous versions of the repository in which changes were made to the R Markdown (analysis/bloodcells_regions.Rmd) and HTML (docs/bloodcells_regions.html) files. If you've configured a remote Git repository (see ?wflow_git_remote), click on the hyperlinks in the table below to view the files as they were in that past version.

File	Version	Author	Date	Message
Rmd	e3da525	zouyuxin	2020-11-22	wflow_publish("analysis/bloodcells_regions.Rmd")
html	61d98d4	zouyuxin	2020-11-22	Build site.
Rmd	885829c	zouyuxin	2020-11-22	wflow_publish("analysis/bloodcells_regions.Rmd")
html	9ab9598	zouyuxin	2020-11-22	Build site.
Rmd	8b1d17e	zouyuxin	2020-11-22	wflow_publish("analysis/bloodcells_regions.Rmd")
html	2b09b42	zouyuxin	2020-11-20	Build site.
Rmd	4044ce5	zouyuxin	2020-11-20	wflow_publish("analysis/bloodcells_regions.Rmd")
html	50a6e41	zouyuxin	2020-11-20	Build site.
Rmd	3bfe683	zouyuxin	2020-11-20	wflow_publish("analysis/bloodcells_regions.Rmd")
html	717d6b1	zouyuxin	2020-11-20	Build site.
Rmd	c2021d0	zouyuxin	2020-11-20	wflow_publish("analysis/bloodcells_regions.Rmd")

There are 248,980 individuals of white British ancestries with 16 blood cells phenotypes. The script to prepare the phenotypes and covariates is get_bloodcells. The filtering steps are also described here

For genotype data, variants with imputation score (INFO) > 0.9, MAF > 1% are included in association studies.

The script to run GWAS is GWAS

For each phenotype, regions for fine-mapping are defined by greedily starting with the most significantly associated SNP, including SNPs within a window of 500kb centered at the SNP, until we include all significant SNPs (p < 5e-8). We merge ovelapping regions. We exclude HLA region (chr6: 25Mb - 36Mb). The steps are

Find the most significantly associated SNP.
Choose region +- 250kb around the SNP.
Find the next most significantly associated SNP ouside the selected regions.
Choose region +- 250kb around the SNP.
Merge regions if they overlap. ...

When we select region across traits, we include all regions from each pheotype and merge overlapping regions. This produces some very large regions with more than 10000 SNPs.

library(data.table)
library(dplyr)
library(kableExtra)
library(knitr)
pheno_names = c("WBC_count", "RBC_count", "Haemoglobin", "MCV", "RDW", "Platelet_count",
                "Plateletcrit", "PDW", "Lymphocyte_perc", "Monocyte_perc",
                "Neutrophill_perc", "Eosinophill_perc", "Basophill_perc",
                "Reticulocyte_perc", "MSCV", "HLR_perc")

trait_regions = list()
for(trait in pheno_names){
  # region = fread(paste0('/gpfs/data/stephens-lab/finemap-uk-biobank/data/raw/BloodCells/regions/', trait, '_regions'))
  region = fread(paste0('~/Desktop/ukb-bloodcells/regions_raw/', trait, '_regions'))
  region = region %>% arrange(desc(logp))
  region_r = c()
  for(i in 1:22){
    region.chr = region %>% filter(CHR == i) %>% arrange(start)
    if(nrow(region.chr) == 0){
      next
    }
    tmp = region.chr %>% group_by(g = cumsum(cummax(lag(end, default = first(end))) < start)) %>%
      summarise(start = first(start), end = max(end), logp = max(logp),.groups = 'drop') %>% 
      mutate(length = end - start) %>%
      mutate(CHR = i) %>% select(CHR, start, end, length, logp)
    region_r = rbind(region_r, tmp)
  }
  trait_regions[[trait]] = region_r
}

Summary of region length for each phenotype before selecting regions across traits:

for(trait in pheno_names){
  tb = rbind(summary(trait_regions[[trait]]$length), summary(trait_regions[[trait]]$logp))
  rownames(tb) = c('region_length', 'region_max_log10p')
  round(tb,3) %>% kbl(caption = paste0(trait, ': ', nrow(trait_regions[[trait]]), ' regions')) %>% kable_styling() %>% print
  cat("\n")
}

<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>WBC_count: 280 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 378854.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 610315.389 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 2301820.000 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.327 </td>
   <td style="text-align:right;"> 8.531e+00 </td>
   <td style="text-align:right;"> 10.807 </td>
   <td style="text-align:right;"> 15.286 </td>
   <td style="text-align:right;"> 14.601 </td>
   <td style="text-align:right;"> 213.266 </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>RBC_count: 310 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 653220.6 </td>
   <td style="text-align:right;"> 755004.000 </td>
   <td style="text-align:right;"> 2940124 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.307e+00 </td>
   <td style="text-align:right;"> 8.779e+00 </td>
   <td style="text-align:right;"> 11.363 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 18.728 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>Haemoglobin: 250 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 301703.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 611650.000 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 3251867.000 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.302 </td>
   <td style="text-align:right;"> 8.249e+00 </td>
   <td style="text-align:right;"> 10.408 </td>
   <td style="text-align:right;"> 17.123 </td>
   <td style="text-align:right;"> 16.425 </td>
   <td style="text-align:right;"> 261.853 </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>MCV: 326 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 670089.4 </td>
   <td style="text-align:right;"> 772562.750 </td>
   <td style="text-align:right;"> 2920824 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.312e+00 </td>
   <td style="text-align:right;"> 9.576e+00 </td>
   <td style="text-align:right;"> 13.151 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 21.175 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>RDW: 275 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 275738.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 650422 </td>
   <td style="text-align:right;"> 500000.00 </td>
   <td style="text-align:right;"> 4089364 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.309 </td>
   <td style="text-align:right;"> 9.108e+00 </td>
   <td style="text-align:right;"> 12.833 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 24.55 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>Platelet_count: 397 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 270883.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 678599.6 </td>
   <td style="text-align:right;"> 766174.000 </td>
   <td style="text-align:right;"> 3646000 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.324 </td>
   <td style="text-align:right;"> 9.018e+00 </td>
   <td style="text-align:right;"> 12.161 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 21.916 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>Plateletcrit: 382 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 270883.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 640938.8 </td>
   <td style="text-align:right;"> 754987.500 </td>
   <td style="text-align:right;"> 3593140 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.317 </td>
   <td style="text-align:right;"> 8.826e+00 </td>
   <td style="text-align:right;"> 11.431 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 19.293 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>PDW: 299 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 268602.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 650832 </td>
   <td style="text-align:right;"> 751519.000 </td>
   <td style="text-align:right;"> 4239970 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.318 </td>
   <td style="text-align:right;"> 8.975e+00 </td>
   <td style="text-align:right;"> 12.738 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 21.975 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>Lymphocyte_perc: 232 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 378854.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 621718.961 </td>
   <td style="text-align:right;"> 753623.00 </td>
   <td style="text-align:right;"> 2239954.000 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.325 </td>
   <td style="text-align:right;"> 9.048e+00 </td>
   <td style="text-align:right;"> 11.351 </td>
   <td style="text-align:right;"> 16.049 </td>
   <td style="text-align:right;"> 16.54 </td>
   <td style="text-align:right;"> 152.602 </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>Monocyte_perc: 262 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 310142.000 </td>
   <td style="text-align:right;"> 5.00e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 659465.2 </td>
   <td style="text-align:right;"> 751945.000 </td>
   <td style="text-align:right;"> 3864804 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.361 </td>
   <td style="text-align:right;"> 8.85e+00 </td>
   <td style="text-align:right;"> 12.112 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 21.705 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>Neutrophill_perc: 218 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 368127.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 612052.165 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 2208916.000 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.303 </td>
   <td style="text-align:right;"> 8.403e+00 </td>
   <td style="text-align:right;"> 10.917 </td>
   <td style="text-align:right;"> 15.978 </td>
   <td style="text-align:right;"> 16.034 </td>
   <td style="text-align:right;"> 193.111 </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>Eosinophill_perc: 277 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 397874.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 656940.123 </td>
   <td style="text-align:right;"> 750663.000 </td>
   <td style="text-align:right;"> 3994646.000 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.328 </td>
   <td style="text-align:right;"> 9.148e+00 </td>
   <td style="text-align:right;"> 13.378 </td>
   <td style="text-align:right;"> 20.676 </td>
   <td style="text-align:right;"> 20.883 </td>
   <td style="text-align:right;"> 212.968 </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>Basophill_perc: 85 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 561068.212 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 1748896.000 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.367e+00 </td>
   <td style="text-align:right;"> 8.567e+00 </td>
   <td style="text-align:right;"> 11.392 </td>
   <td style="text-align:right;"> 17.352 </td>
   <td style="text-align:right;"> 16.739 </td>
   <td style="text-align:right;"> 134.511 </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>Reticulocyte_perc: 237 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 303607.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 653925.3 </td>
   <td style="text-align:right;"> 753414.000 </td>
   <td style="text-align:right;"> 5547602 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.302 </td>
   <td style="text-align:right;"> 8.534e+00 </td>
   <td style="text-align:right;"> 12.249 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 20.093 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>MSCV: 285 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 272793.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 656840.8 </td>
   <td style="text-align:right;"> 754786.000 </td>
   <td style="text-align:right;"> 3254290 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.328 </td>
   <td style="text-align:right;"> 9.262e+00 </td>
   <td style="text-align:right;"> 12.434 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 21.282 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>
<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>HLR_perc: 246 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 381586.000 </td>
   <td style="text-align:right;"> 5.000e+05 </td>
   <td style="text-align:right;"> 500000.000 </td>
   <td style="text-align:right;"> 674728.8 </td>
   <td style="text-align:right;"> 768513.250 </td>
   <td style="text-align:right;"> 5973226 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.304 </td>
   <td style="text-align:right;"> 8.853e+00 </td>
   <td style="text-align:right;"> 12.103 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 19.311 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>

For HLR_perc, the maximum region is at CHR 3 from 46234573 to 52207799, which includes 9572 SNPs.

# gwas_HLR_perc = fread('/gpfs/data/stephens-lab/finemap-uk-biobank/data/raw/BloodCells/gwas/bloodcells_gwas_HLR_perc')
gwas_HLR_perc = fread('~/Desktop/ukb-bloodcells/bloodcells_gwas_HLR_perc')
colnames(gwas_HLR_perc)[1] = 'CHR'
gwas_HLR_perc$P = as.numeric(gwas_HLR_perc$P)
gwas_HLR_perc = gwas_HLR_perc %>% select(CHR, POS, T_STAT, P) %>% mutate(logp = -log10(P))
gwas_HLR_perc.sub = gwas_HLR_perc %>% filter(CHR == 3, POS >= 46234573, POS <=52207799)
plot(gwas_HLR_perc.sub$POS, gwas_HLR_perc.sub$logp, xlab='CHR 3 POS', ylab='-log10(p)', main='HLR_perc')

Version	Author	Date
50a6e41	zouyuxin	2020-11-20
717d6b1	zouyuxin	2020-11-20

For Reticulocyte_perc, the maximum region is at CHR 3 from 48155661 to 53703263, which includes 8888 SNPs.

# gwas_Reticulocyte_perc = fread('/gpfs/data/stephens-lab/finemap-uk-biobank/data/raw/BloodCells/gwas/bloodcells_gwas_Reticulocyte_perc')
gwas_Reticulocyte_perc = fread('~/Desktop/ukb-bloodcells/bloodcells_gwas_Reticulocyte_perc')
colnames(gwas_Reticulocyte_perc)[1] = 'CHR'
gwas_Reticulocyte_perc$P = as.numeric(gwas_Reticulocyte_perc$P)
gwas_Reticulocyte_perc = gwas_Reticulocyte_perc %>% select(CHR, POS, T_STAT, P) %>% mutate(logp = -log10(P))
gwas_Reticulocyte_perc.sub = gwas_Reticulocyte_perc %>% filter(CHR == 3, POS >= 48155661, POS <=53703263)
plot(gwas_Reticulocyte_perc.sub$POS, gwas_Reticulocyte_perc.sub$logp, xlab='CHR 3 POS', ylab='-log10(p)', main='Reticulocyte_perc')

Version	Author	Date
50a6e41	zouyuxin	2020-11-20
717d6b1	zouyuxin	2020-11-20

For PDW, the maximum region is at CHR 8 from 7838230 to 12078200, which includes 16605 SNPs.

# gwas_PDW = fread('/gpfs/data/stephens-lab/finemap-uk-biobank/data/raw/BloodCells/gwas/bloodcells_gwas_PDW')
gwas_PDW = fread('~/Desktop/ukb-bloodcells/bloodcells_gwas_PDW')
colnames(gwas_PDW)[1] = 'CHR'
gwas_PDW$P = as.numeric(gwas_PDW$P)
gwas_PDW = gwas_PDW %>% select(CHR, POS, T_STAT, P) %>% mutate(logp = -log10(P))
gwas_PDW.sub = gwas_PDW %>% filter(CHR == 8, POS >= 7838230, POS <=12078200)
plot(gwas_PDW.sub$POS, gwas_PDW.sub$logp, xlab='CHR 3 POS', ylab='-log10(p)', main='PDW')

Version	Author	Date
50a6e41	zouyuxin	2020-11-20
717d6b1	zouyuxin	2020-11-20

Select regions across phenotype:

tb = bind_rows(trait_regions, .id = "column_label")
res.final = c()
for(i in 1:22){
  tb.chr = tb %>% filter(CHR == i) %>% arrange(start)
  if(nrow(tb.chr) == 0){
    next
  }
  tmp = tb.chr %>% group_by(g = cumsum(cummax(lag(end, default = first(end))) < start)) %>%
    summarise(start = first(start), end = max(end), logp = max(logp), .groups = 'drop') %>% 
    mutate(length = end - start) %>%
    mutate(CHR = i) %>% select(CHR, start, end, length, logp)
  res.final = rbind(res.final, tmp)
}
snpsnum = c()
for(i in 1:nrow(res.final)){
  snpsnum = c(snpsnum, gwas_PDW %>% filter(CHR == res.final$CHR[i],
                                           POS >= res.final$start[i], 
                                           POS <= res.final$end[i]) %>% nrow )
}
res.final$snpsnum = snpsnum

Summary of regions:

tb = rbind(summary(res.final$length), summary(res.final$snpsnum),  summary(res.final$logp))
rownames(tb) = c('region_length', 'region_num_snps', 'region_max_log10p')
round(tb,3) %>% kbl(caption = paste0(nrow(res.final), ' regions')) %>% kable_styling() %>% print

<table class="table" style="margin-left: auto; margin-right: auto;">
<caption>972 regions</caption>
 <thead>
  <tr>
   <th style="text-align:left;">   </th>
   <th style="text-align:right;"> Min. </th>
   <th style="text-align:right;"> 1st Qu. </th>
   <th style="text-align:right;"> Median </th>
   <th style="text-align:right;"> Mean </th>
   <th style="text-align:right;"> 3rd Qu. </th>
   <th style="text-align:right;"> Max. </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:left;"> region_length </td>
   <td style="text-align:right;"> 307855.000 </td>
   <td style="text-align:right;"> 500000.00 </td>
   <td style="text-align:right;"> 658060.500 </td>
   <td style="text-align:right;"> 932085.792 </td>
   <td style="text-align:right;"> 1057863.250 </td>
   <td style="text-align:right;"> 8729501 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_num_snps </td>
   <td style="text-align:right;"> 21.000 </td>
   <td style="text-align:right;"> 1512.50 </td>
   <td style="text-align:right;"> 2015.000 </td>
   <td style="text-align:right;"> 2641.288 </td>
   <td style="text-align:right;"> 3146.250 </td>
   <td style="text-align:right;"> 21219 </td>
  </tr>
  <tr>
   <td style="text-align:left;"> region_max_log10p </td>
   <td style="text-align:right;"> 7.328 </td>
   <td style="text-align:right;"> 9.94 </td>
   <td style="text-align:right;"> 15.255 </td>
   <td style="text-align:right;"> Inf </td>
   <td style="text-align:right;"> 31.907 </td>
   <td style="text-align:right;"> Inf </td>
  </tr>
</tbody>
</table>

There are 972 regoins in total, 68 regions with length greater than 2Mb, 84 regions contain greater than 5000 SNPs.

The region with maximum length and maximum number of SNPs:

gwas_HLR_perc.max = gwas_HLR_perc %>% filter(CHR == 17, POS >= 39754910, POS <=48484411)
plot(gwas_HLR_perc.max$POS, gwas_HLR_perc.max$logp, xlab='CHR 17 POS', ylab='-log10(p)', main='HLR_perc')

Version	Author	Date
9ab9598	zouyuxin	2020-11-22

gwas_Reticulocyte_perc.max = gwas_Reticulocyte_perc %>% filter(CHR == 17, POS >= 39754910, POS <=48484411)
plot(gwas_Reticulocyte_perc.max$POS, gwas_Reticulocyte_perc.max$logp, xlab='CHR 17 POS', ylab='-log10(p)', main='Reticulocyte_perc')

Version	Author	Date
9ab9598	zouyuxin	2020-11-22

gwas_PDW.max = gwas_PDW %>% filter(CHR == 17, POS >= 39754910, POS <=48484411)
plot(gwas_PDW.max$POS, gwas_PDW.max$logp, xlab='CHR 17 POS', ylab='-log10(p)', main='PDW')

Version	Author	Date
9ab9598	zouyuxin	2020-11-22

# gwas_Lymphocyte_perc = fread('/gpfs/data/stephens-lab/finemap-uk-biobank/data/raw/BloodCells/gwas/bloodcells_gwas_Lymphocyte_perc')
gwas_Lymphocyte_perc = fread('~/Desktop/ukb-bloodcells/bloodcells_gwas_Lymphocyte_perc')
colnames(gwas_Lymphocyte_perc)[1] = 'CHR'
gwas_Lymphocyte_perc$P = as.numeric(gwas_Lymphocyte_perc$P)
gwas_Lymphocyte_perc = gwas_Lymphocyte_perc %>% select(CHR, POS, T_STAT, P) %>% mutate(logp = -log10(P))
gwas_Lymphocyte_perc.max = gwas_Lymphocyte_perc %>% filter(CHR == 17, POS >= 39754910, POS <=48484411)
plot(gwas_Lymphocyte_perc.max$POS, gwas_Lymphocyte_perc.max$logp, xlab='CHR 17 POS', ylab='-log10(p)', main='Lymphocyte_perc')

Version	Author	Date
9ab9598	zouyuxin	2020-11-22

sessionInfo()

R version 3.6.3 (2020-02-29)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS  10.16

Matrix products: default
BLAS:   /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] knitr_1.30        kableExtra_1.3.1  dplyr_1.0.2       data.table_1.13.2
[5] workflowr_1.6.2  

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.5        highr_0.8         pillar_1.4.6      compiler_3.6.3   
 [5] later_1.1.0.1     git2r_0.27.1      tools_3.6.3       digest_0.6.27    
 [9] viridisLite_0.3.0 evaluate_0.14     lifecycle_0.2.0   tibble_3.0.4     
[13] pkgconfig_2.0.3   rlang_0.4.8       rstudioapi_0.11   yaml_2.2.1       
[17] xfun_0.19         xml2_1.3.2        stringr_1.4.0     httr_1.4.2       
[21] generics_0.1.0    fs_1.5.0          vctrs_0.3.4       webshot_0.5.2    
[25] rprojroot_1.3-2   tidyselect_1.1.0  glue_1.4.2        R6_2.5.0         
[29] rmarkdown_2.5     purrr_0.3.4       magrittr_1.5      whisker_0.4      
[33] scales_1.1.1      backports_1.2.0   promises_1.1.1    ellipsis_0.3.1   
[37] htmltools_0.5.0   rvest_0.3.6       colorspace_1.4-1  httpuv_1.5.4     
[41] stringi_1.5.3     munsell_0.5.0     crayon_1.3.4

Blood cells regions

Yuxin Zou

11/19/2020