2016-04-07 16 views
2

Bir şekle yaklaşmak için (benim durumumda, bir ayak izi anahattı) kısıtlı bir spline üretmek için bir olasılık arıyorum. Ham veriler olarak, ayak izinin sınırından toplanan birkaç yüz xy koordinat çiftine sahip bir tablom var. Spline sadece veri noktalarına yaklaşmalıdır (spline'ın veri noktalarını geçmesi gerekmez). Spline'ı belli derecelerde düzeltmek istiyorum. Ayrıca, spline'ı sınırlayabilmem gerekiyor: Spline'ın geçmesi gereken bazı kritik veri noktalarını tanımlamak. Kısıtlı B-spline'ları kullanarak bir şekil anahattını tahmin edin

R paketi "koçanı" (Kolejlendirilmiş B-Spline, https://cran.r-project.org/web/packages/cobs/index.html), spline'ı istediği gibi sınırlamak için bir çözüm sunmaya çok yaklaştı. Bununla birlikte, bu paket, belirli bir veri noktası dizisi boyunca yayılmaz; bu da, şeridin bir şeklin sınırını takip etmesini istediğinizde çok önemlidir. X ve y koordinatlarını ayrı ayrı kestirmeye çalıştım, ancak bunları birleştirdikten sonra iki farklı şekil arsada görünür, bu yüzden işe yaramaz (Ya da yanlış bir şeyim var mı?). Bir çözümden haberdar olan var mı?


Güncelleme: Çalışma örneği (dinozor ayak izi anahat)

data.txt:

structure(list(V1 = c(124.9, 86.44, 97.22, 81.34, 49.09, 57.18, 
-77.6, -191.95, -284.67, -383.18, -379.27, -492.85, -547.72, 
-600.67, -713.29, -814.36, -868.27, -926.99, -958.76, -1025.18, 
-1077.16, -1105.07, -1126.25, -1112.77, -1087.74, -989.56, -911.59, 
-859.61, -745.06, -656.5, -682.01, -637.25, -601.71, -539.09, 
-394.79, -219.17, -170.17, -201.48, -122.52, -43.56, 127.97, 
344.42, 539.09, 686.11, 987.63, 1253.31, 1283.15, 1536.32, 1741.14, 
1832.35, 1700.3, 1787.43, 1911.31, 2017.49, 2097.81, 2135.93, 
2093.73, 2066.96, 2063.78, 2022.94, 1978.69, 1919.44, 1904.03, 
1895.37, 1854.22, 1810.23, 1771.09, 1741.48, 1642.45, 1553.96, 
1472.96, 1396.04, 1141.65, 1085.82, 1055.02, 1358.24, 1325.94, 
1031.91, 1287.14, 1265.36, 931.15, 872.12, 811.48, 755.65, 738.32, 
697.41, 682.49, 647.35, 628.25, 620.09, 629.62, 675.22, 709.25, 
718.78, 717.42, 551.09, 535.21, 540.98, 534.73, 546.76, 811.96, 
823.03, 822.07, 607.4, 626.18, 637.73, 659.87, 756.13, 753.72, 
735.91, 720.99, 676.71, 576.6, 508.26, 339.8, 179.53, 121.16, 
45.6, 12.93, -9.87, -12.59, 16, 27.91, 37.78, 49.35, 8.51, 2.72, 
-1.02, 59.22, 58.2, 51.73, 54.45, 0.96, 10.59, 138.62, 149.69, 
144.87, 142.26, 146.34, 125.24, 124.9, 86.44, 97.22, 81.34, 49.09, 
57.18, -77.6, -191.95, -284.67, -383.18, -379.27, -492.85, -547.72, 
-600.67, -713.29, -814.36, -868.27, -926.99, -958.76, -1025.18, 
-1077.16, -1105.07, -1126.25, -1112.77, -1087.74, -989.56, -911.59, 
-859.61, -745.06, -656.5, -682.01, -637.25, -601.71, -539.09, 
-394.79, -219.17, -170.17, -201.48, -122.52, -43.56, 127.97, 
344.42, 539.09, 686.11, 987.63, 1253.31, 1283.15, 1536.32, 1741.14, 
1832.35, 1700.3, 1787.43, 1911.31, 2017.49, 2097.81, 2135.93, 
2093.73, 2066.96, 2063.78, 2022.94, 1978.69, 1919.44, 1904.03, 
1895.37, 1854.22, 1810.23, 1771.09, 1741.48, 1642.45, 1553.96, 
1472.96, 1396.04, 1141.65, 1085.82, 1055.02, 1358.24, 1325.94, 
1031.91, 1287.14, 1265.36, 931.15, 872.12, 811.48, 755.65, 738.32, 
697.41, 682.49, 647.35, 628.25, 620.09, 629.62, 675.22, 709.25, 
718.78, 717.42, 551.09, 535.21, 540.98, 534.73, 546.76, 811.96, 
823.03, 822.07, 607.4, 626.18, 637.73, 659.87, 756.13, 753.72, 
735.91, 720.99, 676.71, 576.6, 508.26, 339.8, 179.53, 121.16, 
45.6, 12.93, -9.87, -12.59, 16, 27.91, 37.78, 49.35, 8.51, 2.72, 
-1.02, 59.22, 58.2, 51.73, 54.45, 0.96, 10.59, 138.62, 149.69, 
144.87, 142.26, 146.34, 125.24), V2 = c(-446.8, -415.83, -394.43, 
-259.19, -104.69, -4.03, 58.59, -80.26, 52.11, -48.33, -142.23, 
-176.89, -233.68, -321.28, -416.57, -457.97, -458.93, -429.09, 
-422.35, -450.27, -431.98, -379.03, -260.63, -123.94, -2.65, 
269.76, 455.55, 548.92, 616.3, 691.38, 756.84, 888.72, 1016.97, 
1157.18, 1198.02, 1101.37, 1025.14, 929.84, 852.25, 766.48, 717.47, 
733.81, 784.18, 835.91, 1225.63, 1198.68, 925.3, 742.4, 814.13, 
732.45, 586.79, 394.84, 212.42, 28.64, -111.58, -337.56, -490.03, 
-526.07, -528.82, -547.2, -551.97, -552.3, -585.51, -551.34, 
-543.16, -526.1, -494.11, -466.88, -355.93, -274.94, -215.04, 
-114.3, -194.21, -103.73, -3.62, 104.2, 230.8, 154.25, 380.55, 
416.62, 260.07, 295.75, 295.75, 251.47, 220.67, 225.96, 180.72, 
121.52, 4.14, -127.23, -176.24, -332.11, -408.35, -494.11, -573.75, 
-582.62, -678.88, -730.38, -788.62, -831.94, -846.38, -895.95, 
-934.46, -968.15, -1033.12, -1097.62, -1150.08, -1157.3, -1254.04, 
-1340.2, -1441.75, -1500.47, -1550.52, -1605.39, -1681.44, -1709.84, 
-1715.22, -1672.34, -1607, -1522.59, -1440.57, -1421.18, -1345.62, 
-1247.95, -1190.77, -1181.58, -1071.65, -1037.62, -1010.39, -998.82, 
-986.57, -937.9, -887.29, -842.05, -831.46, -774.66, -703.91, 
-573.75, -533.59, -448.16, -446.8, -415.83, -394.43, -259.19, 
-104.69, -4.03, 58.59, -80.26, 52.11, -48.33, -142.23, -176.89, 
-233.68, -321.28, -416.57, -457.97, -458.93, -429.09, -422.35, 
-450.27, -431.98, -379.03, -260.63, -123.94, -2.65, 269.76, 455.55, 
548.92, 616.3, 691.38, 756.84, 888.72, 1016.97, 1157.18, 1198.02, 
1101.37, 1025.14, 929.84, 852.25, 766.48, 717.47, 733.81, 784.18, 
835.91, 1225.63, 1198.68, 925.3, 742.4, 814.13, 732.45, 586.79, 
394.84, 212.42, 28.64, -111.58, -337.56, -490.03, -526.07, -528.82, 
-547.2, -551.97, -552.3, -585.51, -551.34, -543.16, -526.1, -494.11, 
-466.88, -355.93, -274.94, -215.04, -114.3, -194.21, -103.73, 
-3.62, 104.2, 230.8, 154.25, 380.55, 416.62, 260.07, 295.75, 
295.75, 251.47, 220.67, 225.96, 180.72, 121.52, 4.14, -127.23, 
-176.24, -332.11, -408.35, -494.11, -573.75, -582.62, -678.88, 
-730.38, -788.62, -831.94, -846.38, -895.95, -934.46, -968.15, 
-1033.12, -1097.62, -1150.08, -1157.3, -1254.04, -1340.2, -1441.75, 
-1500.47, -1550.52, -1605.39, -1681.44, -1709.84, -1715.22, -1672.34, 
-1607, -1522.59, -1440.57, -1421.18, -1345.62, -1247.95, -1190.77, 
-1181.58, -1071.65, -1037.62, -1010.39, -998.82, -986.57, -937.9, 
-887.29, -842.05, -831.46, -774.66, -703.91, -573.75, -533.59, 
-448.16)), .Names = c("V1", "V2"), class = "data.frame", row.names = c(NA, 
-280L)) 
arsa

çıkan bir
require(cobs) 

xy <- dget(data.txt) 

#Cumchord function (from Claude, 2008): Cumulative chordal distance vector 
cumchord<-function(M) 
{cumsum(sqrt(apply((M-rbind(M[1,], 
M[-(dim(M)[1]),]))^2,1,sum)))} 

z <- cumchord(xy) 

#Calculating B-spline for x and y values separately 
x <- cobs(z,xy[,1],nknots=50) 
y <- cobs(z,xy[,2],nknots=50) 

#Plot spline 
plot(xy) 
lines(x$fitted,y$fitted) 

Görüntüenter image description here

+0

Bu ilginç bir soru olduğunu düşünüyorum, ama (a [tekrarlanabilir-örneğin] den http://stackoverflow.com/questions/5963269/how-to-make-a-great-r- yararlanacak tekrarlanabilir-örnek). İnsanların çalışabileceği bir şekil gösteren küçük bir veri paylaşabilir misiniz? – user20650

+0

Sadece merak ediyorum, örneğin, R_ (2008) pp. 208 ile Claude _Morphometrics'deki örneği çoğaltmaya/geliştirmeye mi çalışıyorsunuz? –

+0

@Vincent: İpucu için teşekkürler, tam olarak denediğim şey, ama daha zarif. Yukarıdaki çalışma örneğinde, şimdi Claude 2008'i izledim, ancak "spline" işlevi yerine "cobs" kullanmak için değiştirdim. "Spline" fonksiyonu ile mükemmel bir şekilde çalışır, ancak koçanlarla bu kablolu çift anahatları alırım. Ama neden? – tretelrusch

cevap

2

Yorumlar konusunu takiben, işte bazı grafikler. Tanıdığımdan beri Momocs kullanıyorum ve örnekleri kısaltacağım.

Özetle, sorun, anahatta iki ana hat olduğu.

Ayrıca Julien Claude tarafından spline orijinal kullanımı ve bezier eğrileri ve eliptik Fourier dönüşümleri ile iki ek örnek içerir. 4 tanesi bir taslağı tanımlamak (ve yeniden inşa etmek) için kullanılabilir ve muhtemelen onları burada toplamaya değer.

resim Şimdi kod orijinal şekil ve bu 4 yöntemler panel

toplar. Özellikle uzun değil ama tekrarlı.

# devtools::install_github("vbonhomme/Momocs") 
library(Momocs) # version 1.0.3 

xy <- structure(list(
    V1 = c(124.9, 86.44, 97.22, 81.34, 49.09, 57.18, -77.6, -191.95, -284.67, -383.18, -379.27, -492.85, -547.72, -600.67, -713.29, -814.36, -868.27, -926.99, -958.76, -1025.18, -1077.16, -1105.07, -1126.25, -1112.77, -1087.74, -989.56, -911.59, -859.61, -745.06, -656.5, -682.01, -637.25, -601.71, -539.09, -394.79, -219.17, -170.17, -201.48, -122.52, -43.56, 127.97, 344.42, 539.09, 686.11, 987.63, 1253.31, 1283.15, 1536.32, 1741.14, 1832.35, 1700.3, 1787.43, 1911.31, 2017.49, 2097.81, 2135.93, 2093.73, 2066.96, 2063.78, 2022.94, 1978.69, 1919.44, 1904.03, 1895.37, 1854.22, 1810.23, 1771.09, 1741.48, 1642.45, 1553.96, 1472.96, 1396.04, 1141.65, 1085.82, 1055.02, 1358.24, 1325.94, 1031.91, 1287.14, 1265.36, 931.15, 872.12, 811.48, 755.65, 738.32, 697.41, 682.49, 647.35, 628.25, 620.09, 629.62, 675.22, 709.25, 718.78, 717.42, 551.09, 535.21, 540.98, 534.73, 546.76, 811.96, 823.03, 822.07, 607.4, 626.18, 637.73, 659.87, 756.13, 753.72, 735.91, 720.99, 676.71, 576.6, 508.26, 339.8, 179.53, 121.16, 45.6, 12.93, -9.87, -12.59, 16, 27.91, 37.78, 49.35, 8.51, 2.72, -1.02, 59.22, 58.2, 51.73, 54.45, 0.96, 10.59, 138.62, 149.69, 144.87, 142.26, 146.34, 125.24, 124.9, 86.44, 97.22, 81.34, 49.09, 57.18, -77.6, -191.95, -284.67, -383.18, -379.27, -492.85, -547.72, -600.67, -713.29, -814.36, -868.27, -926.99, -958.76, -1025.18, -1077.16, -1105.07, -1126.25, -1112.77, -1087.74, -989.56, -911.59, -859.61, -745.06, -656.5, -682.01, -637.25, -601.71, -539.09, -394.79, -219.17, -170.17, -201.48, -122.52, -43.56, 127.97, 344.42, 539.09, 686.11, 987.63, 1253.31, 1283.15, 1536.32, 1741.14, 1832.35, 1700.3, 1787.43, 1911.31, 2017.49, 2097.81, 2135.93, 2093.73, 2066.96, 2063.78, 2022.94, 1978.69, 1919.44, 1904.03, 1895.37, 1854.22, 1810.23, 1771.09, 1741.48, 1642.45, 1553.96, 1472.96, 1396.04, 1141.65, 1085.82, 1055.02, 1358.24, 1325.94, 1031.91, 1287.14, 1265.36, 931.15, 872.12, 811.48, 755.65, 738.32, 697.41, 682.49, 647.35, 628.25, 620.09, 629.62, 675.22, 709.25, 718.78, 717.42, 551.09, 535.21, 540.98, 534.73, 546.76, 811.96, 823.03, 822.07, 607.4, 626.18, 637.73, 659.87, 756.13, 753.72, 735.91, 720.99, 676.71, 576.6, 508.26, 339.8, 179.53, 121.16, 45.6, 12.93, -9.87, -12.59, 16, 27.91, 37.78, 49.35, 8.51, 2.72, -1.02, 59.22, 58.2, 51.73, 54.45, 0.96, 10.59, 138.62, 149.69, 144.87, 142.26, 146.34, 125.24), 
    V2 = c(-446.8, -415.83, -394.43, -259.19, -104.69, -4.03, 58.59, -80.26, 52.11, -48.33, -142.23, -176.89, -233.68, -321.28, -416.57, -457.97, -458.93, -429.09, -422.35, -450.27, -431.98, -379.03, -260.63, -123.94, -2.65, 269.76, 455.55, 548.92, 616.3, 691.38, 756.84, 888.72, 1016.97, 1157.18, 1198.02, 1101.37, 1025.14, 929.84, 852.25, 766.48, 717.47, 733.81, 784.18, 835.91, 1225.63, 1198.68, 925.3, 742.4, 814.13, 732.45, 586.79, 394.84, 212.42, 28.64, -111.58, -337.56, -490.03, -526.07, -528.82, -547.2, -551.97, -552.3, -585.51, -551.34, -543.16, -526.1, -494.11, -466.88, -355.93, -274.94, -215.04, -114.3, -194.21, -103.73, -3.62, 104.2, 230.8, 154.25, 380.55, 416.62, 260.07, 295.75, 295.75, 251.47, 220.67, 225.96, 180.72, 121.52, 4.14, -127.23, -176.24, -332.11, -408.35, -494.11, -573.75, -582.62, -678.88, -730.38, -788.62, -831.94, -846.38, -895.95, -934.46, -968.15, -1033.12, -1097.62, -1150.08, -1157.3, -1254.04, -1340.2, -1441.75, -1500.47, -1550.52, -1605.39, -1681.44, -1709.84, -1715.22, -1672.34, -1607, -1522.59, -1440.57, -1421.18, -1345.62, -1247.95, -1190.77, -1181.58, -1071.65, -1037.62, -1010.39, -998.82, -986.57, -937.9, -887.29, -842.05, -831.46, -774.66, -703.91, -573.75, -533.59, -448.16, -446.8, -415.83, -394.43, -259.19, -104.69, -4.03, 58.59, -80.26, 52.11, -48.33, -142.23, -176.89, -233.68, -321.28, -416.57, -457.97, -458.93, -429.09, -422.35, -450.27, -431.98, -379.03, -260.63, -123.94, -2.65, 269.76, 455.55, 548.92, 616.3, 691.38, 756.84, 888.72, 1016.97, 1157.18, 1198.02, 1101.37, 1025.14, 929.84, 852.25, 766.48, 717.47, 733.81, 784.18, 835.91, 1225.63, 1198.68, 925.3, 742.4, 814.13, 732.45, 586.79, 394.84, 212.42, 28.64, -111.58, -337.56, -490.03, -526.07, -528.82, -547.2, -551.97, -552.3, -585.51, -551.34, -543.16, -526.1, -494.11, -466.88, -355.93, -274.94, -215.04, -114.3, -194.21, -103.73, -3.62, 104.2, 230.8, 154.25, 380.55, 416.62, 260.07, 295.75, 295.75, 251.47, 220.67, 225.96, 180.72, 121.52, 4.14, -127.23, -176.24, -332.11, -408.35, -494.11, -573.75, -582.62, -678.88, -730.38, -788.62, -831.94, -846.38, -895.95, -934.46, -968.15, -1033.12, -1097.62, -1150.08, -1157.3, -1254.04, -1340.2, -1441.75, -1500.47, -1550.52, -1605.39, -1681.44, -1709.84, -1715.22, -1672.34, -1607, -1522.59, -1440.57, -1421.18, -1345.62, -1247.95, -1190.77, -1181.58, -1071.65, -1037.62, -1010.39, -998.82, -986.57, -937.9, -887.29, -842.05, -831.46, -774.66, -703.91, -573.75, -533.59, -448.16)), .Names = c("V1", "V2"), class = "data.frame", row.names = c(NA, -280L)) 


### First thing first: double outline --------------------- 
coo_plot(xy) 
ldk_labels(xy) # blurry since superimposed 
coo_plot(xy[1:140, ], lwd=3) # first shape 
coo_draw(xy[-(1:140), ], border="white") # second shape 

# so from here, we will use the first 140th points from xy and name it 'shp' to avoid confusion 
# if you're bored with dinos footprints, you can use beer bottles with shp <- bot[9] for a guinness 

shp <- xy[1:140, ] 

### 1.Natural splines --------------------- 

shp_cumchord <- coo_perimcum(shp) # cumchord equivalent 
shp_spline <- cbind(spline(shp_cumchord, shp[, 1], method="natural", n=120)$y, 
        spline(shp_cumchord, shp[, 2], method="natural", n=120)$y) 

coo_plot(shp, main="natural spline", zoom=1.2) 
coo_draw(shp_spline, border="blue", lwd=2) 

### 2. B-splines with cobs --------------------- 
library(cobs) 

shp_bspline <- cbind(cobs(shp_cumchord, shp[, 1], nknots=50)$fitted, 
        cobs(shp_cumchord, shp[, 2], nknots=50)$fitted) 

coo_plot(shp, main = "bspline", zoom=1.2) 
coo_draw(shp_bspline, border="blue") 

### 3. Bezier curves --------------------- 
# built in function so it's shorter 
shp_bezier <- shp %>% bezier() %>% bezier_i() 

coo_plot(shp, main = "bezier", zoom=1.2) 
coo_draw(shp_bezier, border="blue") 


### 4. elliptic Fourier transforms --------------------- 
# another built in function 

shp_eft <- shp %>% efourier() %>% efourier_i() 

coo_plot(shp, main = "bspline", zoom=1.2) 
coo_draw(shp_eft, border="blue") 

### 5. A panel of original shape and 4 methods --------- 
Out(list(original=shp, 
     nat_spline=shp_spline, bspline=shp_bspline, 
     bezier=shp_bezier, eft=shp_eft)) %>% 
panel(names=TRUE, dim=c(1, 5)) 
+0

Ayrıca, eğer splineların Momocs'a dahil edilmesi gerektiğini düşünüyorsanız, nasıl ve neden tartışmaktan mutluluk duyarım. –

+0

Awesome, kapsamlı genel bakış için ve çift anahat ile hatayı bulmak için teşekkür ederiz! Momocs ile Eliptik Fourier aslında ilk denedim ve oldukça iyi çalıştı. Koçanı denemek istediğim şey, spline'ı sınırlama kabiliyetidir (yani, spline'ın geçmesi gereken noktaları tanımlar). Bu eliptik Fourier dönüşümleri ile de mümkün mü? – tretelrusch

+0

Bunları özellikle kısıtlamak için değil, yer işaretlerini kullanarak hizalamak mümkündür. Ancak, ana hatlar, 7-8 harmonikleri söylese de, (bu noktalardan teknik olarak geçen bir zamanda pikselden daha yakın) çok yakın olacaktır. –