📄 pckkntdm.html
字号:
<html xmlns:mwsh="http://www.mathworks.com/namespace/mcode/v1/syntaxhighlight.dtd"> <head> <meta http-equiv="Content-Type" content="text/html; charset=utf-8"> <!--This HTML is auto-generated from an M-file.To make changes, update the M-file and republish this document. --> <title>How to choose knots when you have to</title> <meta name="generator" content="MATLAB 7.1"> <meta name="date" content="2005-07-27"> <meta name="m-file" content="pckkntdm"> <link rel="stylesheet" type="text/css" href="../../matlab/demos/private/style.css"> </head> <body> <div class="header"> <div class="left"><a href="matlab:edit pckkntdm">Open pckkntdm.m in the Editor</a></div> <div class="right"><a href="matlab:echodemo pckkntdm">Run in the Command Window</a></div> </div> <div class="content"> <h1>How to choose knots when you have to</h1> <introduction> <p>Illustration of the use of OPTKNT and NEWKNT.</p> </introduction> <h2>Contents</h2> <div> <ul> <li><a href="#1">Sample function</a></li> <li><a href="#2">Sample data</a></li> <li><a href="#3">General considerations</a></li> <li><a href="#4">Optimal knots</a></li> <li><a href="#5">Trying OPTKNT</a></li> <li><a href="#7">What happened?</a></li> <li><a href="#8">Knot choice for least squares approximation</a></li> <li><a href="#9">Least squares approximation with uniform knot sequence</a></li> <li><a href="#11">Using NEWKNT to improve the knot distribution for these data</a></li> </ul> </div> <h2>Sample function<a name="1"></a></h2> <p>Here are some sample data, much used for testing spline approximation with variable knots, the so-called Titanium Heat Data, which record some property of titanium measured as a function of temperature. Notice the rather sharp peak. We'll use this function to illustrate some methods for knot selection. </p><pre class="codeinput">[xx,yy] = titanium;frame = [-10 10 -.1 .1]+[min(xx),max(xx),min(yy),max(yy)];plot(xx,yy,<span class="string">'x'</span>), grid <span class="string">off</span>, axis(frame)title(<span class="string">'The titanium heat data show a pronounced peak'</span>)</pre><img vspace="5" hspace="5" src="pckkntdm_01.png"> <h2>Sample data<a name="2"></a></h2> <p>We pick a few data points from these somewhat rough data, since we want to interpolate. The data points picked are marked in the graph. </p><pre class="codeinput">pick = [1 5 11 21 27 29 31 33 35 40 45 49];tau = xx(pick); y = yy(pick);hold <span class="string">on</span>plot(tau,y,<span class="string">'ro'</span>)title(<span class="string">'Titanium heat data with selected points marked'</span>)</pre><img vspace="5" hspace="5" src="pckkntdm_02.png"> <h2>General considerations<a name="3"></a></h2> <p>Since a spline of order k with n+k knots has n degrees of freedom, and we have 12 data sites, tau(1)< ... < tau(12), a fit with a cubic spline, i.e., a fourth order spline, requires a knot sequence t of length 12+4 . </p> <p>Moreover, the knot sequence t must satisfy the Schoenberg-Whitney conditions, i.e., must be such that the i-th data site lies in the support of the i-th B-spline, i.e., </p><pre> t(i) < tau(i) < t(i+k) , all i ,</pre><p>(with equality allowed only in case of a knot of multiplicity k ).</p> <p>One way to choose a knot sequence satisfying all these conditions is as the optimal knots, of Gaffney/Powell and Micchelli/Rivlin/Winograd.</p> <h2>Optimal knots<a name="4"></a></h2> <p>In o p t i m a l spline interpolation, to values at sites tau(1), ..., tau(n) say, one chooses the knots so as to minimize the constant in a standard error formula. Specifically, one chooses the first and the last data site as a k-fold knot. The remaining n-k knots are supplied by OPTKNT. </p> <p>Here is the beginning of the help from OPTKNT:</p><pre class="codeinput"><span class="comment">%OPTKNT Optimal knot distribution.</span><span class="comment">%</span><span class="comment">% OPTKNT(TAU,K) returns an `optimal' knot sequence for interpolation at</span><span class="comment">% data sites TAU(1), ..., TAU(n) by splines of order K.</span><span class="comment">% TAU must be an increasing sequence, but this is not checked.</span><span class="comment">%</span><span class="comment">% OPTKNT(TAU,K,MAXITER) specifies the number MAXITER of iterations to be</span><span class="comment">% tried, the default being 10.</span><span class="comment">%</span><span class="comment">% The interior knots of this knot sequence are the n-K sign-changes in</span><span class="comment">% any absolutely constant function h ~= 0 that satisfies</span><span class="comment">%</span><span class="comment">% integral{ f(x)h(x) : TAU(1) < x < TAU(n) } = 0</span><span class="comment">%</span><span class="comment">% for all splines f of order K with knot sequence TAU.</span></pre><h2>Trying OPTKNT<a name="5"></a></h2> <p>We try this for interpolation on our example, in which we want to interpolate by cubic splines to data (tau(i),y(i)), i=1,...,n : </p><pre class="codeinput">k = 4;osp = spapi( optknt(tau,k), tau,y);fnplt(osp,<span class="string">'r'</span>), hold <span class="string">on</span>, plot(tau,y,<span class="string">'ro'</span>), grid <span class="string">off</span>title(<span class="string">'Spline interpolant to the selected data using optimal knots'</span>)xlabel(<span class="string">'... a bit disconcerting'</span>)</pre><img vspace="5" hspace="5" src="pckkntdm_03.png"> <p>This is a bit disconcerting!</p> <p>Here, marked by stars, are also the (interior) optimal knots:</p><pre class="codeinput">xi = fnbrk(osp,<span class="string">'knots'</span>); xi([1:k end+1-(1:k)]) = [];plot(xi,repmat((frame(3)+frame(4))/2, size(xi)),<span class="string">'*'</span>)xlabel(<span class="string">'... and the optimal knots (*)'</span>)</pre><img vspace="5" hspace="5" src="pckkntdm_04.png"> <h2>What happened?<a name="7"></a></h2> <p>The knot choice for optimal interpolation is designed to make the maximum over a l l functions f of the ratio norm{f - If}/norm{D^k f} of the norm of the interpolation error f - If to the norm of the k-th derivative D^k f of the interpoland as small as possible. Since our data imply that D^k f is rather large, the interpolation error near the flat part of the data is of acceptable size for such an `optimal' scheme. </p> <p>Actually, for these data, the ordinary cubic spline interpolant provided by CSAPI does quite well:</p><pre class="codeinput">cs = csapi(tau,y); fnplt(cs,<span class="string">'g'</span>,2)xlabel(<span class="string">'... and the optimal knots (x), and the interpolant provided by CSAPI.'</span>)hold <span class="string">off</span></pre><img vspace="5" hspace="5" src="pckkntdm_05.png"> <h2>Knot choice for least squares approximation<a name="8"></a></h2> <p>Knots must be selected when doing least-squares approximation by splines. One approach is to use equally spaced knots to begin with, then use NEWKNT with the approximation obtained for a better knot distribution. </p> <p>The next sections illustrate this with the full titanium heat data set.</p><pre class="codeinput">cla, title(<span class="string">''</span>), axis <span class="string">off</span></pre><img vspace="5" hspace="5" src="pckkntdm_06.png"> <h2>Least squares approximation with uniform knot sequence<a name="9"></a></h2> <p>We start with a uniform knot sequence.</p><pre class="codeinput">sp = spap2(augknt(linspace(xx(1), xx(end), 2+fix(length(xx)/4)),k),k, xx,yy);fnplt(sp,<span class="string">'r'</span>), hold <span class="string">on</span>, plot(xx,yy,<span class="string">'x'</span>), axis(frame)title(<span class="string">'Least-squares fit by cubic spline (red) with uniform knot sequence,'</span>)</pre><img vspace="5" hspace="5" src="pckkntdm_07.png"> <p>This is not at all satisfactory. So we use NEWKNT for a spline approximation of the same order and with the same number of polynomial pieces, but better distributed: </p> <h2>Using NEWKNT to improve the knot distribution for these data<a name="11"></a></h2><pre class="codeinput">spgood = spap2( newknt(sp), k, xx,yy);fnplt(spgood,<span class="string">'g'</span>,1.5)xlabel(<span class="string">'... and (green) with a knot sequence chosen from it by NEWKNT.'</span>)hold <span class="string">off</span></pre><img vspace="5" hspace="5" src="pckkntdm_08.png"> <p>This is quite good. Incidentally, one less interior knot would not have sufficed in this case.</p> <p class="footer">Copyright 1987-2005 C. de Boor and The MathWorks, Inc.<br> Published with MATLAB® 7.1<br></p> </div> <!--##### SOURCE BEGIN #####%% How to choose knots when you have to%% Illustration of the use of OPTKNT and NEWKNT.% Copyright 1987-2005 C. de Boor and The MathWorks, Inc.% $Revision: 1.18.4.2 $%% Sample function% Here are some sample data, much used for testing spline approximation% with variable knots, the so-called Titanium Heat Data, which record some% property of titanium measured as a function of temperature. Notice the% rather sharp peak.% We'll use this function to illustrate some methods for knot selection.[xx,yy] = titanium;frame = [-10 10 -.1 .1]+[min(xx),max(xx),min(yy),max(yy)];plot(xx,yy,'x'), grid off, axis(frame)title('The titanium heat data show a pronounced peak')%% Sample data% We pick a few data points from these somewhat rough data, since we% want to interpolate. The data points picked are marked in the graph.pick = [1 5 11 21 27 29 31 33 35 40 45 49];tau = xx(pick); y = yy(pick);hold onplot(tau,y,'ro')title('Titanium heat data with selected points marked')%% General considerations% Since a spline of order k with n+k knots has n degrees of freedom, and% we have 12 data sites, tau(1)< ... < tau(12), a fit with a cubic spline,% i.e., a fourth order spline, requires a knot sequence t of length 12+4 .%% Moreover, the knot sequence t must satisfy the Schoenberg-Whitney% conditions, i.e., must be such that the i-th data site lies in the support% of the i-th B-spline, i.e.,%% t(i) < tau(i) < t(i+k) , all i ,%% (with equality allowed only in case of a knot of multiplicity k ).%% One way to choose a knot sequence satisfying all these conditions is% as the optimal knots, of Gaffney/Powell and Micchelli/Rivlin/Winograd.%% Optimal knots% In o p t i m a l spline interpolation, to values at sites% tau(1), ..., tau(n)% say, one chooses the knots so as to minimize the constant in a standard% error formula. Specifically, one chooses the first and the last data site% as a k-fold knot. The remaining n-k knots are supplied by OPTKNT.%% Here is the beginning of the help from OPTKNT:%OPTKNT Optimal knot distribution.%% OPTKNT(TAU,K) returns an `optimal' knot sequence for interpolation at% data sites TAU(1), ..., TAU(n) by splines of order K.% TAU must be an increasing sequence, but this is not checked.%% OPTKNT(TAU,K,MAXITER) specifies the number MAXITER of iterations to be% tried, the default being 10.%% The interior knots of this knot sequence are the n-K sign-changes in% any absolutely constant function h ~= 0 that satisfies%% integral{ f(x)h(x) : TAU(1) < x < TAU(n) } = 0%% for all splines f of order K with knot sequence TAU.%% Trying OPTKNT% We try this for interpolation on our example, in which we want to% interpolate by cubic splines to data (tau(i),y(i)), i=1,...,n :k = 4;osp = spapi( optknt(tau,k), tau,y);fnplt(osp,'r'), hold on, plot(tau,y,'ro'), grid offtitle('Spline interpolant to the selected data using optimal knots')xlabel('... a bit disconcerting')%%% This is a bit disconcerting!%% Here, marked by stars, are also the (interior) optimal knots:xi = fnbrk(osp,'knots'); xi([1:k end+1-(1:k)]) = [];plot(xi,repmat((frame(3)+frame(4))/2, size(xi)),'*')xlabel('... and the optimal knots (*)')%% What happened?% The knot choice for optimal interpolation is designed to make the% maximum over a l l functions f of the ratio norm{f - If}/norm{D^k f}% of the norm of the interpolation error f - If to the norm of the k-th% derivative D^k f of the interpoland as small as possible. Since our data% imply that D^k f is rather large, the interpolation error near the flat% part of the data is of acceptable size for such an `optimal' scheme.%% Actually, for these data, the ordinary cubic spline interpolant provided by% CSAPI does quite well:cs = csapi(tau,y); fnplt(cs,'g',2)xlabel('... and the optimal knots (x), and the interpolant provided by CSAPI.')hold off%% Knot choice for least squares approximation% Knots must be selected when doing least-squares approximation by splines.% One approach is to use equally spaced knots to begin with, then use NEWKNT% with the approximation obtained for a better knot distribution.%% The next sections illustrate this with the full titanium heat data set.cla, title(''), axis off%% Least squares approximation with uniform knot sequence% We start with a uniform knot sequence.sp = spap2(augknt(linspace(xx(1), xx(end), 2+fix(length(xx)/4)),k),k, xx,yy);fnplt(sp,'r'), hold on, plot(xx,yy,'x'), axis(frame)title('Least-squares fit by cubic spline (red) with uniform knot sequence,')%%% This is not at all satisfactory. So we use NEWKNT for a spline approximation% of the same order and with the same number of polynomial pieces, but better% distributed:%% Using NEWKNT to improve the knot distribution for these dataspgood = spap2( newknt(sp), k, xx,yy);fnplt(spgood,'g',1.5)xlabel('... and (green) with a knot sequence chosen from it by NEWKNT.')hold off%%% This is quite good. Incidentally, one less interior knot would not have% sufficed in this case.displayEndOfDemoMessage(mfilename)##### SOURCE END #####--> </body></html>⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -