The Mystery of the Shape Parameter Lin-Tian Luh Department of Mathematics, Providence University 0 Shalu, Taichung 1 0 email:ltluh@pu.edu.tw 2 phone:(04)26328001 ext. 15126 n fax:(04)26324653 a J January 28, 2010 8 2 ] Abstract. It’s well known that the multiquadrics ( 1) β (c2+ x 2)β, β > 0, and the inverse A − ⌈ ⌉ k k multiquadrics (c2+ x 2)β, β <0, are very useful radialbasis functions for generatingapproximat- N k k ing functions. However the optimal choice of the shape parameter c is unknown. This question has h. beenperplexingmanypeopleformanyyearsandisregardedasoneofthemostimportantquestions t inthetheoryofRBF.Hithertothereareveryfewtheoreticalworksdealingwiththis topic. Instead, a m sporadicresults ofexperimentcanbe seenin the literature. The purposeofthis paper is to uncover its mystery. [ 1 keywords: radial basis function, multiquadric, shape parameter v 7 8 1 Introduction 0 5 Before entering the core of our theory, we would like to make a clarification of the definition. We . 1 are not going to use the simple forms of the multiquadrics and inverse multiquadrics as mentioned 0 in the abstract. Instead, we define 0 1 β v: h(x):=Γ (c2+ x2)β2, β R 2N 0, c>0 (1) −2 | | ∈ \ ≥ i (cid:18) (cid:19) X ,where x istheEuclideannormofxinRn,Γistheclassicalgammafunctionandc,β areconstants. r | | a The function h is called multiquadric or inverse multiquadric, respectively, depending on β > 0 or β <0. The reasonofouradoptingΓisthatitwillmaketheFouriertransformofhandouranalytic work simpler. Ifβ >0,thenh(x)willbeconditionallypositivedefinite[8]oforderm= β where β denotes ⌈2⌉ ⌈2⌉ the smallest integer greater than or equal to β. If β < 0, h(x) is c.p.d.(conditionally positive 2 definite) of order m=0. All these can be found in [8]. For any approximated function f, our interpolating function will be of the form N s(x):= c h(x x )+p(x) (2) i i − i=1 X 1 where p(x) P , the space of polynomials of degree less than or equal to m 1 in Rn, X = m 1 ∈ − − x , ,x is the set of centers(interpolation points). For m = 0, P := 0 . We require that 1 N m 1 s{() i·n·t·erpol}ate f() at data points (x ,f(x )), ,(x ,f(x )). There−fore a{lin}ear system of the 1 1 N N · · ··· form N Q c h(x x )+ b p (x )=f(x ), j =1, ,N (3) i j i i i j j − ··· i=1 i=1 X X and N c p (x )=0 , j =1, ,Q (4) i j i ··· i=1 X , where p , ,p is a basis of P , has to be satisfied. Since h is c.p.d., this requirement will 1 Q m 1 be theor{etica·l·l·y sat}isfied [8]. Howev−er if c is very large, h will be numerically constant, making the linear system (3),(4) numerically unsolvable. Moreover,if c is very large, the coefficient matrix of the linear system will have a very large condition number, making the interpolating function s unreliable when f(x ), ,f(x ) are not accurately evaluated, as pointed out by Madych in [7]. 1 N ··· Our approach for choosing c is based on Theorem2.4 and Corollary2.5 of [4] which we will cite directly but make a slight modification to make them easier to understand. Before introducing the main theorem, we need some fundamental definitions. Let (Rn) denote D the space of complex-valued functions on Rn that are compactly supported and infinitely differen- tiable. For any integer m 1, let ≥ = ϕ (Rn): xαϕ(x)dx=0 for all α 1. n 1 − Definition 1.3 Let n and β be as in (1). The numbers ρ and ∆ are defined as follows. 0 (a) Suppose β 0, ρ=1+ 2 βs+3 and ∆0 = (2m+2+s)(2ρm2m++1+2s)···(2m+3) ⌈2⌉ where m= β . ⌈2⌉ (b) Suppose n 3 β 0, the class of band-limited functions f in L2(Rn) is defined by B = f L2(Rn): fˆ(ξ)=0 if ξ >σ σ { ∈ | | } , where fˆdenotes the Fourier transform of f. Theorem 1.6 Let h be as in (1) with β >0. Any function f in B belongs to and σ h,m C f h m!S(m,n)2−n−1+4βπ−n−14σ1+β4+nec2σc1−β4−n f L2(Rn) (7) k k ≤ k k where c, β are as in (1) apnd S(m,n) is a constant determined by m and n. Proof. We are going to show B by Corollary3.3and Theorem5.2 of [6]. σ h,m ⊆C Letm= β . Thentherequirementa =0forall γ =2minCorollary3.3of[6]isanimmediate ⌈2⌉ γ | | resultof Theorem5.2of [6]. Since B L2(Rn), any f B is a member of where denotes the σ σ ′ ⊆ ∈ S S Schwarzspace. Now, let ρ(ξ) be the Borelmeasurementioned in Corollary3.3of[6]. By (3.9)of[6], it suffices to show that 1/2 m! f := (Dαf)ˆ 2 < k kh  α! k kL2(ρ) ∞ |αX|=m(cid:18) (cid:19)    3 for all f B . We proceed as follows. σ ∈ 1/2 m! (Dαf)ˆ(ξ)2dρ(ξ)  α! | |  |αX|=m ZRn  1/2  m!  = imξαfˆ(ξ)2dρ(ξ)  α! | |  |αX|=m ZRn  1/2  m!  = ξ2α fˆ(ξ)2dρ(ξ)  α! | |  |αX|=m ZRn  1/2  1  1 = (m!)1/2 ξ2α fˆ(ξ)2 dξ (by (3.8) and Theorem5.2 of [6]) |αX|=mα! ZRn | | (2π)2n|ξ|2mhˆ(ξ)  1/2 √m!  1 ξ2α fˆ(ξ)2 1  = | | dξ (by Theorem8.15 of [8]) (2π)n |αX|=mα!ZRn |ξ|2m · 21+β2(|ξc|)−β2−n2Kn+2β(c|ξ|)  ≤ (√2πm)n! ·2−(12+β4)(S(m,n)ZRn|fˆ(ξ)|2· |ξc|nn+2+2ββ · Kn+2β1(c|ξ|)dξ)1/2 where S(m,n) denotes the number of terms in the “ ” X 1/2 m!S(m,n) c−(n4+β) fˆ(ξ)2 ξ n+2β c|ξ| dξ (by Corollary5.12 of [8]) ≤ p(2π)n221+β4 · (ZRn| | | | · pπ2e−c|ξ| ) 1 = p(2mπ)!Sn2(m21+,nβ4) ·c1−(n4+β) rπ2!2 (cid:26)ZRn|fˆ(ξ)p|2|ξ|1+n2+βec|ξ|dξ(cid:27)1/2 m!S(m,n) 2−n−14−β4π−14−nσ1+n4+βec2σc1−(n4+β) f L2(Rn) ≤ · k k < . p ∞ ThusB and(7)follows. ♯ σ h,m ⊆C Ifin(1)β <0,thenorm isdefinedinaslightlydifferentway. Hencewehandleitseparately h k·k and present the following theorem. Theorem 1.7 Let h be as in (1) with β < 0 such that n+β 1 or n+β = 1. Any function f ≥ − in B belongs to and satisfies σ h,m C f h 2−n−1+4βπ−n−14σ1+β4+nec2σc1−(n4+β) f L2(Rn). (8) k k ≤ k k Proof: If β <0, then h is conditionally positive definite of order m =0 [8]. By Theorem5.2 of [6], we know that in Corollary3.3of [6] a =0 for γ =2m. Obviously any f B belongs to since γ σ ′ | | ∈ S f L2(Rn). In order to apply Corollary3.3 of [6] to show that B , it remains to show that σ h,m ∈ ⊆C f B implies fˆ L2(ρ) where dρ(ξ):=r(ξ)dξ = 1 dξ as stated in [6]. ∈ σ ∈ (2π)2nhˆ(ξ) 4 Now, let f B . By (3.9) of [6], it suffices to show that f := fˆ < . We proceed as σ h L2(ρ) ∈ k k k k ∞ follows. 1/2 fˆ = fˆ(ξ)2dρ(ξ) L2(ρ) k k | | (cid:26)ZRn (cid:27) 1/2 = fˆ(ξ)2r(ξ)dξ | | (cid:26)ZRn (cid:27) 1 1/2 = fˆ(ξ)2(2π)−2n dξ where w()=hˆ() by Theorem5.2 of [6] | | · w( ξ) · · (cid:26)ZRn − (cid:27) 1/2 = fˆ(ξ)2(2π)−2n2−1−β2(|ξ|)β+2n 1 dξ by (8.7), p.109 of [8] (ZRn| | c · Kβ+2n(c|ξ|) ) 1/2 2−1−β2−2nπ−2nc−(n+2β) |fˆ(ξ)|2|ξ|n+2β dξ if n+β 1 by Corollary5.12 ≤  ZRn π2 · e√−cc|ξξ|  | |≥ | | of [8] p   1/2 = 2−21−β2−2nπ−2n−12c21−n+2β fˆ(ξ)2 ξ n+2β+1ec|ξ|dξ | | | | (cid:26) ZRn (cid:27) 1/2 2−12−β−2nπ−2n−12c21−n+2βσn+2β+1 fˆ(ξ)2ec|ξ|dξ if n+β+1 0 ≤ | | ≥ (cid:26) ZRn (cid:27) 1/2 2−12−β−2nπ−2n−12σn+2β+1c21−n+2βecσ fˆ(ξ)2dξ ≤ ( Z|ξ|≤σ| | ) 2−14−β−nπ−n−41σn+4β+1c1−(n4+β)ec2σ f L2(Rn) ≤ k k < ♯ ∞ In the preceding deduction we put the restriction n+β 1 on it. This is a drawback because | |≥ the frequently seen case n=1 and β = 1 is not covered. This gap will be discussed shortly. − Corollary 1.8 For any positive number σ and any function f B , under the conditions of Theo- σ ∈ rem1.4, f(x) s(x) m!S(m,n)(2π)−34n√nαnσ1+β4+n ∆0c1+β4−nec2σ(λ)1δ f L2(Rn) (9) | − |≤ k k p p if n+β 1 or n+β = 1. The number S(m,n) is determined by m and n, and is one whenever ≥ − β <0. The other constants are as in Theorem1.4. 2 How to choose c?—a more practical approach In (9) the only things influenced by c are c1+β4−n, ec2σ and (λ)1δ. The number δ is the famous fill distance of the data points, as explained in Corollary2.5 of [4]. The constants β, σ are defined in (1) and Definition1.5 respectively. Theoretically (λ)1δ is very influential. However it contributes little to the error estimate in practice due to two reasons. First, λ is quite large and near one as shown in Theorem1.4. Second, δ can not be very small(near zero). As shown by Madych in 5 [7], whenever δ is too small, the condition number of the linear system (3),(4) will be very large, making its solution meaningless. Note that γ very fast as n . The first examples are n → ∞ → ∞ γ =2, γ =12, γ =78 and γ =632. Thus λ reacts to the change of c very slowly. The number 1 2 3 4 λ will become small enough to make (λ)1δ influential only when c is extremely large. For example, if the dimension n=2 and the domain cube E has side length b =1, then 0 1 2 6Cγ2 λ= 3 (cid:18) (cid:19) where 2 ρ C =max 2ρ′√2e4γ2, 3 , ρ′ = c, γ2 =12. (cid:26) (cid:27) The numberρ isequaltooneorabit largerthanone,depending onβ. ThusC is usuallyverylarge and λ is near one, unless c is very large. Another problem is that the function β h(x):=Γ (c2+ x2)β2 −2 | | (cid:18) (cid:19) willinpracticemakethelinearsystem(3),(4)numericallyunsolvablewhenevercistoolarge. Inour precedingexampleallcenters(interpolationpoints)lieinthecubeofside1. Therefore x x √2. i j | − |≤ If our computer allows seven or fourteen significant digits, then h(xi−xj) will become Γ(−β2)(c2)β2 whenever c 104 or 107. In our example λ>(2)418 and is near one even when c=1012. Therefore ≥ 3 (λ)1δ is only theoretically influential and can be ignored in practice. What’s in practice influential is c1+β4−n and ec2σ. In this section we assume that β, b (the side length of the cube), and δ are fixed and m = 0 max 0, β as in Theorem1.4. The optimal choice or suggested value of c is presented according ⌈2⌉ to thne valuesoof β, n, σ, δ, and b0. Case1 β+1 n 0 Letf B andE be a cube inRn with sidelengthb asinTheorem1.4. Let σ 0 − ≥ ∈ hbe asin(1)withn+β 1orn+β = 1. Let0<δ < b0 befixedinTheorem1.4. Suppose ≥ − 4γn(m+1) β+1 n 0. Then − ≥ (a) the optimal choice of c is to let c = 12ρ√ne2nγnγn(m+1)δ if n = 1 or β = 1, where ρ 6 6 and γ were defined in Definition1.3 and 1.2 respectively; n (b) theoptimalchoiceofc,whenn=1andβ =1,isasfollows. Letc1 =12ρ√ne2nγnγn(m+ 1)δ, c =3b e4, and η = ln23 . Then 0 0 24e4δ (i) choose c=12ρ√ne2nγnγn(m+1)δ if σ2 +η ≥0; (ii) if σ +η <0 and c < 1 3b e4) does not change δ , but 3b0 0 0 0 decreasingc(i.e. c<3b e4) makesδ smaller. After 0<δ < b0 is fixed, c cannotbe less than 0 0 4γn(m+1) 12ρ√ne2nγnγn(m+1)δ because of the restriction δ δ0. ≤ (a) follows from (9) obviously because (λ)1δ can be ignored and the only thing influenced by c is c1+β4−n ec2σ. The smaller c is, the smaller the error bound (9) is. · (b)ismorecomplicated. Althoughwepreviouslymentionedthat(λ)1δ canbeignoredpractically, there is an exception. When n =1 and β = 1, λ is not very large and δ is not required to be very small. Therefore the influence of (λ)1δ has to be taken into consideration. For n=1, β =1, let H(c) := c1+β4−nec2σ(λ)1δ 1 = c41ec2σ 2 12Cδ 3 (cid:18) (cid:19) = ( cc4411eecc22σσ((2233))b2804δec4δ iiff cc≤≥33bb00ee44. H(c) is increasing on [3b e4, ). For c (0, 3b e4], 0 0 ∞ ∈ 1 c H(c) = c41ec2σ 2 24e4δ 3 "(cid:18) (cid:19) # ln2 = c41ec2σeηc where η := 3 <0 24e4δ = c14e(σ2+η)c. Then 1 σ H′(c)=e−43e(σ2+η)c +c +η . 4 2 (cid:20) (cid:16) (cid:17)(cid:21) If σ +η 0, then H (c) > 0 and H is increasing on (0, 3b e4]. In this situation the optimal c is 2 ≥ ′ 0 naturally the minimal possible choice c=12ρ√ne2nγnγn(m+1)δ. Assume σ +η < 0. Then H (c) = 0 iff c = 1 =: c . H (c) < 0 if c > c and H (c) > 0 if 2 ′ 4(σ−+η) ∗ ′ ∗ ′ 2 c0, if f B , then f and σ h,m − ∈ ∈C 1/2 1 1 kfkh ≤(2π)−n2−41 (K0(1)Z|ξ|≤1c |fˆ(ξ)|2dξ+ a0 Z1c<|ξ|≤σ|fˆ(ξ)|2pc|ξ|ec|ξ|dξ) (10) 1 if 1 <σ, where a = √π3−2, and c 0 2Γ(1) 2 1/2 1 f h (2π)−n2−41 fˆ(ξ)2dξ (11) k k ≤ (K0(1)Z|ξ|≤1c | | ) if 1 σ. c ≥ Proof. By (3.9) of [6], f = fˆ and h L2(ρ) k k k k 1/2 fˆ = fˆ(ξ)2dρ(ξ) L2(ρ) k k | | (cid:26)ZRn (cid:27) 8 1/2 = fˆ(ξ)2r(ξ)dξ | | (cid:26)ZRn (cid:27) 1/2 1 = fˆ(ξ)2(2π) 2n dξ where w()=hˆ() by Theorem5.2 of [6] − | | · w( ξ) · · (cid:26)ZRn − (cid:27) 1/2 1 = fˆ(ξ)2(2π)−2n2−12 dξ by (8.7) of [8](p.109) | | · (cξ ) (cid:26)ZRn K0 | | (cid:27) 1/2 fˆ(ξ)2 fˆ(ξ)2 = [(2π)−2n2−12]12 | | dξ+ | | dξ (cξ ) (cξ ) (Z|ξ|≤1c K0 | | Z|ξ|>1c K0 | | ) 1/2 1 1 (2π)−n2−41 fˆ(ξ)2dξ+ fˆ(ξ)2 cξ ec|ξ|dξ ≤ (K0(1)Z|ξ|≤1c | | a0 Z|ξ|>1c | | p | | ) √π by Corollary5.12 of [8] and p.374 of [1] where a = . 0 2Γ(1)√3 2 Since fˆ(ξ)=0 for ξ >σ, our conclusionfollowsimmediately. ♯ | | Lemma 2.2 For any positive number σ and any function f B , under the conditions of Theo- σ ∈ rem1.4 with β = 1 and n=1, − A+B 1/2 |f(x)−s(x)|≤√2π ∆0(λ)1δ c (12) (cid:26) (cid:27) p where A = K01(1) |ξ|≤1c |fˆ(ξ)|2dξ, B = 2√√3Γπ(12) 1c<|ξ|≤σ|fˆ(ξ)|2 c|ξ|ec|ξ|dξ if 1c ≤ σ and B = 0 if 1 σ. c ≥ R R p Proof. This is just an immediate result of Theorem1.4 and Lemma2.1. ♯ With Lemma2.1 and 2.2 we can now begin to analyze the appropriate choice of c. Case3 n=1,β = 1 Let f B and E be a cube in Rn with side length b as in Theorem1.4. σ 0 − ∈ Let h be as in (1) with n = 1, β = 1. If 0 < δ < min 1 , b0 in Theorem1.4 is fixed, under − 24σe4 8 the conditions of Theorem1.4, the suggested choice for c is to let c= 1. (cid:8) (cid:9) σ Reason: Aspointedoutinthe beginningofthissection,the number(λ)1δ in(12)canbe essentially ignored. Hencethe onlythingin(12)influencedbychangingcis A+B. Note thatifc (0, 1], then c ∈ σ B =0. We define H (c)= A if c (0, 1] and H (c)= A+B if c [1, ). Then 1 c ∈ σ 2 c ∈ σ ∞ 1 σ σ 1 H =σA= fˆ(ξ)2dξ = f 2 =H . 1 σ (1) | | (1)k kL2(Rn) 2 σ (cid:18) (cid:19) K0 Z|ξ|≤σ K0 (cid:18) (cid:19) Now, for c (0, 1], ∈ σ 1 lim H (c) = lim fˆ(ξ)2dξ c→0+ 1 c→0+ cK0(1)Z|ξ|≤1c | | 1 = lim f 2 c 0+ c 0(1)k kL2(Rn) → K = . ∞ 9 For c [1, ), ∈ σ ∞ 1 1 2√3Γ(1) lim H (c) = lim fˆ(ξ)2dξ+ 2 fˆ(ξ)2 cξ ecξdξ c→∞ 2 c→∞ c (K0(1)Z|ξ|≤1c | | √π Z1c<|ξ|≤σ| | p | | | | ) lim 1 kfk2L2(Rn) + 2√3Γ(12) √cσecσ f 2 ≤ c→∞ c ( K0(1) √π · k kL2(Rn)) = . ∞ By the structure of the integrands, we know that in fact lim H (c)= , not just . c 2 Therefore, the upper bound of f(x) s(x) in (12) goes →to∞infinity as∞c 0+ or c≤∞ . This | − | → → ∞ phenomenon together with the agreement of H and H at 1 suggests that c = 1 is a safe and 1 2 σ σ natural choice. Remark: (a) In Theorem1.4 we require that δ δ where δ 0+ as c 0+. If δ is fixed, 0 0 ≤ → → c cannot approach0 arbitrarily. Note that for n=1, β = 1, − c if c 3b e4 δ0 = b204e4 if c≤3b0e4. (cid:26) 8 ≥ 0 The requirement δ c is equivalent to 24e4δ c. Therefore our exploration is in fact restricted ≤ 24e4 ≤ in the range c [24e4δ, 3b e4] [3b e4, ) = [24e4δ, ). This is a drawback. However the 0 0 ∈ ∪ ∞ ∞ assumption δ 36√2(210)4 >36√21012 ≥ which will make the linear system (3),(4) numerically unsolvable because c is too large. Since γ very fast as n , the minimum requirement for c will become extremely stringent for n → ∞ → ∞ high dimensions. On the other hand, as is well known by RBF people, ill-conditioning will happen whenever the fill distance is small. Therefore the restriction δ δ is a big problem. Madych’s 0 ≤ experiment [7] shows that when the domain size b =8, andn=1, the conditionnumber is already 0 toolargeifδ =0.125. Ofcourseδ canbeincreasedbyincreasingc. Howeverthelinearsystem(3),(4) 0 will become numerically unsolvable when c is too large. Moreover,the conditionnumber also grows 10

The Mystery of the Shape Parameter

Similar The Mystery of the Shape Parameter