R-LCSS计算相似性

使用R，对轨迹相似度进行度量

G1：                                                                         SUID：

由编号经纬度三列组成      由编号，该编号对应轨迹的总长度组成

两者的编号必须保持一致。

将一段轨迹中GPS点与另一段轨迹进行匹配，需要将一段轨迹中的点与另一段轨迹中的点一一匹配，其时间复杂度为O（ ）。

其中轨迹中的任意两个点匹配公式如下：

公式（5-1）中，p_n^i表示轨迹i中的第n个点，dist(p_n^i,p_m^j )表示点p_n^i和点p_m^j之间的空间距离，即两个GPS点间的大圆距离，其中ɠ表示两点匹配的阈值。本次实验考虑到道路的宽度等因素，综合判断，设定ɠ值为20m，当ɠ值大于20m时，simPnt〖  〖(p〗_n^i  ,p〗_m^j)值为0，当ɠ小于20m时，其取值为0到1之间，作为一个相似性度量。

在本实验中， 是由 函数累加而来。其中需要解决的问题是，轨迹点 可能有多个点与之匹配，例如，点 通过 的验证与 匹配，点 也通过 的验证与 匹配，在这样的情况下，就要选出一个最佳匹配项^[7]。LCS中，显然是不能够有重复的点，所以本次ScoreLCS模型的作用有两个：一方面是找出最大相似性系数的匹配点，另一方面是能够在避免同一个一点被匹配两次，这样在LCS中就会出现重复的匹配点了。具体的解决方法下矩阵：

LCS轨迹按行挑选，第一行选择出[ , ]点为最佳匹配后，只需要将j轨迹对应的第二列相似性系数全部赋值为0，随后再继续按行的最大值挑选匹配点，这样就不会产生点的重复现象。

library("RODBC")
G1=read.csv("G1.txt", sep = " ",header=F)
SUID=read.csv("SUID.txt", sep = " ",header=F)


Dim=nrow(SUID)

ScoreLCS=matrix(0,Dim,Dim)
rownames(ScoreLCS)=c(1:Dim)
colnames(ScoreLCS)=c(1:Dim)

SimSeq=matrix(0,Dim,Dim)
rownames(SimSeq)=c(1:Dim)
colnames(SimSeq)=c(1:Dim)

D=matrix(0,Dim,Dim)
rownames(D)=c(1:Dim)
colnames(D)=c(1:Dim)

Math.PI = 3.1415926535897931

#算法
for ( i in 1:Dim)
{
   #从一天的数据中选取出两个SUID的数据
   SUID_VALUE1=SUID[i,1]
   
   point1LAT=G1[which(G1[,1]==SUID_VALUE1),2]
   point1LON=G1[which(G1[,1]==SUID_VALUE1),3]
   if(i==Dim)
   {
      break()
   }
   
   for (j in (i+1):Dim)
   {
      SUID_VALUE2=SUID[j,1]
      point2LAT=G1[which(G1[,1]==SUID_VALUE2),2]
      point2LON=G1[which(G1[,1]==SUID_VALUE2),3]
      
      p1 = length(point1LAT)
      p2 = length(point2LAT)
      
      #simPnt计算
      simPnt=matrix(0,p1,p2)
      
      for (k in 1:p1)
      {
         for (l in 1:p2)
         {
            Distance = 2 * asin(sqrt((sin((point1LAT[k] - point2LAT[l]) / 360 * Math.PI)*sin((point1LAT[k] - point2LAT[l]) / 360 * Math.PI)) 
                                     + cos(point1LAT[k] / 180 * Math.PI) * cos(point2LAT[l] / 180 * Math.PI) 
                                     * (sin((point1LON[k] - point2LON[l]) / 360 * Math.PI)*sin((point1LON[k] - point2LON[l]) / 360 * Math.PI)))) * 6378137
            if (Distance >= 20)
            {
               simPnt[k, l] = 0
            }else
            {
               simPnt[k, l] = 1 - Distance / 20
            }
         }
      }
      #ScoreLCS计算
      if(p1>p2)
      {
         LCS= rep(0,p1)
      }else
      {
         LCS= rep(0,p2)
      }
      coun=1
      for(k in 1:p1)
      {
         m=0
         pointI=0
         pointJ=0
         for (l in 1:p2)
         {
            if (simPnt[k, l] > m)
            {
               m = simPnt[k, l]
               pointI=k
               pointJ=l
            }
         }
         if (m > 0)
         {
            LCS[coun]=(pointI)
            coun=coun+1
            ScoreLCS[i,j] =ScoreLCS[i,j]+ simPnt[pointI, pointJ]
            for (l in 1: p1)
            {
               simPnt[l, pointJ] = 0
            }
         }          
      }
      #根据LenSeq函数计算SimSeq相似性度量
      SUMDistance=0
      LCS=LCS[which(LCS!=0)]
      if (length(LCS) <= 1)
      {
         SUMDistance = 0
      }else
      {
         for (l in 2:length(LCS))
         {
            
            SUMDistance =SUMDistance+ 2 * asin(sqrt((sin((point1LAT[LCS[l]] - point1LAT[LCS[l-1]])/360 * Math.PI)*sin((point1LAT[LCS[l]] - point1LAT[LCS[l-1]])/360 * Math.PI))
                                                    +cos(point1LAT[LCS[l]] / 180 * Math.PI)
                                                    *cos(point1LAT[LCS[l-1]] / 180 * Math.PI)*(sin((point1LON[LCS[l]] - point1LON[LCS[l-1]])/360 * Math.PI)*sin((point1LON[LCS[l]] - point1LON[LCS[l-1]]) /360 * Math.PI))
            ))*6378137
           
         }
         
      }
      
      if (SUMDistance == 0)
      {
         SimSeq[i,j] = 0
      }else
      {
         if (SUID[i, 1] != 0)
         {
            SimSeq[i,j] = SUMDistance / SUID[i, 2]
         }else
         {
            SimSeq[i,j] = 0
         }
         
      }
      D[i,j] = 1 - SimSeq[i,j];
   }
   print(i)
}