Как загрузить счет из существующих матриц в пакете diffHi c - PullRequest
Новая
       55

Как загрузить счет из существующих матриц в пакете diffHi c

0 голосов
/ 28 апреля 2020

У меня есть необработанные матрицы взаимодействия Hi C (с использованием фрагментов рестрикционных ферментов), рассчитанные с использованием инструмента HiCExplorer. Я пытаюсь использовать их для обнаружения дифференциальных контактов с помощью пакета diffHic R. У меня есть матрицы в формате h5, но экспортированные матрицы TXT (через более прохладную версию). Т.е. у меня есть следующее

  1. Матрица фрагмента кровати (char, start, end): 
chr2L   0   5663
chr2L   5663    6452
chr2L   6452    7307
Матрица подсчета (idx1, idx2, count; где idx1 и idx2 относятся к индексу фрагмента в матрице фрагментов): 
0   0   1
0   1   1
0   2   6

Я намеревался последовать совету diffHi c виньетка (https://bioconductor.org/packages/release/bioc/vignettes/diffHic/inst/doc/diffHicUsersGuide.pdf раздел 3.6) и загрузка моих данных в несколько ContactMatrix (по 1 на образец) из InteractionSet, но я застрял при создании такого ContactMatrix.

Что я сделал: 

> library(rtracklayer)
> library(tidyverse)
# load fragments (bed like matrix)
> restrSites = rtracklayer::import( file("RestrSites_bins.txt"), format = "bed")

# read in the count matrix (second matrix), increase index number by 1 since R is 1-based
> M <- read_delim("mysecondmatrix.txt", delim = "\t", col_names = c("idx1", "idx2", "cnt") , col_types = > cols( idx1 = col_integer(), idx2 = col_integer(), cnt = col_integer() ) ) %>% 
>   dplyr::mutate( start=as.integer(idx1+1), end=as.integer(idx2+1) ) %>% 
>   dplyr::select(start, end, cnt) 

# create a GInteractions
> gi <- GInteractions(M$start, M$end, restrSites)

# convert to ContactMatrix --version1 : use boolean to exclude all fragments from "weird" chr
> keep = as.logical( seqnames(regions(gi)) %in% c("chr2L","chr2R","chr3L","chr3R", "chrX", "chr4") )
> cm <- InteractionSet::inflate(gi, keep, keep, fill=M$cnt[keep])

Error in out.mat[(ac2[relevantA] - 1L) * nR + ar1[relevantA]] <- fill[relevantA] : 
  NAs are not allowed in subscripted assignments
In addition: Warning message:
In (ac2[relevantA] - 1L) * nR : NAs produced by integer overflow

# convert to ContactMatrix --version2: keep all
> cm <- InteractionSet::inflate(gi, NULL, NULL, fill=M$cnt[keep])

Error in out.mat[(ac2[relevantA] - 1L) * nR + ar1[relevantA]] <- fill[relevantA] : 
  NAs are not allowed in subscripted assignments
In addition: Warning message:
In (ac2[relevantA] - 1L) * nR : NAs produced by integer overflow

, если бы кто-нибудь мог указать мне на ошибку или предложить другой способ загрузки данных в пакет diffHi c, это было бы здорово (я могу также конвертируйте h5 в формат hi c, bedpe, ...).

...