ข้ามไปที่เนื้อหาหลัก

RANdom SAmple Consensus (RANSAC) algorithm

Introduction

fig. 1 picture from 'Multiview geometry in Computer Vision' book.

If you use orthogonal regression ( minimizes the sum of squared perpendicular distances -- LMS ), there will be a problem if you have an outliers (see fig 1a).

RANSAC algorithm

RANSAC algorithm will cope with this problem by discarding outliers.

Slide from 25th year of RANSAC, Philip Torr slides has very clear picture of the algorithm.






support = number of points that lie within a distance threshold

points within the threshold distance of a line with most support are the inliers.

If a point is an outliers, a line will not have so much support. ( see fig 1b from mvg book above )

Explain the algorithm

- First we randomly pick two red point and estimate m, c for y=mx+c ( this is easy, right )

- for consider if a point is a inlier
for every yellow point (x, y)
If | y - (m*x + c) | < t, number of inlier need to justify model

-- re-estimate m, c again by using all the inlier

sample code ?

a simple python example code for RANSAC is available here.
or more complicate code here. all is the same code as a psudocode from wiki.

or you can use my porting octave code here ( NOTE if you are using matlab, this code should be modified )

# Run RANSAC to see if we can recover the line (y=x) from the data.
function [model_params, model_error, model_inliers] = ransac(data,n,k,t,d)
  % k -- the number of iterations.
  % t -- the threshold for deciding when a data point fits a model (i.e. is an inlier).
  % d -- the number of inliers needed to justify the model.

# k=100,d=20
  model_params = [0.0, 0.0] # slope and intercept
  model_error = 1e6
  model_inliers = []

  for i=1:k # k iterations
      fprintf(1, "new iteraion\n\n")
      inliers = choice(data, n)
   
      # Fit a line (or generally any model) using lstsq.
      a = [inliers(:,1), ones(1, n)']
      b = inliers(:,2)

      # linear least square -- rewrite y = mx+c in y = A*p form
      # A = [x 1]      
      # p = [m, c]
      p = a \ b
    
      compatible = []
      for j=1:size(data, 1)
          fprintf(1, '--')
          pt = data(j,:)
          if abs(pt(2) - dot(p, [pt(1), 1.0])) < t # | y - (m*x + c) | < t
              compatible = [compatible; pt]
          end  
      end  

      if size(compatible, 1) > d
          # The current model is good enough so we should recompute it using all compatible points.
          a = [compatible(:,1), ones(1, size(compatible, 1))']
          b = compatible(:,2)
          p = a \ b
       
          # if residuals < model_error
              model_params = p
              # model_error = residuals
              model_inliers = compatible
          # end
      end
  end

function result=choice(seq, n)
  indx = 1:size(seq)-1
  shuffle(indx)
  result = seq(indx(1:n), :)


adaptively find k (the number of iterations)
k = log(1-p)/log(1-(1-ε)^n)

proof ( ref : wiki ) :

p
= prob(RANSAC algorithm in some iteration selects only inliers )

w = prob( choosing an inlier each time a single point is selected ) = number of inliers / number of all points

A common case is that w is not well known beforehand, but some rough value can be given.

n = Sample Size

wn = prob( n points are inliers )
1 − wn = prob( at least one of the n points is an outlier ) [ bad model will be estimated from this point set ]
( 1 − wk = prob( algorithm never selects a set of n points which all are inliers )
( 1 − wn )k = 1 − p
take log and we will get
k = log(1-p)/log(1-w^n)

It should be noted that this result assumes that the n data points are selected independently, that is, a point which has been selected once is replaced and can be selected again in the same iteration.

Apply to Homography
When adapt this to Homography, estimating line is not a visual line. The error is computing from transfer error d.

This RANSAC for Homography steps are adapt from Alexei (Alyosha) Efros's slide.

RANSAC loop:
1. Select 4 feature pairs (at random)
[ n in the code, or s in the book is fixed to 4 ]
2. Compute homography H (exact)
3. If d(x’, H x) largest set of inliers
[ so the number of inliers needed to justify the model ( d in the code or capital T in the book ) is not used here. ]
5. Re-compute least-squares H estimate on all of the inliers

ความคิดเห็น

โพสต์ยอดนิยมจากบล็อกนี้

ปีอธิกสุรทิน อธิกมาส อธิกวาร และ การทดปฏิทิน

 ปีที่แบ่งตามรอบพระอาทิตย์ (สุริยคติ) มี 2 แบบ คือ 1. ปกติสุรทิน คือ มี 365 วัน 2. อธิกสุรทิน คือ มี 366 วัน (เดือนกุมภาฯมี 29 วัน) การคำนวน ปีอธิกสุรทิน ตรงนี้บางท่านจะจำได้แต่เพียง ว่า หาร 4 ซึ่งไม่ใช่แค่นั้นครับ ความจริงแล้ว จะมีสูตรคำนวณที่ถูกต้องคือ ให้เอา ค.ศ.ตั้ง แล้วเอา 4 หาร หากหารลงตัวก็ใช่ ยกเว้น 100 หารลงตัว แต่หาก 400 ลงตัวก็ให้นับเป็นอธิกสุรทินด้วย (เช่น ปี 1900 ไม่เป็นอธิกสุรทิน แต่ปี 2000, 2004 เป็นปีอธิกสุรทิน) ปีที่แบ่งตามรอบพระจันทร์ (จันทรคติ) แบ่งเป็น 3 คือ 1. ปกติมาส-ปกติวาร  (บางที่เขียนย่อ เป็น ปกติมาส-วาร)       คือ ปีที่เป็นปกติ มีเดือนคู่ ข้างขึ้น 15 วัน ข้างแรม 15 วัน       และมีเดือนคี่ ข้างขึ้น 15 วัน ข้างแรม 14 วัน       รวมวันใน 1 ปี เป็น (30*6+29*6) = 354 วัน 2. ปกติมาส-ปีอธิกวาร (บางที่เรียกเป็น อธิกวาร)       คือ ปีที่เป็นปกติ แต่เดือน 7 จะมีข้างแรม 15 วัน        รวมวันใน 1 ปี เป็น 354+1 = 355 วัน 3. ปีอธิกมาส-ปกติวาร (บางที่จะเรียกเป็น อธิกมาส)       คือ ปีที่มีเดือนแปดเพิ่มอีกเดือน หรือที่เรียกกันว่า มีเดือนแปดสองหน       รวมวันใน 1 ปี เป็น 354+30

วิธีการไป อย. กระทรวงสาธารณสุขจากหัวลำโพง

ทางไป : รถไฟฟ้า MRT หัวลำโพง ไปลงที่ สถานี กระทรวงสาธารณสุข  ถ้ากดที่ตู้ต้องเปลี่ยนไปหน้าจอสายสีม่วง สนน ราคา 48 53 บาท ต่อมอไซด์ ถ้าไป อย. 20 บาท จากหน้าทางเข้า  ถ้าฝนตกแนะนำให้โบกแท็กซี่จากข้างหน้า ข้างในหาแท็กซี่ยากมาก ถ้าจะเดินประมาณ 2.4 km ให้ระวังหลงเข้าไปรพ ศรีธัญญา รพ ศรีธัญญาพื้นที่ข้างในใหญ่มาก และเหมือนจะล้อมด้วยคลอง เหมือนจะมีทางออกแค่ทางที่เข้าไปนั่นแหละ ทางกลับ : รถเมล์ 97 จาก อย. ตรงข้ามประกันสังคม ทางที่ 1 : ถ้าจะใกล้ลงหน้าปากซอยขึ้นสายสีม่วงที่สถานีกระทรวงสาธารณสุขที่เดิม ทางที่ 2 : ผ่าน ท่าน้ำนนท์​ กลับเรือได้ ทางที่ 3 : ผ่านหน้าพระจอมพระนครเหนือด้วยนะ ผ่าน สถานีรถไฟฟ้า MRT บางซื่อ  ( จาก อย. ไป MRT บางซื่อ 17 บาท,  จาก MRT บางซื่อ ไป MRT หัวลำโพง 44 บาท นั่งกลับได้ 2 ทาง ทางหัวลำโพง กับ ไปเปลี่ยนที่ท่าพระ ไม่รู้ว่าทางไหนเร็วกว่ากัน ) ทางที่ 4 :  ผ่าน สะพานควาย  ทางที่ 5 :  นั่งถึงอนุสาวรีย์ชัยสมรภูมิได้ ค่ารถเมล์ 21 บาท ค่ารถไฟฟ้าไป BTS สะพานตากสิน 47 บาท 

สอบสัมภาษณ์ MBA คำถามและการเตรียมตัว

 * “แนะนำตนเอง” การแนะนำตนเองไม่ใช่แค่บอกชื่อ-นามสกุล ตำแหน่งงาน สถานที่ทำงาน หรือ ประวัติการศึกษาเท่านั้น ข้อมูลเหล่านี้ต้องพูดถึง แต่ไม่ใช่ประเด็นสำคัญ ส่วนที่สำคัญในการแนะนำตนเองก็คือต้องขายความเป็นตัวตนของเรา ความสามารถของเรา และ/หรือวัตถุประสงค์ในการเลือกเรียนหลักสูตรนี้  พยายามตอบคำถามให้สอดคล้องกับ MBA ไม่ต้องนาน ประมาณ 2–3 นาที เน้นเนื้อ ไม่เน้นน้ำ ซ้อมพูดเยอะๆ ถือว่าเป็น First Impression * ทำไมจึงเลือกสมัครเข้าเรียนหลักสูตรนี้  ทำไมถึงมาเรียน MBA ทำไมอยากเรียน MBA ทำไม อยากเรียนตอนนี้  * ทำไม ต้องเรียน MBA ที่นี่ -- ลองศึกษา Program ของมหาลัยที่จะไปดูน้าว่ามหาลัยมีอะไรเด่น * คิดว่าถ้าเรียน MBA จะมี Challenge อะไรบ้าง * สนใจโปรแกรมอะไรบ้าง * หลังเรียนจบอยากทำอะไร * ต้องการอะไรจากหลักสูตรนี้  เรียนแล้วคิดว่าจะได้อะไร เอาไปใช้อะไรในชีวิต * ทำไมไม่เรียนสาขาอื่น ถ้าอายุงานถึงเรียนอย่างอื่นได้ * ในองค์กรที่ทำงานอยู่สามารถเติบโตได้ถึงตำแหน่งไหน * Performance ปัจจุบันเป้นยังไง  * ดูดีอยู่แล้ว แล้วมาเรียน MBA ทำไม เพราะงานที่ทำอยู่ก็มีโกาสก้าวหน้าในสายอาชีพบริหารอยู่แล้ว * ไม่ได้เรียนม

วิธีใช้ ubuntu ต่อ อินเทอร์เน็ตทรู ( true ) โดยโมเด็ม billion bipac 7000 usb adsl modem

ก็อปไฟล์ cxacru-fw.bin ไปที่ /lib/firmware ไฟล์ cxacru-fw.bin download ได้ที่นี่ ก็อปไฟล์ br2684ctl ไปที่ /usr/sbin ไฟล์ br2684ctl download ได้ที่นี่ $ sudo pppoeconf nextๆ ไปเรื่อยๆ ใส่ username, password ของทรู ตามปกติ แล้วเขียนไฟล์ดังนี้ true.sh #!/bin/sh modprobe cxacru modprobe br2684 sudo /usr/sbin/br2684ctl -b -c 0 -a 0.100 # Communicating over ATM 0.0.100, encapsulation: LLC sudo ifconfig nas0 up pon dsl-provider # Plugin rp-pppoe.so loaded เสร็จแล้วสั่ง รัน shell script $ . ./true.sh คราวต่อไปรัน . ./true.sh อย่างเดียวก็ได้แล้วๆ reference : siamgeek บทความอื่นๆเกี่ยวกับ ubuntu

เรื่อง matrix ที่อาจจะลืมกันไปแล้ว

Rank ของ matrix Rank ของ matrix A คือ จำนวน independent columns (หรือ rows) ของ A นั่นคือ square matrix จะ full rank ถ้า ทุกคอลัมน์ independent กัน เมื่อ full rank, det จะ = 0 วิธีหา rank อาจหาได้โดย [U, W, V] = svd(A) แล้วดูว่า rank คือ จำนวน residual ของ W ที่ไม่เป็น 0 full rank = singular matrix = หา inverse ได้ สมบัติของ rank 1. rank(AB) <= min(rank(A), rank(B)) ย้ำว่า < หรือ = นะ ดูสมบัติอื่นๆ ได้จาก wiki Null Matrix เมตริกซ์ศูนย์ (Zero Matrix หรือ Null Matrix ) คือ เมตริกซ์ที่มีสมาชิกทุกตัวเป็นศูนย์หมด Orthogonal Matrix Cramer's rule Ax = b Cramer's rule ใช้ได้เมื่อ A เป็น square matrix เท่านั้น กรณีที่มีจำนวน สมการ มากกว่าจำนวน ตัวแปร ( A mxn เมื่อ m > n ) หรือเราเรียกว่า over parameter เราไม่สามารถหา inverse ของ rectangular matrix ได้ ให้ไปใช้ psudoinverse แทน x = A + b หรือ หรือ หรือ ไปใช้ SVD แก้สมการซะ คำตอบคือ last col of v ! Gaussian elimination method ใช้แก้สมการ เช่นเดียวกับ กฏของ คราเมอร์ วิธีคิดหลักๆ คือ ทำให้สามเหลี่ยมล่างเป็น 0 ให้หมด โดย

เทคนิคคิดเลขเร็วโดยใช้ วิธีคิด แบบ เวทคณิต ( Vedic Mathematics example )

จากที่สงสัยเรื่อง ลูกคิด ของ จินตคณิต ที่ลองไปค้นดู ปรากฎว่า เจอ เวทคณิต ซึ่งเขาบอกว่า อยู่ในคัมภีร์พระเวท ลองอ่านดูแล้ว รู้สึกว่าฝึกสมอง ก็ทำให้คิดเลขเร็วดี เลยสรุปมาให้ ตามนี้ Tutorial 1 การลบเลข ALL FROM 9 AND THE LAST FROM 10 ทุกตัวลบจาก 9 และตัวสุดท้ายลบจาก 10 เช่น 1000 - 357 = 643 10,000 - 1,049 = 8951 ถ้า 1,000 - 83 ให้มองว่ามี 0 อยู่ข้างหน้า เป็น 1,000 - 083 = 917 ฝึกบ่อยๆ ก็คล่อง แล้วก็ไม่ต้องใช้เครื่องคิดเลขด้วย ลองทำดูสิ 1) 1000 - 777 = 2) 1000 - 283 = 3) 1000 - 505 = 4) 10,000 - 2345 = 5) 10,000 - 9876 = 6) 10,000 - 1011 = 7) 100 - 57 = 8) 1000 - 57 = 9) 10,000 - 321 = 10) 10,000 - 38 = 3,000 - 467 ก็ทำเหมือนกัน โดยลบตัวแรกสุดของ 3,000 ไป 1 จากนั้นก็ทำเหมือนเดิม จะได้ว่า 3,000 - 467 = 2,533 Tutorial 2 VERTICALLY AND CROSSWISE สำหรับตัวเลขที่น้อยกว่าฐานนิดหน่อย ลอง 88x98 88 น้อยกว่า 100 อยู่ 12 98 น้อยกว่า 100 อยู่ 2 12x2 = 24 88-2 หรือ 98-12 ได้ 86 ดังนั้นตอบ 8,624 ดูอีกตัวอย่าง หรือ ลองทำนี่ดู 1) 87 x 98 = 2) 88 x

เลขฐานสอง ติดลบ เรื่องที่อาจจะลืมกันไปแล้ว

คอมพิวเตอร์ใช้การเปิดปิด หลอดสุญญากาศ ดังนั้นค่าที่เป็นไปได้คือ 0 กับ 1 ไม่มีติดลบ จึงกำหนดให้ใช้ 2's complement มากำหนดเลขลบ วิธีทำคือ เปลี่ยนเลข 1 เป็น 0 เปลี่ยนเลข 0 เป็น 1 แล้ว บวกหนึ่ง เช่น 1 คือ 00000001 เปลี่ยนเป็น 11111110 บวก 1 ได้ 11111111 บิตที่อยู่หน้าสุดจะบอกว่าเป็นเลขบวกหรือลบ ( 0 = +, 1 = -) พิสูจน์ จาก สมการคณิตศาสตร์​ 1 + (-1) = 0 00000001 + ???????? = 0 00000001 + (11111110 + 000000001 ) = 0 นั่นเอง วิธีที่ง่ายกว่านั้นในการทำ 2's complement คือ 1. หา 1 ตัวสุดท้าย 010100 1 2. invert ตัวหน้า 1 ทั้งหมด 101011 1 สำหรับคนที่ลืมไปแล้ว 1's complement คือเปลี่ยนเลข 1 เป็น 0 เปลี่ยนเลข 0 เป็น 1 ตามปกติ เช่น ~1 1 = 00000001 ~1 = 11111110 ซึ่งมีค่าเท่ากับ -2 ที่มา : วิชาการดอทคอม , wikipedia

ส่งไปรษณีย์ทีละมากๆ ที่ช่องไปรษณีย์สำหรับธุรกิจ

  ถ้าเราส่งไปรษณีย์ทีละ 10 กล่องขึ้นไป สามารถไปส่งโดยใช้ช่องทางธุรกิจได้ โดยต้องกรอกใบรับฝากรวม ( Receipt for bulk Posting ) เป็นลิสต์รายการให้เขาไปด้วย โดยกรอกพัสดุแต่ละรายการ และ ไปยื่นให้เขาพร้อมกับพัสดุที่จะส่ง วิธีกรอก คือ ให้กรอกพัสดุแบบเดียวกันไว้แผ่นเดียวกัน  เช่น พัสดุ10 กล่อง กล่องขนาดเท่ากัน น้ำหนักเท่ากันหมด กรอกไว้ 1 แผ่น ถ้าน้ำหนักต่างกัน ขนาดกล่องต่างกัน กรอกแยกแผ่นไว้ดีที่สุด ซึ่งใบนี้สามารถไปขอได้ที่ไปรษณีย์ฝ่ายธุรกิจ สามารถนำมาทำใส่ A4 ก็ได้ ขอบคุณคุณพี่ amarin.ch ที่ไปรษณีย์กลาง ( BANGKOK G.P.O. ) มากๆ นะครับ สำหรับคำแนะนำ ขอบคุณที่ช่วยคีย์ให้ทีละรายการสำหรับมือใหม่ที่ยังไม่รู้ว่ามีใบรับฝากรวมอย่างผมด้วยครับ คราวหน้าผมจะทำใบรับฝากรวมไปครับ

ข้อแตกต่างระหว่าง Mahalanobis distance กับ Euclidean Distance : ทำไม และ เมื่อไหร่ ต้องใช้ Mahalanobis distance ใน data mining

Euclidean Distance นิยาม EuclideanDistance = sqrt(sum( (A - B) .^ 2 )) โชว์เหนือ เขียนแบบ linear algebra EuclideanDistance = norm(A - B) ข้อเสียของ Euclidean distance 1. sensitive to scales ของตัวแปร ในกรณี geometric ตัวแปรทุกตัวมีหน่วยเดียวกันหมด คือ ระยะทาง แต่เมื่อพิจารณาตัวแปรที่มีข้อมูลหลายชนิดพร้อมๆกัน เช่น ใน data mining เราอาจจะพิจารณา อายุ, ความสูง, น้ำหนัก ฯลฯ พร้อมๆกันหมด สเกลมันเอามาเปรียบเทียบกันไม่ได้ 2. Euclidean distance ใช้กับตัวแปรที่ correlated กันไม่ได้ เช่น สมมติว่าเรามี data set 5 ตัวแปร ที่ซึ่งค่าของตัวแปรหนึ่งเหมือนกับอีกตัวแปรหนึ่งเด๊ะๆ ( กรณีนี้เหมือนเด๊ะ เลยเป็น completely correlated ) Euclidean distance จะคำนวณโดย weight ข้อมูลที่ซ้ำกันมากขึ้น ทำให้มีปัญหา Mahalanobis distance นิยาม เมื่อ S คือ covariance matrix และ x, y มี distribution เดียวกัน Mahalanobis distance มันพิจารณ่า covariance matrix ไปด้วย เลยขจัดปัญหาเรื่อง scale และ correlation ที่ Euclidean Distance มีได้ ใน MATLAB ใช้ฟังก์ชั่น mahal() หรือ pdist() ดูตัวอย่าง mahaldist.m ของคุณ Peter J. Ackl

อยู่เหงาๆ เราไปเที่ยว - ไหว้พระขอพร ศาลเจ้าแม่ทับทิม (อาม่า), เจริญกรุง, กรุงเทพ; 天后聖母廟, 石龙軍路, 曼谷, 泰国; Thap Thim Chinese Goddess Shrine, Chareon Krung 63 Road, Bangkok, Thailand

天后聖母廟, 石龙軍路, 曼谷, 泰国 ไหว้ศาลเจ้าแม่ทับทิม ขอให้การค้าเจริญรุ่งเรือง ตำนานเจ้าแม่ทับทิมเกิดที่ตำบลตุ้ยบ๊วย เขตบ่นเซียว เกาะไหหลำ มีผู้เฒ่าแซ่พัว เป็นผู้มีความซื่อสัตย์สุจริต ทำงานขยันขันแข็ง ครั้งหนึ่งแกออกไปหาปลา โดยผูกแหเป็นช้อนดักปลา เวลาผ่านไปแกยังหาปลาไม่ได้ คืนนั้นก็ประสบความล้มเหลว เมื่อช้อนแหขึ้นมาทีไรก้อมีแต่ท่อนไม้ ด้วยความโมโหแกเลยขว้างท่อนไม้นั้นออกไปให้ไกล แต่แล้วเมื่อช้อนแหขึ้นมาใหม่ก็ปรากฏท่อนไม้ท่อนเดิมอีก ต่อจากนั้นแกก็ขว้างท่อนไม้ขึ้นฝั่ง และแกก็ฉุก คิดว่าแปลกที่ท่อนไม้ธรรมดาจะสามารถลอยทวนน้ำได้ คงจะเป็นสิ่งวิเศษ และแกก็ได้นำท่อนไม้นั้นขึ้นฝั่ง และเพ่งมองท่อนไม้นั้นพร้อมกับอธิษฐานว่า หากท่อนไม้นี้มีความศักดิ์สิทธิ์ขอให้คืนนี้จับปลาได้มาก เมื่อพ้นจากความจนแล้ว เมื่อขึ้นฝั่งจะนำท่อนไม้นี้แกะสลักเป็นเทวรูปศักดิ์สิทธิและสักการะบูชาเช้าวันไม่ให้ขาด เมื่ออธิษฐานจบแกเอาท่อนไม้นั้นวางบนหัวเรือ ปรากฏว่าช้อนเพียงสองถึงสามครั้งก็ได้ปลาตัวโตเต็มเรือ จึงนำปลาขึ้นฝั่งวันนั้นปลาของแกขายได้ราคา เพราะชาวประมงคนอื่นจับได้น้อยแกจึงมีเงินจับจ่ายใช้สอย และทุกครั้งที่แกออกหาปลา