2

我正在尝试根据不同患者的医生就诊序列创建马尔可夫转换矩阵。在我的马尔可夫模型中,状态是不同的医生,而联系是患者的访问。患者可以留在同一个提供者处,也可以过渡到另一个提供者进行下一次访问。使用该信息我需要创建一个转换矩阵。

这是excel中的一部分数据。数据包括对近 100 家不同提供商的 3 万多次访问。

这是excel中数据的一部分。 数据

如何使用此 excel 数据(或 csv)并创建一个马尔可夫转换矩阵作为访问次数,例如:....

我需要的矩阵如下所示:

在此处输入图像描述

如何使用 R 将我的数据转换为转换矩阵?

我对 R 很陌生,真的需要帮助。

谢谢

4

2 回答 2

1

我想在不使用 data.table 的情况下比较我的方法,发现它快 45 倍(并且可能更容易理解)。

首先,我从接受的答案中对 data.table 解决方案进行计时:

rm(list=ls())
library(readxl)
library(data.table)

############## Using data.table method() ######################
data <- setDT(read_excel("Book2.xlsx"))[!is.na(PatId)]
data[ , (names(data)) := lapply(.SD, as.integer)]
provs <- data[ , sort(unique(SeenByProv))]
nprov <- length(provs)
markov <- matrix(nrow = nprov, ncol = nprov, dimnames = list(provs, provs))

system.time(      ## Timing the main loop
  for (pr in provs){
    markov[as.character(pr), ] <-
      data[ , {nxt <- SeenByProv[which(SeenByProv == pr) + 1L]
      .(prov = provs, count =
          sapply(provs, function(pr2) sum(nxt == pr2, na.rm = TRUE)))}, by = PatId
      ][, sum(count), by = prov]$V1
  }
)
#   user  system elapsed 
#  3.128   0.000   3.135 
table(markov)
#markov
#   0    1    2    3    4    5    6    7    8    9   10   11   13   22  140 
#3003  308   89   34   14   11    6    4    1    3    4    1    1    1    1 

接下来仅使用基本 R 调用:

############## Using all base R calls method() ###################
tm_matrix<-matrix(0, nrow = nprov, ncol = nprov, dimnames = list(provs, provs))
d<-read_excel("Book2.xlsx")
d<-d[!is.na(d$PatId),] # Note: Data is already ordered by PatId, DaysOfStudy

baseR<-function(tm_matrix){
  d1<-cbind(d[-nrow(d),-3],d[-1,-3]); # Form the transitions and drop the DaysofStudy
  colnames(d1)<-c("SeenByProv","PatId","NextProv","PatId2");
  d1<-d1[d1$PatId==d1$PatId2,];       # Drop those transition between different patients
  d1$SeenByProv<-as.character(d1$SeenByProv); # transform to strings to use as rownames
  d1$NextProv  <-as.character(d1$NextProv);   # and column names
  for (i in 1:nrow(d1)){                      # Fill in the transition matrix
    tm_matrix[d1$SeenByProv[i],d1$NextProv[i]]<-tm_matrix[d1$SeenByProv[i],d1$NextProv[i]]+1
  };
  return(tm_matrix)
}
system.time(tm_matrix<-baseR(tm_matrix))
#   user  system elapsed 
#  0.072   0.000   0.072 

table(tm_matrix)
#tm_matrix
#   0    1    2    3    4    5    6    7    8    9   10   11   13   22  140 
#3003  308   89   34   14   11    6    4    1    3    4    1    1    1    1 

all.equal(markov,tm_matrix)
#[1] TRUE

我的 base-R 方法快 3.135/0.072 = 43.54

于 2016-01-06T00:00:31.550 回答
1

这是一种适用于您的示例数据的方法。

我将用于readxl获取数据并data.table对其进行操作。

读取数据:

library(readxl)
library(data.table)

data <- setDT(read_excel("~/Desktop/Book2.xlsx"))[!is.na(PatId)]

#read_excel doesn't have the option to specify integers... silly...
data[ , (names(data)) := lapply(.SD, as.integer)]

预分配转移矩阵:

provs <- data[ , sort(unique(SeenByProv))]
nprov <- length(provs)

markov <- matrix(nrow = nprov, ncol = nprov,
                 dimnames = list(provs, provs))

逐行分配

for (pr in provs){
  markov[as.character(pr), ] <-
    data[ , {nxt <- SeenByProv[which(SeenByProv == pr) + 1L]
    .(prov = provs, count = 
        sapply(provs, function(pr2) sum(nxt == pr2, na.rm = TRUE)))}, by = PatId
    ][, sum(count), by = prov]$V1
}

这可能会在一些地方加速,但它确实有效。

于 2015-12-31T21:46:19.087 回答