2011-04-15 51 views
2

我想從這些數據中選擇一組距離值爲< = 40和激活值< = 15的點。他們必須是連續的,並且必須至少有3個。我開始使用的代碼是rle()命令。我無法弄清楚如何將這些信息連接回數據 - 將這兩個條件> = 3的那些「運行」拉出來。我使用了哪個命令和cbind來查看並檢查rle是否正常工作。那麼我如何獲得這些組的相關信息?然後讓事情變得更加複雜......我需要一次完成超過2天的信息。這意味着,從下面的示例數據中可以看出,我有一個日期/時間,我需要在48小時的窗口內找到所有> = 3個位置的組。窗戶是排他性的,在午夜開始,48小時後在午夜結束。我在下面提供的樣本數據不足以滿足最後一部分的要求,但如果有必要,我願意將其過去。使用rle在R中選擇一組連續數據

run <- rle(Dist <= 40 & Activ <= 15) 

tempdist <- Dist <= 40 
tempactive <- Activ <= 15 
which(tempdist == TRUE & tempactive == TRUE) 
combine <- cbind(tempdist, tempactive) 

structure(list(Fix = 1:15, Date = structure(c(1L, 1L, 1L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("2010.06.23", 
"2010.06.24", "2010.06.25", "2010.06.26", "2010.06.27", "2010.06.28", 
"2010.06.29", "2010.06.30"), class = "factor"), Time = structure(c(86L, 
89L, 95L, 6L, 11L, 69L, 106L, 107L, 113L, 121L, 123L, 128L, 136L, 
14L, 22L), .Label = c("0:00:36", "0:00:39", "0:00:42", "0:00:48", 
"0:01:04", "0:01:28", "0:01:35", "1:00:15", "1:00:19", "1:00:49", 
"1:01:05", "1:01:29", "10:00:13", "10:00:35", "10:00:36", "10:00:42", 
"10:01:05", "10:01:11", "11:00:15", "11:00:20", "11:00:35", "11:00:42", 
"11:00:57", "11:01:05", "12:00:13", "12:00:19", "12:00:36", "12:00:42", 
"12:01:39", "12:03:03", "13:00:13", "13:00:32", "13:00:36", "13:00:47", 
"13:01:05", "13:01:52", "14:00:18", "14:00:36", "14:00:48", "14:01:04", 
"14:03:04", "15:00:35", "15:00:36", "15:01:04", "15:01:28", "15:01:36", 
"16:00:18", "16:00:19", "16:00:20", "16:00:36", "16:00:45", "16:03:05", 
"17:00:35", "17:00:42", "17:01:05", "17:01:11", "18:00:13", "18:00:35", 
"18:00:36", "18:00:37", "18:00:48", "18:01:04", "19:00:13", "19:00:27", 
"19:00:30", "19:00:35", "19:00:42", "19:00:43", "2:00:18", "2:00:36", 
"2:00:47", "2:00:54", "2:00:59", "2:01:30", "2:02:37", "20:00:13", 
"20:00:35", "20:00:36", "20:01:05", "21:00:14", "21:00:19", "21:00:20", 
"21:00:21", "21:00:35", "21:00:48", "21:01:08", "21:02:23", "22:00:35", 
"22:00:36", "22:00:38", "22:01:36", "22:03:04", "22:03:05", "23:00:15", 
"23:00:16", "23:00:20", "23:00:36", "23:01:05", "23:02:05", "23:03:03", 
"3:00:35", "3:00:36", "3:00:47", "3:00:48", "3:00:59", "3:01:58", 
"4:00:13", "4:00:17", "4:00:19", "4:00:35", "4:00:39", "5:00:16", 
"5:00:35", "5:00:36", "5:01:23", "6:00:13", "6:00:16", "6:00:18", 
"6:00:19", "6:00:35", "6:00:36", "7:00:18", "7:00:19", "7:00:36", 
"7:01:05", "7:01:22", "8:00:18", "8:00:35", "8:00:36", "8:00:42", 
"8:01:05", "9:00:13", "9:00:19", "9:00:35", "9:00:36", "9:00:53", 
"9:01:36"), class = "factor"), X = c(NA, 6351162.121, 6351137.038, 
6351180.514, 6351175.284, 6351189.867, 6351168.96, 6351230.592, 
6351187.825, 6351182.809, 6351179.366, 6351181.313, 6351157.538, 
6351180.688, 6351178.51), Y = c(NA, 436904.9145, 436928.697, 
436925.2068, 436946.703, 436961.3917, 436967.0372, 436934.0528, 
436978.3563, 436949.7116, 436963.6412, 436967.4078, 436949.5642, 
436964.5655, 436961.398), Temp = c(14L, 19L, 16L, 14L, 18L, 20L, 
18L, 14L, 13L, 15L, 17L, 13L, 16L, 18L, 20L), Activ = c(0.5, 
0.5, 57.5, 21.5, 0.5, 1, 48, 62.5, 42.5, 2, 9.5, 0.5, 0.5, 0.5, 
0.5), DateTime = structure(list(sec = c(8, 36, 16, 28, 5, 18, 
58, 13, 35, 36, 19, 35, 53, 35, 42), min = c(1L, 0L, 0L, 1L, 
1L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L), hour = c(21L, 22L, 
23L, 0L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L), mday = c(23L, 
23L, 23L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 
24L), mon = c(5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 
5L, 5L, 5L), year = c(110L, 110L, 110L, 110L, 110L, 110L, 110L, 
110L, 110L, 110L, 110L, 110L, 110L, 110L, 110L), wday = c(3L, 
3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L), yday = c(173L, 
173L, 173L, 174L, 174L, 174L, 174L, 174L, 174L, 174L, 174L, 174L, 
174L, 174L, 174L), isdst = c(0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 
0L, 0L, 0L, 0L, 0L, 0L, 0L)), .Names = c("sec", "min", "hour", 
"mday", "mon", "year", "wday", "yday", "isdst"), class = c("POSIXlt", 
"POSIXt"), tzone = "UTC"), DistX = c(NA, -25.0830000005662, 43.4760000007227, 
-5.23000000044703, 14.5829999996349, -20.9069999996573, 61.6320000002161, 
-42.7669999999925, -5.01599999982864, -3.4429999999702, 1.94699999969453, 
-23.7750000003725, 23.1500000003725, -2.17800000030547, -7.12799999956042 
), DistY = c(NA, 23.7824999999721, -3.49020000000019, 21.4961999999941, 
14.6886999999988, 5.64550000004238, -32.9844000000157, 44.30349999998, 
-28.6447000000044, 13.9296000000322, 3.76659999997355, -17.8435999999638, 
15.0013000000035, -3.16750000003958, 0.131600000022445), Dist = c(NA, 
34.5653612056503, 43.6158694984158, 22.1232799205819, 20.6983525112342, 
21.6558149058434, 69.9033194303938, 61.5777265027509, 29.0805621350161, 
14.3487980388844, 4.24005714096063, 29.72616164217, 27.5855306584331, 
3.84405258179195, 7.12921472209522), LnDist = c(NA, 3.54285205938314, 
3.77542106345205, 3.09663044416928, 3.0300541082747, 3.07527400526048, 
4.24711313638682, 4.12030022240512, 3.37006998311865, 2.66366617834105, 
1.44457674579472, 3.39202752144892, 3.31729138374919, 1.34652716994060, 
1.96420109120712), TimeDif = c(59.4666666666667, 59.6666666666667, 
61.2, 59.6166666666667, 59.2166666666667, 61.6666666666667, 58.25, 
60.3666666666667, 60.0166666666667, 59.7166666666667, 60.2666666666667, 
60.3, 59.7, 60.1166666666667, 60), Velocity = c(NA, 0.579307729703637, 
0.712677606183265, 0.371092198835593, 0.349535927575022, 0.351175376851514, 
1.20005698592951, 1.02006173113337, 0.484541440738952, 0.240281295655334, 
0.0703549304362937, 0.492971171511941, 0.462069190258511, 0.063943209, 
0.118820245368254), Heading = c(NA, 5.47117955593762, 1.65090330804596, 
6.04452351734706, 0.781787182793282, 4.97612800901505, 2.06219261575797, 
5.51543193207954, 3.31494594980981, 6.04087053685035, 0.477085335105382, 
4.06855550443118, 0.995825324541324, 3.74395368878987, 4.7308492849138 
)), .Names = c("Fix", "Date", "Time", "X", "Y", "Temp", "Activ", 
"DateTime", "DistX", "DistY", "Dist", "LnDist", "TimeDif", "Velocity", 
"Heading"), row.names = c(NA, 15L), class = "data.frame") 

回答

1

我相信有更有效的方法來做到這一點,但它是星期五,這種方法的工作原理。 rle()返回一個包含長度和值的列表對象。正如jimmyb提到的那樣,要把從rle()返回的內容轉化爲對你的data.frame有用的東西。 rep()是做這件事的武器。

#Your rle code 
run <- rle(dat$Dist <= 40 & dat$Activ <= 15) 

#Add columns for the length of the run and the value of that run 
dat$run <- rep(run$lengths, run$lengths) 
dat$cond <- rep(run$values, run$lengths) 

#Subset the data based on the length and value 
subset(dat, run >= 3 & cond) 

根據樣本數據規定,這將返回第10 - 15行。這是您的預期嗎?

+0

是的,它確實返回了正確的信息行......但是......我想要的是1)在一個時間窗內將數據按組2聚類)。所以,它提出正確的信息。也許我應該在示例中添加更多行數據,因爲我仍然需要在24或48小時窗口中查找和識別這些羣集。 – Kerry 2011-04-18 04:26:03

1

使用RLE這樣做的一個很好的方式是通過以下方式:

v = c(rep(1, 3), rep(2, 4), 4, rep(6, 2)) 
r = rle(v) 

現在,您可以創建矢量索引到原來的結構。關鍵是rep功能。

tapply(v, rep(r$lengths, r$lengths), function(a) { 
    # a is a vector of all runs of a certain length. 
}) 

tapply(v, rep(1:length(r$lengths), r$lengths), function(a) { 
    # a is a vector of each run. 
    if (length(a) > 4) 
    # do something 
    else 
    # do something else 
})