Как вычислить различия между последовательностями, когда последовательности содержат пробелы?
Я хочу кластеризовать последовательности с оптимальным соответствием TraMineR::seqdist() из данных, которые содержат пропуски, то есть последовательности, содержащие пропуски.
library(TraMineR)
data(ex1)
sum(is.na(ex1))
# [1] 38
sq <- seqdef(ex1[1:13])
sq
# Sequence
# s1 *-*-*-A-A-A-A-A-A-A-A-A-A
# s2 D-D-D-B-B-B-B-B-B-B
# s3 *-D-D-D-D-D-D-D-D-D-D
# s4 A-A-*-*-B-B-B-B-D-D
# s5 A-*-A-A-A-A-*-A-A-A
# s6 *-*-*-C-C-C-C-C-C-C
# s7 *-*-*-*-*-*-*-*-*-*-*-*-*
sm <- seqsubm(sq, method='TRATE')
round(sm,digits=3)
# A-> B-> C-> D->
# A-> 0 2.000 2 2.000
# B-> 2 0.000 2 1.823
# C-> 2 2.000 0 2.000
# D-> 2 1.823 2 0.000
Когда я бегу seqdist()
dist.om <- seqdist(sq, method="OM", indel=1, sm=sm)
Я получаю
Error: 'with.missing' must be TRUE when 'seqdata' or 'refseq' contains missing values
но когда я установил опцию with.missing=TRUE Я получаю
[>] including missing values as an additional state
[>] 7 sequences with 5 distinct states
[>] checking 'sm' (one value for each state, triangle inequality)
Error: [!] size of substitution cost matrix must be 5x5
Итак, как мы можем вычислить различия между последовательностями, используя seqdist() с выходом seqsubm() правильный путь, когда данные содержат пропуски, то есть последовательности содержат пробелы?
Примечание: я не очень уверен, имеет ли это смысл. Пока я просто исключаю наблюдения с пропусками, но из-за своих данных я теряю много наблюдений по этому поводу. Поэтому было бы полезно узнать, есть ли такая возможность.
1 ответ
Существуют разные стратегии для вычисления расстояний, когда у вас есть пробелы.
1) Первое решение состоит в том, чтобы рассматривать отсутствующее состояние как дополнительное состояние. Это то, что seqdist делает, когда вы установите with.missing=TRUE, В этом случае sm Матрица должна содержать затраты на замену любого состояния отсутствующим состоянием. С помощью seqsubm вам просто нужно указать with.missing=TRUE для этой функции тоже. По умолчанию затраты на замену "отсутствующими" устанавливаются как фиксированные значения. miss.cost (2 по умолчанию).
sm <- seqsubm(sq, method='TRATE', with.missing=TRUE)
round(sm,digits=3)
# A-> B-> C-> D-> *->
# A-> 0 2.000 2 2.000 2
# B-> 2 0.000 2 1.823 2
# C-> 2 2.000 0 2.000 2
# D-> 2 1.823 2 0.000 2
# *-> 2 2.000 2 2.000 0
Чтобы получить стоимость замещения для "пропавших без вести" на основе вероятностей перехода
sm <- seqsubm(sq, method='TRATE', with.missing=TRUE, miss.cost.fixed=FALSE)
round(sm,digits=3)
# A-> B-> C-> D-> *->
# A-> 0.000 2.000 2.000 2.000 1.703
# B-> 2.000 0.000 2.000 1.823 1.957
# C-> 2.000 2.000 0.000 2.000 1.957
# D-> 2.000 1.823 2.000 0.000 1.957
# *-> 1.703 1.957 1.957 1.957 0.000
Используя последний sm, мы получаем расстояния между последовательностями
dist.om <- seqdist(sq, method="OM", indel=1, sm=sm, with.missing=TRUE)
round(dist.om, digits=2)
# [,1] [,2] [,3] [,4] [,5] [,6] [,7]
# [1,] 0.00 22.87 21.91 18.41 6.41 17.00 17.03
# [2,] 22.87 0.00 13.76 11.56 19.91 19.87 22.57
# [3,] 21.91 13.76 0.00 14.25 18.96 18.91 21.57
# [4,] 18.41 11.56 14.25 0.00 13.70 15.70 18.14
# [5,] 6.41 19.91 18.96 13.70 0.00 15.70 16.62
# [6,] 17.00 19.87 18.91 15.70 15.70 0.00 16.70
# [7,] 17.03 22.57 21.57 18.14 16.62 16.70 0.00
Конечно, последовательности будут тогда близки друг к другу только потому, что они имеют много пропущенных состояний (*). Поэтому вы можете захотеть сохранить только те последовательности, в которых, например, пропущено менее 10 % элементов.
2) Второе решение состоит в том, чтобы удалить пробелы, которые вы делаете в seqdef, (Однако обратите внимание, что это меняет выравнивание.)
## Here, we drop seq 7 that contains only missing values
sq <- seqdef(ex1[-7,1:13], left='DEL', gaps='DEL')
sq
# Sequence
# s1 A-A-A-A-A-A-A-A-A-A
# s2 D-D-D-B-B-B-B-B-B-B
# s3 D-D-D-D-D-D-D-D-D-D
# s4 A-A-B-B-B-B-D-D
# s5 A-A-A-A-A-A-A-A
# s6 C-C-C-C-C-C-C
sm <- seqsubm(sq, method='TRATE')
round(sm,digits=3)
# A-> B-> C-> D->
# A-> 0.000 1.944 2 2.000
# B-> 1.944 0.000 2 1.823
# C-> 2.000 2.000 0 2.000
# D-> 2.000 1.823 2 0.000
dist.om <- seqdist(sq, method="OM", indel=1, sm=sm)
round(dist.om, digits=2)
# [,1] [,2] [,3] [,4] [,5] [,6]
# [1,] 0.00 19.61 20.00 13.78 2.00 17
# [2,] 19.61 0.00 12.76 9.59 17.61 17
# [3,] 20.00 12.76 0.00 13.29 18.00 17
# [4,] 13.78 9.59 13.29 0.00 11.78 15
# [5,] 2.00 17.61 18.00 11.78 0.00 15
# [6,] 17.00 17.00 17.00 15.00 15.00 0