Trong nhiều lĩnh vực ứng dụng thông thường ta phải đối mặt với nhiệm vụ mô tả các số liệu nhận được ( đo được chẳng hạn) bằng hàm số giải tích. Có hai cách tiếp cận với vấn đề này Khi nội suy các số liệu được coi là chính xác và điều mong muốn là tìm cách mô tả những gì xảy ra giữa các điểm số liệu. Cách tiếp cận này được xét ở phần sau. Trong mục này sẽ xét phương pháp vẽ đường cong còn gọi là phương pháp hồi quy, tức là tìm một đường cong trơn nhẵn phù hợp tối đa ( khớp nhất ) với các dữ liệu nhưng không nhất thiết phải đi qua mọi điểm dữ liệu. Hình dưới đây minh hoạ cho cả hai cách tiếp cận. Điểm dấu ‘o’ là các điểm dữ liệu; đường liền nét nối chúng miêu tả cách nội suy tuyến tính, còn đường chấm chấm là đường khớp nhất với dữ liệu.
1. Khớp đường cong là gì?
Khớp đường cong bao gồm việc trả lời hai câu hỏi cơ bản : Thế nào là khớp nhất ? và có thể dùng loại đường cong nào ?. “ Khớp nhất “ có thể hiểu theo 7 rất nhiều cách và cũng có vô số các loại đường cong. Vì vậy chúng ta sẽ phải làm như thế nào ? Có thể hiểu “khớp nhất” là khi tổng sai số bình phương tại các điểm số liệu được tối thiểu hoá và sử dụng đường cong được giới hạn bởi đa thức thì việc khớp đường cong là dễ dàng và cho kết quả tốt.
Về mặt toán học điều này được gọi là làm khớp đường cong bình phương tối thiểu với đa thức. Nếu cách phát biểu làm cho bạn bối rối thì hãy xem lại hình trên. Khoảng cách thẳng đứng giữa đường cong chấm chấm và các điểm số liệu là sai số tại điểm này. Lấy bình phương khoảng cách này tại tất cả các điểm số liệu và cộng chúng lại với nhau ta có tổng sai số bình phương. Đường cong chấm chấm là đường làm cho tổng này bé nhất, tức là “khớp nhất”. Thuật ngữ “ bình phương tối thiểu “ là cách gọi tắt của “ tối thiểu hoá tổng sai số bình phương “.
Trong MatLab hàm polyfit giải quyết vấn đề khớp đường cong bình phương tối thiểu. Để minh hoạ cho cách dùng hàm này chúng ta bắt đầu bằng số liệu sau:
>> x=[0 .1 .2 .3 .4 .5 .6 .7 .8 .9 1];
>> y=[-.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2];Để sử dụng polyfit ta cần cấp cho nó các số liệu trên và bậc của đa thức mà ta muốn nó khớp nhất với dữ liệu.
Nếu ta chọn n=1 là bậc của đa thức thì ta sẽ được đường xấp xỉ là đường thẳng và được gọi là hồi quy tuyến tính. Mặt khác nếu n=2 ta sẽ có xấp xỉ là đa thức bậc 2. Ta chọn đa thức này
>> n=2;
>> p= polyfit(x,y,n)
p= -9.8108 20.1293 -0.0317Đầu ra của polyfit là véc-tơ dòng chứa các hệ số đa thức nên ta có lời giải là:
y=-9.8108x3+20.1293x-0.0317.Để so sánh lời giải khớp nhất với các điểm dữ liệu, ta sẽ vẽ cả hai : Đầu tiên là bước chuẩn bị số liệu :
>> xi=linspace(0,1,100);%tạo dữ liệu cho trục x để vẽ đa thức
>> z=polyval(p,xi) % tính đa thức p tại các điểm xả.
và sau đó bằng lệnh :
>> plot(x,y,'o',x,y,xi,z,':')ta sẽ ghi các số liệu gốc x và y bằng điểm dấu ‘o’, rồi vẽ chúng bằng các đường thẳng nối giữa các điểm số liệu và cuối cùng là vẽ số liệu của đa thức xấp xỉ xi và z bằng đường chấm chấm.
>> xlabel('x'), ylabel('y=f(x)')
>> title('Xap xi bang duong cong bac 2’)đặt nhãn cho đồ thị. Kết quả là ta sẽ nhận được lại đô thị ở đầu mục này.
Việc chọn bậc đa thức thường là hơi tuỳ ý. Cần phải có 2 điểm để xác định đường thẳng hay đa thức bậc nhất. Mất 3 điểm để xác định đa thức bậc 2. Theo tiến trình này cần phải có n+1 điểm dữ liệu để xác định đơn trị đa thức bậc n. Như vậy trong trường hợp trên đây khi có 11 điểm chúng ta có thể chọn tới đa thức bậc 10. Tuy nhiên do tính chất nghèo nàn về mặt số của các đa thức bậc cao nên không nên chọn đa thức bậc cao hơn mức cần thiết. Ngoài ra khi số bậc đa thức tăng lên, đường xấp xỉ trở nên gồ ghề hơn vì các đa thức bậc cao có thể bị vi phân nhiều lần trước khi chúng trở thành không. Ví dụ chọn đa thức bậc 10 sau :
>> pp= polyfit(x,y,10)
>> format short e % thay đổi dạng biểu diễn số
>> pp.'
% cho hiện các hệ số đa thức dưới dạng cột
ans =
-4.6436e+005
2.2965e+006
-4.8773e+006
5.8233e+006
-4.2948e+006
2.0211e+006
-6.0322e+005
1.0896e+005
-1.0626e+004
4.3599e+002
-4.4700e-001Chú ý là kích thước của các hệ số đa thức trong trường hợp này so với trường hợp bậc 2 trước đây là khớp nhau. Để ý thêm về sự khác nhau bảy bậc giữa hệ số nhỏ nhất (-4.4700e-01) và hệ số lớn nhất (5.8233e+06). Vẽ đồ thị của đa thức này và so với số liệu gốc cùng đường cong bậc 2 xem chúng có khớp nhau không ?
>>zz=polyval(pp,xi);% tính giá trị đa thức bậc 10
>> plot(x,y,'o',xi,z,':',xi,zz)
>> xlabel('x'), ylabel('y=f(x)')
>> title("Xap xi bang duong cong bac 2 và duong cong bac 10')
Trong đô thị, dữ liệu gốc được đánh dấu ‘o, đường cong bậc 2 là đường chấm chấm, còn đường bậc 10 là đường liền nét. Hãy để ý đến các gợn sóng xuất hiện giữa các điểm dữ liệu tại bên phải và bên trái của đường cong bậc 10. Dựa trên đồ thị này ta thấy rõ là triết lý càng nhiều càng tốt không được áp dụng ở đây.
2. Nội suy bậc 1
Mời các bạn đón xem tiếp phần 2.