如何使用R基于线性度分离实验扫频起始时间的索引段

我手里有一组实验中多个Protocol的单次扫频起始时间向量，希望通过起始时间的线性度区分不同的Protocol。绘制这个向量能直观看到哪些扫频是连续的，但不知道具体怎么通过线性度完成分离。

starting_times = c(1518.280, 1523.622, 1529.188, 1534.527, 1539.858, 1545.006, 1550.458, 1555.838, 1561.153, 1566.463, 1571.848, 1577.106, 1582.271, 1587.658, 1592.874, 1598.086, 1603.334, 1608.481, 1613.953, 1619.115, 1673.661, 1695.512, 1716.557, 1856.711, 1866.470, 1869.777, 1873.147, 1886.839, 1890.145, 1893.404, 1896.853, 1900.199, 1903.585, 1921.432, 1931.714, 1937.140, 1942.540, 1947.849, 1953.022, 1958.291, 1963.643, 1968.793, 2008.844, 2020.818, 2029.011, 2044.400, 2053.175, 2077.344)

plot(starting_times, type = "b", pch = 16, col = "steelblue", main = "扫频起始时间序列", xlab = "扫频序号", ylab = "起始时间")

分离方法 #

方法1：基于时间差突变分组 #

同一Protocol下的连续扫频，时间间隔相对稳定；不同Protocol切换时，会出现明显的时间差突变。可以通过计算相邻时间的差值，设定阈值划分分组：

# 计算相邻时间差
time_diff <- diff(starting_times)
# 取差值的95分位数作为突变阈值（可根据实际调整）
threshold <- quantile(time_diff, 0.95)
# 找出突变点位置
break_points <- which(time_diff > threshold)
# 生成分组标签
group_labels <- cumsum(c(1, time_diff > threshold))
# 整理带分组信息的数据框
starting_times_df <- data.frame(
  scan_index = seq_along(starting_times),
  start_time = starting_times,
  protocol_group = factor(group_labels)
)

# 可视化分组结果
plot(starting_times_df$scan_index, starting_times_df$start_time, 
     type = "b", pch = 16, col = starting_times_df$protocol_group,
     main = "按时间差突变分组后的扫频序列", xlab = "扫频序号", ylab = "起始时间")
legend("topleft", legend = levels(starting_times_df$protocol_group), 
       col = 1:nlevels(starting_times_df$protocol_group), pch = 16, title = "Protocol分组")

方法2：基于线性拟合残差分组 #

同一Protocol的起始时间近似线性（连续扫频按固定间隔执行），可以用滑动窗口线性拟合，通过残差大小判断是否属于同一组：

library(zoo)

# 定义滑动窗口大小（根据连续扫频数量调整）
window_size <- 5
# 计算每个窗口的线性拟合残差标准差
resid_sd <- rollapply(starting_times, width = window_size, 
                      FUN = function(x) {
                        fit <- lm(x ~ seq_along(x))
                        sd(residuals(fit))
                      }, align = "center")
# 补全首尾缺失值
resid_sd <- c(rep(resid_sd[1], floor(window_size/2)), 
              resid_sd, 
              rep(resid_sd[length(resid_sd)], floor(window_size/2)))
# 设定残差阈值划分分组
resid_threshold <- quantile(resid_sd, 0.9)
# 生成分组标签
group_labels_resid <- cumsum(c(1, diff(c(0, resid_sd > resid_threshold)) == 1))
starting_times_df$protocol_group_resid <- factor(group_labels_resid)

# 可视化残差分组结果
plot(starting_times_df$scan_index, starting_times_df$start_time, 
     type = "b", pch = 16, col = starting_times_df$protocol_group_resid,
     main = "按线性拟合残差分组后的扫频序列", xlab = "扫频序号", ylab = "起始时间")
legend("topleft", legend = levels(starting_times_df$protocol_group_resid), 
       col = 1:nlevels(starting_times_df$protocol_group_resid), pch = 16, title = "Protocol分组")

方法3：聚类分析 #

如果时间间隔变化复杂，可通过K-means聚类对时间差特征自动划分分组：

set.seed(123)
# 用时间差作为特征进行K-means聚类，聚类数根据实际情况调整
kmeans_clust <- kmeans(time_diff, centers = 5)
# 扩展聚类标签到每个扫频点
group_labels_clust <- c(1, kmeans_clust$cluster)
starting_times_df$protocol_group_clust <- factor(group_labels_clust)

# 可视化聚类结果
plot(starting_times_df$scan_index, starting_times_df$start_time, 
     type = "b", pch = 16, col = starting_times_df$protocol_group_clust,
     main = "K-means聚类分组后的扫频序列", xlab = "扫频序号", ylab = "起始时间")
legend("topleft", legend = levels(starting_times_df$protocol_group_clust), 
       col = 1:nlevels(starting_times_df$protocol_group_clust), pch = 16, title = "Protocol分组")

内容的提问来源于stack exchange，提问作者nrsc