Как рассчитать расстояния между генами, используя несколько хромосом

Question

Я хочу выяснить, объединяются ли определенные гены вместе. Теперь у меня уже есть список генов, а также места их старта и остановки, и я уже знаю, как рассчитать расстояние между этими генами. проблема в том, что я не знаю, как принять во внимание переход в хромосомах.

Невозможно измерить расстояние между геном в хромосоме 1 и геном в хромосоме 2.

Я думал о расчете расстояния следующим образом: начальное местоположение гена 2 - конечное местоположение гена 1. Затем у вас есть расстояние между этими генами.

Но как мне объяснить это: когда вы достигаете следующей хромосомы, код R захватывает начальное местоположение гена на хромосоме 2, но остановочное местоположение гена на хромосоме 1, и это невозможно (для мое исследование, по крайней мере).

Так что мне интересно, как объяснить это в R. Мне просто нужно каким-то образом пропустить гены, если они находятся на разных хромосомах.

Надеюсь, вы, ребята, сможете мне помочь.

о коде ниже: три вектора - это просто векторы начального и конечного местоположения и хромосомы. они все равной длины. хромосомы - это вектор, содержащий номер хромосомы для каждого гена

start_vector <- as.vector(sorted_coords$start_position)
end_vector <- as.vector(sorted_coords$end_position)
chromosomes <- as.vector(sorted_coords$chromosome_name)

chromosomes[is.na(chromosomes)] <- 24

count = 0
for(i in 1:length(chromosomes)){
  if(count != chromosomes[i]){
    start <- i - 1
    end <- i + 1

    start_vector <- start_vector[-start]
    end_vector <- end_vector[-end]
    count <- count + 1  
    }

}

Я ожидаю вектор расстояний всех генов, исключая расстояния генов, которые лежат на разных хромосомах.

Answer 1

library(tidyverse) # for all the tidyverse goodies
library(reshape2) # For the melt function

Поскольку вы не предоставили воспроизводимый пример, я позволил себе создать свой собственный игрушечный фрейм данных, как показано ниже. Он имеет только 2 хромосомы, но этот метод должен быть применим к произвольному количеству хромосом и генов.

sorted_coords <- tibble(start_position = abs(rnorm(10)*10),
                        end_position = abs(rnorm(10)*10),
                        chromosome_name = c(rep(1,5),rep(2,5)))

Редактировать: ОП уточнил, что они хотели найти расстояние только до гена, а не до любого другого гена. Метод для выполнения последней части находится внизу, потому что я нашел это интересным. Новое решение здесь:

sorted_coords %>% 
  group_by(chromosome_name) %>%  
  arrange(chromosome_name, start_position) %>% 
  mutate(distance = start_position - lag(end_position, n = 1, default = 0))

1. Мы группируем по хромосомам, чтобы не делать ошибочных расчетов между хромосомами.

2. В конце сортируем по названию хромосомы для сортировки. Мы организуем по стартовой позиции, чтобы гены были в правильном порядке.

3. Мы рассчитываем расстояние как предложено. Начальная позиция текущей строки - конечная позиция предыдущей строки. Мы указываем (хотя это по умолчанию), что мы смотрим на строку раньше, и что если нет строки до того, как значение конечной позиции по умолчанию равно 0.

---

Старый ответ

Если вы хотите сравнить каждый ген с каждым другим геном, самый быстрый способ сделать это - создать матрицу. Как вы указали, мы хотим вычесть начало гена 1 до конца гена 2. Мне это не кажется правильным, но прошло уже много времени с тех пор, как я сделал биохимию :). Поскольку вам нужен один список пар, мы можем свернуть его (функция плавления).

Код ниже немного неясен для понимания, поэтому давайте разберем его.

sorted_coords %>% 
  group_by(chromosome_name) %>% 
  do( outer(.$start_position, .$end_position) %>% 
        melt() %>% 
        setNames(c("rows", "columns", "distance")))

1. Мы берем наш фрейм данных и группируем его по каждой хромосоме, как вам нужно.
2. Команда do позволяет нам выполнять сложные операции. Команда group_by гарантирует, что все, что мы делаем, по существу является изолированным для каждой хромосомы.
3. Внешняя функция создает нашу матрицу для нас. . - это передаваемый нами фрейм данных (подмножество определенной хромосомы). Мы передаем два столбца, для которых нам нужно найти разницу.
4. Функция melt превращает матрицу в фрейм данных для нас, указывая, какие два гена она использовала для вычисления разницы. Они численно перечислены, и вы можете вернуться и сравнить его. Я бы предложил использовать порядок, чтобы установить порядок, чтобы вы могли легко ссылаться на него.
5. setNames просто устанавливает имена столбцов на что-то более разрешимое.

Это должно быть намного быстрее, чем запуск цикла for для всех ваших процессов. Если вы предоставите больше информации, я, скорее всего, уточню ответ.