利用互信息量法确定延迟时间T
###确定延迟时间tau——互信息量法
MI<-function(data, tau_max, n){
  I_sq = rep(NA,30)  # 建立一个空数组，用来存放每个τ值的互信息
  data_length = length(data) # 计算时间序列的长度
  for (tau in 31:tau_max){
    S = data[1:(data_length-tau)]
    Q = data[(tau+1):data_length]
    as1 = min(S)
    bq = min(Q)
    delts = (max(S)-as1)/n
    deltq = (max(Q)-bq)/n
    N_sq = matrix(0,ncol=n,nrow = n)
    for (i in 1:n){
      for (j in 1:n){
        for (k in 1:(data_length-tau)){
          as_k = (S[k]-as1)/delts
          bq_k = (Q[k]-bq)/deltq
          if (as_k >= i-1 && as_k < i && bq_k >=j-1 && bq_k < j){
            N_sq[i, j] = N_sq[i, j]+1
          }
        }
      }
    }
    Ntotal = sum(N_sq)
    Ps<-rep(0,n)
    Pq<-rep(0,n)
    for (i in 1:n) {
      Ps[i] = sum(N_sq[,i])/Ntotal   # 计算位于一维s格子内的概率
      Pq[i] = sum(N_sq[i,])/Ntotal  # 计算位于一维q格子内的概率
    }
    Psq = N_sq/Ntotal         # 计算位于二维格子(ii,jj)内概率
    H_s = 0  # 计算S的熵
    H_q = 0  # 计算q的熵
    for (i in 1:n){
      if (Ps[i] != 0)
      {H_s = H_s - Ps[i] * log(Ps[i])}
      if (Pq[i] != 0)
      {H_q = H_q - Pq[i] * log(Pq[i])}
    }
   
    H_sq = 0  # 计算(s,q)的联合熵
    for (i in 1:n){
      for (j in 1:n){
        if (Psq[i, j] != 0)
        {H_sq = H_sq - Psq[i, j] * log(Psq[i, j])}
      }
    }
   
    I_sq[tau-30] = H_s+H_q-H_sq  # 计算tau下的互信息函数
  }
  return(I_sq)
  }
选择互信息函数的第一个极小值点对应的T作为最优延迟时间

R里面nonlinearTseries包里面的timeLag函数可以用于选择延迟时间tau
参数technique = "acf"时，表示用自相关函数法，此时参数selection.method选择默认
参数technique = "ami",表示用互信息量法，此时参数selection.method = "first.minimum"

嵌入维数m可用函数estimateEmbeddingDim()函数，此时使用Cao方法确定维数的
estimateEmbeddingDim(data, time.lag = tau,do.plot = F,max.embedding.dim = 20)
tseriesChaos包中的false.nearest()函数使用伪近邻方法确定嵌入维数
false.nearest(data, m = max.m, d = tau, t = 150, eps = sd(data)/10)