Я пишу симуляцию, состоящую из сетки n x n ячеек. В разное время моделирования ячейка случайным образом выбирается для «разделения». Когда клетка делится, она умирает и создает две дочерние клетки. Одна дочь заменяет исходную ячейку, а другая дочь случайным образом заменяет одного из восьми соседей по сетке.
Сетка кодируется фреймом данных с n ^ 2 строками в начале, по одной строке для каждой ячейки (каждая ячейка имеет рождение_время = 0, смерть_время = 50 и родитель = 0 в начале). По ходу моделирования для каждого события деления добавляются две строки, представляющие дочерние клетки, и обновляются времена смерти родителя (и соседнего предшественника). Дочерям назначается время рождения!=0, время_смерти=50 и родитель (примеры см. ниже).
После того, как симуляция выполнялась в течение определенного периода времени (50 в приведенных ниже примерах), я беру выборку ячеек, имеющих одинаковую координату x. Для этих ячеек я хотел бы использовать историческую информацию, закодированную в моем фрейме данных сетки, чтобы найти время их слияния, то есть время смерти всех ячеек, которые являются предками двух или более ячеек в окончательной выборке. Я ищу функцию для выполнения этого в R (или помощь в построении алгоритма, который я мог бы кодировать в R самостоятельно).
Ниже приведены три примера, которые, я надеюсь, прояснят мои требования:
Тест1:
> grid1
cellID x_coordinate y_coordinate onEdge parent birth_time death_time
1 1 1 1 1 0 0 50
2 2 2 1 1 0 0 50
3 3 3 1 1 0 0 2
4 4 4 1 1 0 0 50
5 5 5 1 1 0 0 50
6 6 1 2 1 0 0 50
7 7 2 2 0 0 0 50
8 8 3 2 0 0 0 2
9 9 4 2 0 0 0 50
10 10 5 2 1 0 0 50
11 11 1 3 1 0 0 50
12 12 2 3 0 0 0 50
13 13 3 3 0 0 0 12
14 14 4 3 0 0 0 50
15 15 5 3 1 0 0 50
16 16 1 4 1 0 0 50
17 17 2 4 0 0 0 50
18 18 3 4 0 0 0 21
19 19 4 4 0 0 0 50
20 20 5 4 1 0 0 50
21 21 1 5 1 0 0 50
22 22 2 5 1 0 0 50
23 23 3 5 1 0 0 50
24 24 4 5 1 0 0 50
25 25 5 5 1 0 0 50
26 26 3 2 0 8 2 12
27 27 3 1 1 8 2 50
28 28 3 2 0 26 12 33
29 29 3 3 0 26 12 21
30 30 3 3 0 29 21 33
31 31 3 4 0 29 21 45
32 32 3 3 0 30 33 45
33 33 3 2 0 30 33 50
34 34 3 4 0 31 45 50
35 35 3 3 0 31 45 50
Я выбираю склепы, которые существуют в конечное время (50) и имеют x-координату = 3. Обратите внимание, что хотя в этом тестовом примере я выбираю все 5 крипт, в реальной симуляции будет выбрано подмножество.
> sample1
cellID x_coordinate y_coordinate onEdge parent birth_time death_time
23 23 3 5 1 0 0 50
27 27 3 1 1 8 2 50
33 33 3 2 0 30 33 50
34 34 3 4 0 31 45 50
35 35 3 3 0 31 45 50
В этом примере ячейка в (3,5) не связана с другими (за исключением псевдородительского узла всех ячеек (0). Все остальные четыре ячейки связаны, и есть 3 события деления, которые являются информативными для филогении.
> res1
cellID x_coordinate y_coordinate onEdge parent birth_time death_time
1 8 3 2 0 0 0 2
3 29 3 3 0 26 12 21
5 31 3 4 0 29 21 45
В приведенном ниже дереве показаны отношения, которые я пытаюсь зафиксировать.
Вот два других примера: Тест2:
> grid2
cellID x_coordinate y_coordinate onEdge parent birth_time death_time
1 1 1 1 1 0 0 50
2 2 2 1 1 0 0 2
3 3 3 1 1 0 0 50
4 4 4 1 1 0 0 45
5 5 5 1 1 0 0 50
6 6 1 2 1 0 0 50
7 7 2 2 0 0 0 2
8 8 3 2 0 0 0 45
9 9 4 2 0 0 0 21
10 10 5 2 1 0 0 21
11 11 1 3 1 0 0 50
12 12 2 3 0 0 0 50
13 13 3 3 0 0 0 33
14 14 4 3 0 0 0 50
15 15 5 3 1 0 0 50
16 16 1 4 1 0 0 50
17 17 2 4 0 0 0 33
18 18 3 4 0 0 0 12
19 19 4 4 0 0 0 50
20 20 5 4 1 0 0 50
21 21 1 5 1 0 0 50
22 22 2 5 1 0 0 50
23 23 3 5 1 0 0 50
24 24 4 5 1 0 0 12
25 25 5 5 1 0 0 50
26 26 2 2 0 7 2 50
27 27 2 1 1 7 2 50
28 28 3 4 0 18 12 50
29 29 4 5 1 18 12 50
30 30 4 2 0 9 21 50
31 31 5 2 1 9 21 50
32 32 2 4 0 17 33 50
33 33 3 3 0 17 33 50
34 34 3 2 0 8 45 50
35 35 4 1 1 8 45 50
> sample2
cellID x_coordinate y_coordinate onEdge parent birth_time death_time
3 3 3 1 1 0 0 50
23 23 3 5 1 0 0 50
28 28 3 4 0 18 12 50
33 33 3 3 0 17 33 50
34 34 3 2 0 8 45 50
Ячейки в примере2 совершенно не связаны между собой (их самым последним общим предком является псевдоузел 0). Функция не должна ничего возвращать (или только время 0).
Тест3:
> grid3
cellID x_coordinate y_coordinate onEdge parent birth_time death_time
1 1 1 1 1 0 0 50
2 2 2 1 1 0 0 50
3 3 3 1 1 0 0 50
4 4 4 1 1 0 0 50
5 5 5 1 1 0 0 50
6 6 1 2 1 0 0 50
7 7 2 2 0 0 0 31
8 8 3 2 0 0 0 34
9 9 4 2 0 0 0 37
10 10 5 2 1 0 0 50
11 11 1 3 1 0 0 50
12 12 2 3 0 0 0 22
13 13 3 3 0 0 0 8
14 14 4 3 0 0 0 8
15 15 5 3 1 0 0 6
16 16 1 4 1 0 0 50
17 17 2 4 0 0 0 2
18 18 3 4 0 0 0 2
19 19 4 4 0 0 0 3
20 20 5 4 1 0 0 50
21 21 1 5 1 0 0 50
22 22 2 5 1 0 0 50
23 23 3 5 1 0 0 50
24 24 4 5 1 0 0 50
25 25 5 5 1 0 0 50
26 26 2 4 0 17 2 50
27 27 3 4 0 17 2 3
28 28 3 4 0 27 3 45
29 29 4 4 0 27 3 6
30 30 4 4 0 29 6 50
31 31 5 3 1 29 6 50
32 32 4 3 0 14 8 50
33 33 3 3 0 14 8 22
34 34 3 3 0 33 22 45
35 35 2 3 0 33 22 31
36 36 2 3 0 35 31 50
37 37 2 2 0 35 31 34
38 38 2 2 0 37 34 50
39 39 3 2 0 37 34 37
40 40 3 2 0 39 37 49
41 41 4 2 0 39 37 50
42 42 3 3 0 34 45 49
43 43 3 4 0 34 45 50
44 44 3 3 0 42 49 50
45 45 3 2 0 42 49 50
> sample3 <- subset(grid3, x_coordinate==3 & death_time==50)
> sample3
cellID x_coordinate y_coordinate onEdge parent birth_time death_time
3 3 3 1 1 0 0 50
23 23 3 5 1 0 0 50
43 43 3 4 0 34 45 50
44 44 3 3 0 42 49 50
45 45 3 2 0 42 49 50
В этой сетке много событий, перекрывающих координату x 3, но только два из них являются информативными:
> res3
cellID x_coordinate y_coordinate onEdge parent birth_time death_time
1 42 3 3 0 34 45 49
2 34 3 3 0 33 22 45
Если кому-то это покажется полезным, вот мой полугрубый рисунок состояния каждой сетки в каждый момент времени (игнорируйте две верхние строки):
Спасибо большое за помощь!
Ваша проблема сложна для понимания и я не до конца понимаю, что вам нужно и почему для результата выбрана каждая строка данных. Моя функция проверяет предков в каждом поколении тех, кто выживает по соседству, и возвращает их информацию. Возможно, это послужит ориентиром в решении вашей проблемы.
find.elders = function(x, dead, dat){
locals = dat[dat$x_coordinate == x & dat$death_time != dead,]
survivors = dat[dat$x_coordinate == x & dat$death_time == dead,]
anc = survivors$parent
res = NULL
while(any(anc != 0)){
anc = anc[anc > 0]
cat("Ancestors:", anc, "\n")
res = c(res, which(locals$parent %in% anc))
survivors = locals[locals$cellID %in% anc,]
anc = survivors$parent
}
#res = c(res, which(locals$parent %in% anc))
locals[res,]
}
find.elders(3, 50, grid1)
Для будущих читателей я понимаю, что вопрос довольно сложен и труден для понимания. Прошу прощения за это.
В итоге я решил свою проблему, добавив в кадры данных сетки атрибут «pathString» в форме 0/1, 0/1/27 и т. д., где для каждой ячейки ячейка хранит всех своих предков, а также «саму себя».
Затем я мог бы использовать функцию as.Node()
в пакете data.tree
в R, чтобы преобразовать мою сетку в объект дерева, который впоследствии можно преобразовать в объект phylo с помощью функции as.phylo()
в ape.
После того, как выбранные ячейки сохранены в виде дерева, существующие функции в ape
и ggtree
упрощают все остальное.
См. data(acme)
пакета data.tree
и пример #Tree здесь: https://rdrr.io/cran/data.tree/man/as.Node.data.frame.html
Если клеточную родословную можно переформатировать в филогению (класс
ape::phylo), function
getMRCA` вернет предков иnode.height
возраст предка. Однако у меня проблемы с чтением ваших данных как филогении.