Вычисление частоты генов определенных c областей домена с использованием R - PullRequest
Новая
       138

Вычисление частоты генов определенных c областей домена с использованием R

0 голосов
/ 04 мая 2020

Я изучаю топологию генома у растений, чтобы понять организацию / комплементацию ядерного генома.

Я хотел бы узнать быстрые методы для расчета плотности генов по конкретным c интересующим областям генома, чтобы построить график частот генов на комп на каждый бин +/- 50Kb от начала и конца области ,

У кого-нибудь есть идеи, чтобы решить эту проблему? Заранее спасибо

1 Ответ

0 голосов
/ 04 мая 2020

Если вы хотите, чтобы у меня была публикация в Biostar с дополнительной информацией о том, что я делаю:

публикация в BioStar с цифрами, примерными для желаемого результата

Так что для Плотность генов, я имею в виду количество генов, которые перекрывают мои регионы. Но я нормализовал это число с моим размером региона, потому что они имеют разный размер.

здесь у вас есть мой объект dataframe региона: 

> head(regionFileCustom)
            chr startCust  endCust                     id score strand
1 CMiso1.1chr01    -50000  3550000  Domain-1CMiso1.1chr01  1000      +
2 CMiso1.1chr01  14450000 15300000 Domain-17CMiso1.1chr01  1000      +
3 CMiso1.1chr01  15250000 15500000 Domain-18CMiso1.1chr01  1000      +
4 CMiso1.1chr01  15600000 16650000 Domain-19CMiso1.1chr01  1000      +
5 CMiso1.1chr01  16600000 17050000 Domain-20CMiso1.1chr01  1000      +
6 CMiso1.1chr01  26950000 27450000 Domain-34CMiso1.1chr01  1000      +

затем мой объект гена dataframe: 

> head(features)
            chr start   end                    id score strand
1 CMiso1.1chr01   955 24594 CMiso1.1chr01g0058061     .      +
2 CMiso1.1chr01 24595 25591 CMiso1.1chr01g0058071     .      -
3 CMiso1.1chr01 25797 30580 CMiso1.1chr01g0058081     .      +
4 CMiso1.1chr01 31006 35465 CMiso1.1chr01g0058091     .      -
5 CMiso1.1chr01 41322 41398 CMiso1.1chr01g0058101     .      -
6 CMiso1.1chr01 42694 46218 CMiso1.1chr01g0058111     .      -

Итак, для подсчета я попробую эту собственную функцию сборки: 

# function to compute the gene density per bin and the gene frequency (gene per bin per nbr of bin per genome region)

width = 201L
step=201L
n = 100
minoverlap=21L
maxgap=-1L
ignore.strand=TRUE # I don't know yet if I have to take in count the strand for this part of the analysis.

# For test windowBin.gr = windowBinVarName

computeFeatureFreq <- function(features, regionFile, outDFName){
  # convert feature as GRange
  features.gr = makeGRangesFromDataFrame(features)
  # convert regionFile as Grange
  regionFile.gr = makeGRangesFromDataFrame(regionFile)
  # Counting gene density per genome region
  # feature.density.windowBin = countOverlaps(windowBin.gr,features.gr,ignore.strand=ignore.strand, minoverlap=minoverlap, maxgap=maxgap)
  # feature.density.windowBin = countOverlaps(features.gr, windowBin.gr,ignore.strand=ignore.strand, minoverlap=minoverlap, maxgap=maxgap)
  feature.density.region = countOverlaps(features.gr, regionFile.gr)
  # retreive coordonates for each windowBin
  regionFile.df=as.data.frame(regionFile.gr)
  feature.density.region.df = data.frame(chr=regionFile.df$seqnames, start=regionFile.df$start, end=regionFile.df$end, count=feature.density.region)
  return(feature.density.region.df)
}

Если что-то не понятно, скажите мне, спасибо

...