我有土壤水分亏缺（SMD）表与170列（每列是一个月）和103937行，我想计算p和q如下面的等式，我写了一个代码，但它不能从第四行;它说：

向量中的错误（长度= 21，na.rm = TRUE是未使用的参数）

数据中有很多NA，我不想包括：方程式

p=1-(m/m+b)....1
    q=C/(m+b)......2

其中m和b是最干燥和最潮湿条件下累积SMD与1个月至18个月不同持续时间之间线性回归的斜率和截距（图9） . 对于每个网格单元，要在干燥条件下评估m和b，首先选择历史中最干燥的月份（最低SMD）并绘制一个月的持续时间 . 然后计算每两个相邻月份的SMD运行总和，并选择最低累积SMD，持续两个月 . 重复相同的过程直到持续18个月并且选择最高累积SMD用于潮湿条件 . 然后使用线性回归拟合这些图并识别斜率，-m和截距b . C来自干旱专论的最佳拟合线，其范围从-100到100，然后按比例缩放以适合PDSI类别（-4,4）的范围 . 代码如下：

SM=read.table('SMD.csv',header=T,sep=',')
df=data.frame(data[3:21])#subset from 3 to 21 column; i have 2000 column and 103937rows.
matrix=data.matrix(df)
x=t(t(c(matrix[3:21])))
dry=vector(Length=21, na.rm=TRUE)
wet=vector(length= 21,na.rm=TRUE )
 slope_dry= vector(length= 103937,na.rm=TRUE )
 slope_wet= vector(length= 103937,na.rm=TRUE )
 inter_dry= vector(length= 103937,na.rm=TRUE )
 inter_wet= vector(length= 103937,na.rm=TRUE )

 for (a in 1:103937){
    for (i in 1:103937) {
       sum_SMD=vector(length=nrow(matrix)-i+1)
       for (j in 1 : (nrow(matrix)-i+1)) {
           for(b in j :(j+i-1))
       sum_SMD[j]<-sum_SMD[j]+SMD[b,a]
       }
       dry[i]<-min(sum_SMD)
       wet[i]<-max(sum_SMD)

       }
       model_dry<-lm(dry~x)
          slope_dry[a]<-coefficients(model_dry)[2]
          inter_dry[a]<-coefficients(model_dry)[1]
       model_wet<-lm(wet~x)
          slope_wet[a]<-coefficients(model_wet)[2]
          inter_wet[a]<-coefficients(model_wet)[1]
       }

  c_dry=slope_dry/25
   #c_dry=-4
   p_dry=1-slope_dry/(slope_dry+inter_dry)
   q_dry=c_dry/(slope_dry+inter_dry)

   #c_wet=4
  c_wet=slope_wet/25
   p_wet=1-slope_wet/(slope_wet+inter_wet)
   q_wet=c_wet/(slope_wet+inter_wet)
r