Classes of elementary function solutions to the CEV model. I
aa r X i v : . [ q -f i n . M F ] A p r Classes of elementary function solutionsto the CEV model. I.
Evangelos Melas a, a Department of Economics, Unit of Mathematics and Informatics,Athens, Greece
Abstract
In the equity markets the stock price volatility increases as thestock price declines. The classical Black − Scholes − Merton (BSM) op-tion pricing model does not reconcile with this association. Cox in-troduced the constant elasticity of variance (CEV) model in 1975, inorder to capture this inverse relationship between the stock price andits volatility. An important parameter in the model is the parameter β , the elasticity of volatility. The CEV model subsumes some of theprevious option pricing models. For β = 0, β = − /
2, and β = − − rootmodel of Cox and Ross, and the Bachelier model. Both in the caseof the BSM model and in the case of the CEV model it has becometraditional to begin a discussion of option pricing by starting with thevanilla European calls and puts. The pricing formulas for these finan-cial instruments give concrete information about the pricing of optionsonly after the employment of some intermediate approximation scheme.However, there are simpler solutions to both models than those per-taining to the standard calls and puts. Mathematically, it makes senseto investigate the simpler cases first. Furthermore we do not allow our-selves to be drawn into any rash generalizations or inferences from thevanilla European case by prematurely focusing on those cases and weobtain concrete information for the pricing of options without needingto introduce any intermediate approximation schemes. In the case ofBSM model simpler solutions are the log and power solutions. Thesecontracts, despite the simplicity of their mathematical description, areattracting increasing attention as a trading instrument. Similar simplesolutions have not been studied so far in a systematic fashion for theCEV model. We use Kovacic’s algorithm to derive, for all half − integervalues of β , all solutions “in quadratures” of the CEV ordinary differ-ential equation. These solutions give rise, by separation of variables, [email protected] o simple solutions to the CEV partial differential equation. In partic-ular, when β = ..., − , − , − , − , , , , , ..., we obtain four classesof of denumerably infinite elementary function solutions, when β = − and β = we obtain two classes of of denumerably infinite elementaryfunction solutions, whereas, when β = 0 we find two elementary func-tion solutions. In the derived solutions we have also dispensed withthe unnecessary assumption made in the the BSM model asserting thatthe underlying asset pays no dividends during the life of the option. Approximately 45 years ago Black and Scholes [1, 2], and independently,Merton [3] (see also [4]), by making a number of crucial assumptions, de-veloped the most widely used model (hereinafter referred to as the “BSMmodel”) in the pricing of financial options. The BSM model states that byconstantly adjusting the proportions of stocks and options in a portfolio,the investor can create a riskless hedge portfolio, where all market risks areeliminated. They were led to a PDE, the Black − Scholes − Merton partialdifferential equation (PDE) (hereinafter referred to as the “BSM PDE”),which governs the price of the option over time.Black and Scholes, and subsequently Merton, derived [2, 3] their pricingformula for European calls and puts, by transforming the BSM PDE intothe “heat equation”. The “heat equation” is a well known parabolic PDE,which has been studied extensively by physicists, and describes the evolutionover time of the distribution of heat in a certain region of space under giveninitial and boundary conditions.The pricing formula derived by Black and Scholes, and subsequentlyMerton, contains the Gaussian probability density function (see e.g. [5]),and therefore, some kind of approximation scheme is needed, in order toextract any concrete information from this formula for the pricing of options.The model also assumes that volatility, a measure of the estimation ofthe future variability for the asset underlying the option contract, remainsconstant over the option’s life, which is not the case because volatility fluc-tuates with the level of supply and demand.The constant volatility hypothesis in BSM model often leads to resultswhich are inconsistent with market data. To improve the discrepancy, theConstant Elasticity of Variance model diffusion process (hereinafter referredto as the “CEV model”) was proposed in [6] (in [7] various jump processeswere incorporated into the model) to model the heteroscedasticity and theleverage effect in returns of common stocks. An important parameter in the2odel is the elasticity of volatility β which controls the relationship betweenvolatility and price.The pricing formula under the CEV diffusion for European options wasgiven for β < β > − neutral pricing theory [6, 8].A breakthrough was made in [9] by Schroder who expressed the pricingformula, for all values of β , in terms of the non − central Chi − squared dis-tribution, which facilitates the computations significantly. The CEV modelhas been further investigated in [10, 11, 12].There exists an extensive literature with different approximation schemes(see, for instance, [13], [14], [15], [16], [17], and [18]) for the the efficient eval-uation of the complementary non − central Chi − square distribution function,and only after the application of such an approximation scheme concrete in-formation can be obtained from the pricing formulas of Cox [6] and Emanueland MacBeth [8] for the pricing of European options under CEV diffusion. Both in the case of the BSM model and in the case of the CEV modelit has become rather traditional to begin a discussion of option pricing bystarting with the vanilla European calls and puts. Other instruments aregenerally labeled as “exotic”. However, there are simpler solutions to bothmodels than those pertaining to the standard calls and puts.Mathematically, it makes sense to investigate the simpler cases first. Wecan get a feeling for the behaviour of options without needing to introduceany intermediate approximation schemes. Rather more importantly, we donot allow ourselves to be drawn into any rash generalizations or inferencesfrom the vanilla European case by prematurely focusing on those cases.In the case of BSM model simpler solutions are the log and power so-lutions. These contracts, despite the simplicity of their mathematical de-scription, are certainly not devoid of financial interest [19, 20, 21]. On thecontrary such contracts are attracting increasing attention as a trading in-strument.For instance, we note in passing, that in 1994, Anthony Neuberger, byusing the log solution, introduced [19] a simple, yet a very peculiar prod-uct, the Log contract. This contract would forever change how we look atvolatility and lay the foundation for the introduction of variance swaps inthe late 1990s and the new volatility index, VIX by CBOE in 2003.Though Neuberger may not have realized it then but what he had done3as to initiate the process by which volatility would no longer remain amathematical and an abstract concept but become a tangible asset.Similar simple solutions have not been studied so far in a systematicfashion for the CEV model. This is the task we undertake in this paper. Wefind classes of simpler solutions to the CEV model. In the next subsectionwe define what we mean by simpler solutions and how we obtain them.
By using Kovacic’s algorithm [22], we find classes of simpler solutionsto the CEV model, for all half − integer values of β. Simpler solutions heremeans solutions expressed in terms of elementary functions. Most impor-tantly we do not need any intermediate approximation scheme in order toobtain concrete information from these solutions for the pricing of the fi-nancial instruments they describe under appropriate initial and boundaryconditions. This is in contradistinction to the pricing formulas for Europeancalls and puts, described in Section 3, which give information for the pric-ing of options only after the employment of an appropriate approximationscheme.In particular by using Kovacic’s algorithm [22], we obtain classes of Li-ouvillian solutions to the CEV ODE, and we so derive, with separation ofvariables, solutions to the CEV PDE, in terms of elementary functions, forall half − integer values of β. Kovacic’s algorithm will find all possible Liouvillian solutions (i.e., essen-tially, all solutions in terms of quadratures) of linear second order homoge-neous ODEs with complex rational function coefficients. Hereafter, solvableby quadratures means that we consider a differential equation solved whenwe are left only with evaluating an antiderivative.The Liouvillian solutions we find to the CEV ODE, are given, for allhalf − integer values of β , in section 9, and are essentially products of poly-nomials, truncated confluent hypergeometric functions of the first kind, withexponential functions and powers of the independent variable. These Liou-villian solutions to the CEV ODE induce, by separation of variables, solu-tions in terms of elementary functions, to the CEV PDE.These simpler solutions to the CEV PDE are amenable to analytic ma-nipulation and easy to use, and thus, are to be juxtaposed with the solutionsto the CEV PDE given in terms of power series which are described in sec-tion 3 in the first category of the existing literature on the pricing of optionsunder the CEV model.This paper is organised as follows: In section 2 the bare essentials of the4EV model are given and the CEV model is compared with the BSM model.In section 3 a very brief review of the literature on the pricing of optionsunder the CEV diffusion is given. In section 4 a remark is made on themeaning of the term “closed − form” solutions which is used in the literatureon the pricing of options, either with the BSM model or with the CEVmodel. In section 5 the BSM PDE, and the associated ODE, are derived.In section 6 the CEV PDE, and the associated ODE, are derived. In section7 we explain that Kovacic’s algorithm is an application of Picard − Vessiottheory to linear second order homogeneous ODEs and we trace the originof Picard − Vessiot theory to the Galois theory of polynomials. In section8 Kovacic’s algorithm is applied to the CEV ODE. In section 9 classes ofLiouvillian solutions to the CEV ODE are obtained, and associated solutionsto the CEV PDE are derived, for all half − integer values of β . In section 10we outline prospects for future research. Finally, in the Appendix we givean outline of Kovacic’s algorithm. The ability to build a portfolio of the kind envisaged by Black, Scholesand Merton [1, 2, 3] relies on the assumptions of continuous trading andcontinuous sample paths of the asset price. Other assumptions made in thederivation of the BSM model are: • The asset price follows a lognormal random walk. • The volatility of the underlying asset is constant over the life of theoption contract. • The risk − free interest rate r is constant and known. • The underlying asset pays no dividends during the life of the option. • There are no arbitrage opportunities. • Short selling is allowed (full use of proceeds from the sale is permitted). • Fractional shares of the underlying asset may be traded.The BSM model also applies to any financial instrument the future ofwhich is uncertain at the present time. Basically it has to do with the pricingof options, but anything vaguely connected such as corporate debt is equallygrist for its mill. 5olatility, a measure of how much a stock can be expected to move inthe near − term, is probably the most important single input to any optionpricing model. In BSM model it is assumed that volatility is a constant overtime. This means that the variance of the return is constant over the life ofthe option contract and is known to market participants. While volatility can be relatively constant in very short term, it is neverconstant in longer term. Therefore, a remedy for this shortcoming of Black-Scholes − Merton model is needed. The CEV model introduced by Cox [6] isan example of such a remedy.The CEV spot price model is a one − dimensional diffusion model withthe instantaneous volatility specified to be a power function of the under-lying spot price, σ ( S ) = aS β . It was introduced by Cox as one of the earlyalternatives to the geometric Brownian motion, assumed in the BSM model,to model asset prices (see also [7] where alternative stochastic processes forthe asset price are introduced).The CEV stochastic process is closely related to Bessel processes andis analytically tractable. In fact there are analytic forms of option pricingformulas for the CEV diffusion; for β < β > − central Chi − square distribution function.For β = 0 the CEV model reduces to the constant volatility geometricBrownian motion process employed in the Black, Scholes and Merton model.When β = −
1, the volatility specification is that of Bachelier (the asset pricehas the constant diffusion coefficient, while the logarithm of the asset pricehas the a/S volatility). We mention in passing that Bachelier’s [23] work hadfallen into oblivion, at least in the financial circles, for an extensive period,and was rediscovered by mathematical economists such as Paul Samuelsonin the 1960s. Recently it has become again popular since it assumes thatinterest rates can be negative. For β = − / − root model of Cox and Ross [7].Thus Cox rather than assuming constant volatility, expressed the volatil-ity as a function of the price of the underlying asset. This is precisely themain advantage of the CEV model: The volatility of the underlying as-set is linked to its price level, thus exhibiting an implied volatility smile(or implied volatility skew) that is a convex and monotonically decreasingfunction of exercise price, similar to the volatility smile curves observed in6ractice (see, for example, [24]). The CEV framework is also consistent withthe so − called leverage effect (i.e., the existence of a negative correlation be-tween stock returns and realized stock volatility) as documented for instancein [25]. In [26], MacBeth and Merville, compared the performances of the BSMand CEV models through simulations and real examples. Their resultsshow that the CEV model has a better performance, which underpricesin − the − money call options and overprices out − of − the money call options.In [27] Dias and Nunes showed that a firm that uses the standard geomet-ric Brownian motion assumption is exposed to significant errors of analysiswhich may lead to non − optimal investment and disinvestment decisions.Given that the lognormal assumption with constant volatility does not cap-ture the implied volatility smile effect observed across a wide range of mar-kets and underlying assets Dias and Nunes used instead the CEV diffusionprocess and gave analytical solutions for perpetual American − style call andput options under the CEV diffusion. Their results strongly highlight thecase for moving beyond the simplistic real options models based on the log-normal assumption to more realistic models incorporating volatility smileeffects. The CEV diffusion process has been extensively used to obtain thesolutions of several flnancial option pricing problems. In particular, theCEV call option pricing formula for valuing European options has beeninitially expressed in terms of the standard complementary gamma distri-bution function by Cox [6] for β <
0, and by Emanuel and MacBeth [8] for β >
0. Schroder [9] has subsequently extended the CEV model by express-ing the corresponding formula in terms of the complementary non − centralChi − square distribution function.There exists an extensive literature devoted to the efficient computationof the complementary non − central Chi − square distribution function, withseveral alternative representations available (see, for instance, [9], [13], [14],[15], [16], [17], and [18]).Things are more complicated in the case of exotic options such as look-back and barrier options. Davydov and Linetsky [28] evaluated European7ookback options under the CEV process with a model based on the nu-merical inversion of the Laplace transform of the option price. To evaluatelookback options Linetsky [29] used spectral theory. In [30], Davydov andLinetsky, priced single − barrier and double − barrier options under the CEVdiffusion process, using again spectral theory.In [31], Costabile, used the binomial process to approximate the CEVprocess and priced lookback options. Boyle, Tian and Imai [32] tackledthe same problem using Monte Carlo simulations. In [33], Boyle and Tian,priced single − barrier and double − barrier options under the CEV diffusionprocess in the numerical trinomial lattice framework.Therefore the work which has been done so far on the pricing of optionsunder the CEV diffusion can be conveniently divided into two categories:1. In the first category the pricing formula is given in closed − form [6,8, 28], which, after Schroder’s work [9], is expressed in terms of thecomplementary non − central Chi − square distribution function. Thereexists an extensive literature ( [13], [14], [15], [16], [17], and [18]) de-voted to the efficient computation of the complementary non − centralChi − square distribution function, with several alternative representa-tions and associated approximation schemes available. Pricing (exotic)options by employing spectral theory (see, for instance, [29] and [30])also employs functions represented by power series, and falls also intothis category.2. In the second category the pricing of (exotic) options is attained bydiscrete approximations of the CEV process: By a discrete approxi-mation of the CEV process using the binomial tree method [31], thetrinomial tree method [33], the Monte Carlo method [32].Needless to say the above review of the literature is extremely far frombeing exhaustive and it only highlights two main trends in the current re-search for pricing options under the CEV diffusion. − form” solutions It is appropriate at this point to make a remark regarding the meaningof the term “closed − form” solutions which is used in the literature on thepricing of options, either with the BSM model or with the CEV model.This term does not imply, as one might anticipate, that the solutions athand give answers for the pricing of options without the implementation of8ome kind of approximation scheme. On the contrary, as we explain below,both in the case of the BSM model and in the case of the CEV model wehave to resort to some kind of approximation in order to value the optionsby using the pricing formulas.Black and Scholes found [2] the celebrated formulas for the pricing ofthe vanilla European calls and puts by transforming the BSM PDE into the“heat equation”.This PDE, which characterises the propagation of heat in a continuousmedium, has been extensively studied in physics. Its fundamental solution,subject to appropriate boundary conditions, is the Gaussian density function(e.g. [34], p. 81). Black and Scholes found the formulas for the pricingof options by transforming this solution back to the initial dependent andindependent variables appearing in the BSM PDE.The formulas so derived by Black and Scholes for the pricing of optionsare characterised in the literature as “closed − form” solutions. However, wehave to bear in mind, that in order to extract any concrete information fromthem, we need to compute the integrals N ( d i ) = √ π R d i −∞ e − s ds, where d i are appropriate limits, numerically by quadrature rules such as the Simson’srule or Gaussian rule.Moreover, in the case of the CEV model the pricing formulas whichhave been derived by Cox [6] for β <
0, by Emanuel and MacBeth [8] for β >
0, and by Schroder [9], give answers for the pricing of options only afterthe implementation of some kind of approximation scheme. For instance,as we noted in subsection 3, in the case of Schroder’s pricing formulas,which are expressed in terms of the complementary non − central Chi − squaredistribution function, there is a whole ongoing literature ( [13], [14], [15], [16],[17], and [18], to mention only a very small fraction of it,) of approximationsschemes for its efficient and fast evaluation.These classes of solutions are to be juxtaposed with the classes of so-lutions to the CEV model we derive in this paper when 2 β belongs to theintegers. The classes of solutions we derive in this paper do not involve anyintegrals or functions expressed in terms of series and we do not need anyintermediate approximation schemes in order to use these solutions for thepricing of various financial instruments. Thus the solutions we derive in thispaper complement the solutions of the first category, described in subsection3, for the pricing of options under the CEV diffusion.9 Derivation of the BSM PDE and BSM ODE
Black, Scholes and Merton made the fundamental observation [1, 2, 3]that if one could perfectly hedge an option, then one could price it aswell. The reason being that a perfectly hedged portfolio has no uncertainty,and hence has a risk − free rate of return given by the spot interest rate r.In summary, their derivation goes as follows: Consider a portfolio Πwhich contains one option and − ∆ units of the underlying asset. The valueof the portfolio is Π = V − ∆ · S t , (1)where ∆ is to be determined, S t is the asset value at time t , and V = V ( S, t )is the value of the option.It is assumed that the asset value has lognormal dynamics, i.e., it satisfiesthe SDE dS t = rS t dt + σS t dW t , (2)where dS t = S t + dt − S t is the change in asset value from t to t + dt , σ is thestandard deviation per unit time (volatility) of the underlying asset value,and W t is a Wiener process. Equation (2) says that the percentage changein asset value from t to t + dt is normally distributed with mean µdt andvariance σ dt. We assume that V ∈ C , ( R × [0 , T ]) , so by applying Ito Lemma we havethat the change in the value of the portfolio is given by: d Π = dV − ∆ · dS t = ∂V∂t dt + ∂V∂S dS t + ∂ V∂S ( dS t ) − ∆ · dS t = (cid:18) ∂V∂t + 12 σ S ∂ V∂S (cid:19) dt + (cid:18) ∂V∂S − ∆ (cid:19) dS t . (3)If we choose ∆ = ∂V∂S then the change d Π in the value of the portfolio Πis no longer sensitive to random changes in the value S t of the underlyingasset, and we obtain: d Π = (cid:18) ∂V∂t + 12 σ S ∂ V∂S (cid:19) dt. (4)Thus the choice ∆ = ∂V∂S yields a perfectly hedged portfolio, i.e., a portfoliowith no uncertainty. 10ince there is no uncertainty left in the portfolio, its value, by the prin-ciple of no arbitrage, has to be the same as if being on a bank account withthe risk − free interest rate r. If the return of Π were larger than r then onewould simply take a loan at a risk − free rate r and buy the portfolio Π toobtain a guaranteed profit. Conversely, If the return of Π were smaller than r then one would simply short the portfolio and invest it in the bank.However this is ruled out by the principle of no arbitrage and thereforewe must have: d Π = r Π dt = r ( V − ∆ · S t ) dt. (5)By combining Equations (4) and (5) we arrive at ∂V∂t + 12 σ S ∂ V∂S + rS ∂V∂S − rV = 0 . (6)This is the BSM PDE.For the BSM PDE, which is a backward parabolic equation, we mustspecify final and boundary conditions, for otherwise the PDE does not havea unique solution. For instance for a vanilla European call c ( S, t ) withexercise price E and expiry date T , the final condition is just its payoff at Tc ( S, T ) = max ( S − E, , for all S ≥ . (7)The asset − price boundary conditions are applied at zero asset price, S = 0 , and as S → + ∞ . At S = 0 we have c (0 , t ) = 0 , for all t ≥ . (8)The second boundary condition, as S → + ∞ , reads c ( S, t ) ∼ S − Ee − r ( T − t ) , as S → + ∞ , for all t ≥ . (9)Moreover, for a vanilla European put option p ( S, t ) , the final conditionis the payoff at T p ( S, T ) = max ( E − S, , for all S ≥ . (10)The asset − price boundary conditions are applied again at zero asset price, S = 0 , and as S → + ∞ . At S = 0, assuming that interest rates are constant,we have p (0 , t ) = Ee − rt , for all t ≥ . (11)11s S → + ∞ , the option is unlikely to be exercised and so for t > , we have p ( S, t ) → , as S → + ∞ , for all t ≥ . (12)With these final conditions and asset − price boundary conditions Blackand Scholes [1, 2], and independently, Merton [3] derived their pricing for-mulas for European calls and puts. Their pricing formulas are essentiallythe so called fundamental solution of the heat equation (see e.g. [34], p. 81),which is the Gaussian density function.We note in passing that the mathematical analysis of American optionsis more complicated than that of European options. It is almost always im-possible to find a useful explicit solution to any given free boundary problem,and so a primary aim is to construct efficient and robust numerical methodsfor their computation.In contrast to ODEs there is no unified theory of PDEs. Some equationshave their own theories, while others have no theory at all. The solutionsof a PDE form in general an infinite − dimensional space and there is a widerange of methods to probe this solution space.This also applies to the BSM PDE (6). This is a second − order linearhomogeneous separable backward parabolic PDE. A method to solve a linearhomogeneous separable PDE is the method of separation of variables. Inthis method such a PDE involving n independent variables is converted inton ODEs.In particular the BSM PDE (6) with the separation of variables V ( S, t ) = A ( S ) B ( t ) (13)is reduced to the following pair of ODEs12 σ S d A d S + rS d A d S − ( r + λ ) A = 0 , (14)d B d t + λB = 0 , (15)where λ is an arbitrary constant. We name ODE (14) the BSM ODE. Thegeneral solutions to the ODEs (14) and (15) are given respectively by A ( S ) = αS ξ + ζ σ + βS ξ − ζ σ , (16)where , ξ = σ − r , ζ = p ( σ − r ) + 8 σ ( r + λ ) , (17)and ,B ( t ) = γe − λt , (18)12here α, β, and γ are arbitrary constants.From Equations (13), (16), and (18) the following solution to the BSMPDE (Equation (6)) is readily obtained V ( S, t ) = e − λt ( aS ξ + ζ σ + bS ξ − ζ σ ) , (19)where a, b arbitrary constants and ξ, ζ are specified in Equation (17).It is noteworthy that both the fundamental solution of the heat equationwhich gives rise to the pricing formulas of Black, Scholes and Merton andthe solutions to the CEV model we derive in this paper are obtained byexploiting symmetries, albeit symmetries of different objects; we elaboratemore on this issue in section 7.The fundamental solution of the heat equation is one of the most funda-mental formulas in all of mathematics. It is the famous Gaussian or normaldensity of probability theory. It describes an initial point source of heat ata given position.It can be derived in a number of different ways. Namely, it can be derivedby Fourier transform ( see e.g. [35] p.98, p.99) and it can also be derivedby using the Lie symmetries of the heat equation. In fact by using the Liesymmetries of the heat equation the fundamental solution can be derived intwo distinct ways.One derivation results from the invariance of the heat equation underthe one − parameter Lie group of dilations ( [34] p.91, p.92 and [36] p.119).Another derivation ([36] p.119, p.120) results from a symmetry of the heatequation which cannot be anticipated from basic physical principles and canonly be obtained by applying Lie theory of symmetries to the heat equation.Remarkably the last derivation is the simplest of the three.The simpler log and power solutions of the BSM PDE, which not only arethey not devoid of financial interest but also they are attracting increasingattention as a trading instrument, are derived with separation of variables(see e.g. [35], p. 127, p. 132). No resort to symmetry principles is needed.In order to derive classes of simpler solutions to the CEV model in thispaper we follow an hybrid method: We apply symmetry principles not tothe CEV PDE itself, but to the CEV ODE, which is derived from the CEVPDE by separation of variables and so we obtain classes of simpler solutionsto the CEV ODE first. Then the classes of simpler solutions to the CEVODE easily induce classes of simpler solutions to the CEV PDE.The derivation of the CEV ODE and the CEV PDE is similar to thederivation of the BSM ODE and the BSM PDE, and it is outlined in Sub-section 6.1. 13 Derivation of the CEV PDE and CEV ODE
Cox [6] considered the CEV diffusion process governed by the SDE when q = 0 dS t = ( r − q ) S t dt + αS β +1 t dW t , (20)where q is a dividend yield parameter and the instantaneous volatility isspecified to be a power function of the underlying spot price, σ ( S ) = αS β , α is the volatility scale parameter. An easy calculation shows that β is theelasticity of the volatility σ ( S ).The elasticity parameter β is the central feature of the model, and con-trols the relationship between the volatility and price of the underlying asset.The dividend yield parameter q has been introduced to dispense with theunnecessary assumption made in the BSM model asserting that the under-lying asset pays no dividends during the life of the option. dS t , r, W t remainas defined in Section 5.The CEV diffusion model seems to be a natural extension of the BSMmodel. In fact as it is pointed out in Subsection 6 the CEV model subsumessome of the previous option pricing models. For β = 0, β = − /
2, and β = − − root model [7] of Cox and Ross, and the Bachelier model [23].In Subsection 2.2 is made clear that the CEV model performs betterthan the BSM model [26, 27]. Moreover, in [37], Beckers considered theCEV model and its implications for option pricing on the basis of empiricalstudies and concluded that the CEV class could be a better descriptor ofthe actual stock price than the traditionally used lognormal model. The derivation of the CEV ODE and the CEV PDE is similar to thederivation of the BSM ODE and the BSM PDE. We give it in the detailneeded in order to clarify the similarities and the dissimilarities between thetwo derivations.Consider a portfolio Π which contains one option and − ∆ units of theunderlying asset. The value of the portfolio isΠ = V − ∆ · S t , (21)where ∆ is to be determined, S t is the asset value at time t , and V = V ( S, t )is the value of the option. It is assumed that the asset value satisfies theSDE (Equation (20)) dS t = ( r − q ) S t dt + αS β +1 t .
14e assume that V ∈ C , ( R × [0 , T ]) , so by applying Ito Lemma we havethat the change in the value of the portfolio is given by: d Π = dV − ∆ · dS t = ∂V∂t dt + ∂V∂S dS t + ∂ V∂S ( dS t ) − ∆ · dS t = (cid:18) ∂V∂t + 12 α S β +2 ∂ V∂S (cid:19) dt + (cid:18) ∂V∂S − ∆ (cid:19) dS t . (22)The choice ∆ = ∂V∂S yields a perfectly hedged portfolio, i.e., a portfoliowith no uncertainty. With this choice the change d Π in the value of theportfolio Π is given by: d Π = (cid:18) ∂V∂t + 12 α S β +2 ∂ V∂S (cid:19) dt. (23)Since there is no uncertainty left in the portfolio, the change d Π in itsvalue Π has to be equal to d Π = ( r − q )Π dt = ( r − q )( V − ∆ · S t ) dt, (24)where r is the risk − free interest rate and q is a dividend yield parameter.By combining Equations (23) and (24) we arrive at ∂V∂t + 12 α S β +2 ∂ V∂S + ( r − q ) S ∂V∂S − ( r − q ) V = 0 . (25)This is the CEV PDE.For the CEV PDE, which is a backward parabolic equation, we mustspecify final and boundary conditions, for otherwise the PDE does not havea unique solution. For instance for vanilla European calls and puts withexercise price E and expiry date T the final and boundary conditions arethose specified in Section 5.A method to solve a linear homogeneous separable PDE is the method ofseparation of variables. In particular the CEV PDE (25) with the separationof variables V ( S, t ) = C ( S ) D ( t ) (26)is reduced to the following pair of ODEs12 α S β +2 d C d S + ( r − q ) S d C d S − ( r − q + λ ) C = 0 , (27)d D d t + λD = 0 , (28)15here λ is an arbitrary constant. We name ODE (27) the CEV ODE.The general solution to the ODE (28) is given by D ( t ) = δe − λt , (29)where δ is an arbitrary constant. The general solution to the ODE (27) isgiven for example in [28] and it reads C ( S ) = εC ( S ) + ζC ( S ) , (30)where ε, ζ are arbitrary constants and C ( S ) , C ( S ) , are given respectivelyby C ( S ) = S β + e ǫ x M k,m ( x ), β < S β + e ǫ x W k,m ( x ), β > C ( S ) = S β + e ǫ x W k,m ( x ), β < S β + e ǫ x M k,m ( x ), β > M k,m ( x ) and W k,m ( x ) are the Whittaker functions (see e.g. [38],Chapter 13, p.505), and, x = | r − q | α | β | S − β , ǫ = sign (( r − q ) β ) , (33) k = ǫ (cid:18)
12 + 14 | β | (cid:19) − r − q + λ | ( r − q ) β | , m = 14 | β | . (34)Whittaker’s functions M k,m ( x ) and W k,m ( x ) are defined in terms of theconfluent hypergeometric functions of the first and the second kind F and U respectively (see e.g. [38], Chapter 13, p.504 − M k,m ( x ) = e − x x m + F (cid:18) m − k + 12 , m ; x (cid:19) , (35) W k,m ( x ) = e − x x m + U (cid:18) m − k + 12 , m ; x (cid:19) . (36)From Equations (26), (29), (30), (31)and (32) the following solution tothe CEV PDE (Equation (25)) is readily obtained V ( S, t ) = e − λt ( cC ( S ) + dC ( S )) , (37)16here c, d are arbitrary constants.We note that when β = 0 the solution space of the CEV ODE (Equation(27)) is spanned by the Liouvillian solutions S ξ ζ σ and S ξ − ζ σ , where , (38) ξ = σ − r − q ) , ζ = p ( σ − r − q )) + 8 σ ( r − q + λ ) . (39)Therefore, in this case the general solution of the CEV ODE is Liouvillianand it is given by A ( S ) = eS ξ ζ σ + f S ξ − ζ σ , (40)where e, f are arbitrary constants. The solution A ( S ), as expected, isnothing but the solution A ( S ) (Equation (16)) where r has been replacedby r − q. From Equations (26), (29), and (40), we conclude that when β = 0 thefollowing solution to the CEV PDE (Equation (25)) is readily obtained V ( S, t ) = e − λt ( gS ξ ζ σ + hS ξ − ζ σ , (41)where g, h are arbitrary constants. To the best of our knowledge Kovacic’s algorithm has not been applied sofar to problems appearing in the field of Economics (the only reference whereit appears that the authors apply Kovacic’s algorithm to problems appearingin the field of Economics is [39], but in this reference the authors do notuse Kovacic’s algorithm to derive new solutions to economic models butrather they comment on the possible applications of Kovacic’s algorithm tocertain financial models). This is to be juxtaposed with with the numerousapplications of Lie theory of continuous symmetries to problems arising inthe field of Economics (see e.g. [40] and references therein).Kovacic’s algorithm is the outcome of Picard − Vessiot theory when thisis applied to linear second order homogeneous ODEs with rational functioncoefficients. Interestingly enough Picard − Vessiot theory and Lie theory havesprung from the same origin, Galois theory of polynomials. In subsection7.1 we highlight the main results of these theories and how these are relatedto the main result of their origin, Galois theory of polynomials.By so doing we bring to the fore one of the major connections betweendifferential equations and abstract algebra − one that is not commonly em-phasized. This connection is related to the solvability of the equation. In17he case of abstract algebra, we have the result that the polynomial equa-tion is solvable by radicals if and only if the Galois group of the polynomialis a solvable group. The differential equation result is that if the equationadmits a solvable group, then it is solvable by quadratures. − Vessiot theory
The idea behind the Galois theory of polynomials [41, 42, 43] is to associateto a polynomial a group, the Galois group, which is the group of symmetriesof the roots that preserve all the algebraic relations among these roots, anddeduce properties of the roots from properties of this group.For instance, a solvable Galois group implies that the roots can be ex-pressed in terms of radicals, and the fact that the symmetric groups S n for n ≥ p ( x ) = 0 be an irreducible polynomial equation with coefficients ina field F . Every element of the Galois group of p ( x ) permutes the roots of p ( x ) and leaves invariant all the algebraic relations satisfied by the roots.This allows to describe the Galois group of p ( x ) as a group of automor-phisms of a field extension: Let A be the group of automorphisms of thefield extension of F which is formed from F by adjoining to it the roots ofthe polynomial equation p ( x ) = 0 . Galois group is the subgroup of A whichleaves fixed pointwise the elements of F .Galois theory motivated two major developments in the theory of dif-ferential equations: Lie group theory [44, 45] for continuous symmetries ofdifferential equations, and, Picard − Vessiot theory [46, 47, 48] of the dif-ferential field extensions generated by the solutions of a linear differentialequation.Picard − Vessiot theory is closer in spirit to Galois theory than Lie theory.In Picard − Vessiot theory one replaces fields by differential fields: fields witha derivation D . Just as adjoining a root of a polynomial equation to a fieldgives an extension of fields, adjoining a root of a differential equation to afield gives an extension of differential fields.As in Galois theory, one can form the differential Galois group of anirreducible linear homogeneous ordinary differential equation E with coef-ficients in a field F ( x ). The differential Galois group of E is the group ofautomorphisms of the extension of the differential field F ( x ) which leavesthe elements of F ( x ) fixed.A main result of Picard − Vessiot theory is that a linear homogeneousdifferential equation can be solved by quadratures if and only if its differen-18ial Galois group is solvable. Picard − Vessiot theory falls into the realm ofdifferential algebra.Lie’s motivation in constructing his theory was the desire to extend thetheory developed by Galois (and Abel [49]) to differential equations. Ifthe discrete invariance group of an algebraic equation could be exploited togenerate algorithms to solve the algebraic equation by radicals, might it bepossible that the continuous invariance group of a differential equation couldbe exploited to solve the differential equation by quadratures?In fact Lie showed in 1894 [45] that as an algebraic equation of degree nis solvable by radicals if its Galois group is an n th − order solvable group, an n th − order ordinary differential equations can be integrated by quadraturesif it admits a solvable n − parameter symmetry group.Two are the key ideas which lie at the foundations of Lie theory andenable its development: The first key idea is Lie’s great advance to replacethe complicated, nonlinear invariance condition of an analytic function undera one − parameter Lie group of transformations by a vastly more useful linearinfinitesimal condition and to recognize that if an analytic function satisfiesthe infinitesimal condition then it also satisfies the finite condition and viceversa.The second key idea in Lie theory is his somewhat unique view of differen-tial equations: ODEs and PDEs can be viewed as locally analytic functionsin a space whose coordinates are independent variables, dependent variables,and the various derivatives of one with respect to the other.This implies in particular that we can derive one − parameter Lie groupswhich leave ODEs or PDEs invariant by applying the first key idea. It turnsout (see e.g [50] p. 129) that this can be achieved in a two − step process. Inthe first step we derive the vector fields whose integral curves are the orbitsof the one − parameter Lie groups which leave the ODEs or PDEs invariant.In the second step we use these vector fields to derive the one − parameterLie groups which leave the ODEs or PDEs invariant.This makes apparent that Lie theory falls into the realm of local differen-tial geometry. Sophus Lie was at heart a geometer, and it was through thislens that he viewed much of his work. On the other hand Picard − Vessiottheory falls into the realm of differential algebra. Lie theory applies to anydifferential equation, whereas Picard − Vessiot theory applies only to linearhomogeneous ODEs of n th order.One expects that in the case of linear homogeneous ODEs of n th orderthe two theories should be related to each other since a main result in boththeories is that if a linear ODE of n th order admits a solvable group, thenit is solvable by quadratures. However the links between Lie theory and19icard − Vessiot theory remained hidden for a long time, mostly becauseof the apparent walls that separate the mathematical disciplines of localdifferential geometry and differential algebra.In fact evidence was given [51] to support the common wisdom thatthe two theories are not related to each other. It came as a surprise whenIbragimov found a bridge [52] between Lie symmetries and Galois groups:He constructed the Galois groups for several simple algebraic equations byfirst calculating their Lie symmetries and then restricting the symmetrygroup to the roots of the equation in question. Thereafter more papershave appeared [53] which study the interplay and connections between thedifferential Galois group and Lie symmetries of linear homogeneous ODEsof n th order.We have emphasized that a main result in Picard − Vessiot theory is thatif a linear homogeneous ODE of n th order admits a solvable group, then it issolvable by quadratures. In fact Picard and Vessiot proceeded even further[46, 47, 48] and stated sufficient and necessary conditions for the existenceof Liouvillian solutions to a linear homogeneous ODE of n th order.Roughly speaking Liouvillian solutions are solutions in quadratures whichcan be expressed in terms of exponentials, integrals and algebraic functions;a more precise definition of Liouvillian solutions is given in subsection 7.2.A formal modern proof of the criteria given by Picard and Vessiot for theexistence of Liouvillian solutions was given by Kolchin [54, 55]. These re-sults are such that lead to several algorithms [56, 57] to decide if a linearhomogeneous ODE of n th order has a Liouvillian solution.Picard − Vessiot − Kolchin theory and the ensuing algorithms to decide ifa linear homogeneous ODE of n th order admits Liouvillian solutions canbe simplified for 2 nd order linear homogeneous ODEs because of severalfacts that are summarized in [58] (Chapter 4.3.4). The resulting algorithmis essentially the algorithm presented by Kovacic in [22]. Kovacic’s algo-rithm [22] predated and motivated much of the work on Liouvillian solu-tions of general linear homogeneous ODEs of n th order. A beautiful accountof Picard − Vessiot − Kolchin theory and of Kovacic’s algorithm is given inSinger’s lectures [59].
Kovacic’s algorithm [22] finds a “closed − form” solution of the differentialequation y ′′ + ay ′ + by = 0 (42)20here a and b are rational functions of a complex variable x , provided a“closed − form” solution exists. The algorithm is so arranged that if no solu-tion is found, then no solution can exist. The “closed − form” solution meansa Liouvillian solution, i.e. one that can be expressed in terms of algebraicfunctions, exponentials and indefinite integrals. (As functions of a complexvariable are considered, trigonometric functions need not be mentioned ex-plicitly, as they can be written in terms of exponentials. Logarithms areindefinite integrals and hence are allowed).In more concrete terms, a Liouvillian function, is a function of one com-plex variable, which is the composition of a finite number of arithmeticoperations (+ , − , × , ÷ ), exponentials, constants, solutions of algebraic equa-tions, and antiderivatives. It follows directly from the definition that the setof Liouvillian functions is closed under arithmetic operations, composition,and integration. It is also closed under differentiation. It is not closed underlimits and infinite sums.All elementary functions are Liouvillian. Examples of well − knownnfunctions which are Liouvillian but not elementary are the nonelementaryintegrals, for example: The error function, the Fresnel integrals. All
Li-ouvillian solutions are solutions of algebraic differential equations, but notconversely.Examples of functions which are solutions of algebraic differential equa-tions but not Liouvillian include: the Bessel functions (except special cases),and the hypergeometric functions (except special cases). Such a special caseare the truncated confluent hypergeometric functions which are going to beof particular importance in this study. More generally, all functions whichare represented as power series ( not truncated ) are not
Liouvillian. A moreprecise definition involves the notion of Liouvillian field [22] and it is notgoing to be given here.Let η be a (non − zero) Liouvillian solution of the differential equation(42). It follows that every solution of this differential equation is Liouvillian.Indeed the method of reduction of order produces a second solution, namely η Z e − R a η . (43)This second solution is evidently Liouvillian and the two solutions are lin-early independent. Thus any solution, being a linear combination of thesetwo, is Liouvillian.A well − known change of dependent variable may be used to eliminate21he term involving y ′ from the differential equation (42). Let z = e R a y. (44)Then Equation (42) yields z ′′ + (cid:18) b − a − a ′ (cid:19) z = 0 . (45)Equation (45) still has rational function coefficients and evidently (seeEquation (44)) y is Liouvillian if and only if z is Liouvillian. Thus nogenerality is lost by assuming that the term involving y ′ is missing from thedifferential Equation (42). Before giving the main result obtained by Kovacic[22] we first introduce some notation and some terminology. C denotes thecomplex numbers and C ( x ) the rational functions over C . A function ω of x is called an algebraic function of degree k, where k is a positive integer,when ω solves an irreducible algebraic equationΠ( ω, x ) = k X i=0 P i ( x )(k − i)! ω i = 0 (46)where P i ( x ) are rational functions of x . Let ν ∈ C ( x ) (to avoid triviality, ν ∈| C ). Then the following holds ([22]) Theorem 1
Equation z ′′ = νz, , ν ∈ C ( x ) (47) has a Liouvillian solution if and only if it has a solution of the form z = e R ω d x (48) where ω is an algebraic function of x of degree 1,2,4,6 or 12. The search of Kovacic’s algorithm for Π( ω, x ) is based on the knowledgeof the poles of ν and consists in constructing and testing a finite number ofpossible candidates for Π( ω, x ). If no Π( ω, x ) is found then the differentialEquation (47) has no Liouvillian solutions. If such a Π( ω, x ) is found and ω is a solution of the Equation (46) then the function η = e R ω d x is a Liouvilliansolution of (47). If ν ( x ) = s ( x ) t ( x ) , (49)22ith s, t ∈ C [ x ], ( C [ x ] denotes the polynomials over C), relatively prime,then the poles of ν are the zeros of t ( x ) and the order of the pole is themultiplicity of the zero of t . The order of ν at ∞ , o( ∞ ), is defined aso( ∞ ) = max(0 , o s − d o t ) (50)where d o s and d o t denote the leading degree of s and t respectively (fora justification of this definition see [60] page 10; Kovacic originally gave adifferent definition, see [22] page 8). In the Appendix we give an outlineof Kovacic’s algorithm. This is the outline of an improved version of thealgorithm given by Duval and Loday − Richaud [60].Now, besides the original formulation of this algorithm [22] we have itsseveral versions and improvements [60], [61], [62] and extensions to higherorder equations [63], [64], [65]. The formulation of the Kovacic algorithmgiven in [62] is alternative to its original form [22] and to that presentedpreviously. It seems that it is much more convenient for computer imple-mentation and it has been implemented in Maple. However, for differentialequations with simple structure of singularities, and depending on parame-ters it seems that the previous form of the algorithm, which is the form givenin [60], is well suited. Moreover, the original formulation of the algorithm[22] consists in fact of three separated algorithms each of them repeatingsimilar steps. In [60] one can find a modification of the original formulationunifying and improving these three algorithms in one. This form is veryconvenient for applications and it is the one employed here.
We apply Kovacic’s algorithm to the CEV ODE (Equation (27))12 α S β +2 d C d S + ( r − q ) S d C d S − ( r − q + λ ) C = 0Dividing both sides of Equation (27) by the coefficient of d C d S we bring it tothe form of Equation (42)d C d S + 2( r − q ) α S β +1 d C d S − r − q + λ ) α S β +2 C = 0 . (51)A well − known change of dependent variable (Equation (44)) may beused to eliminate the term involving d C d S from the differential equation (51).Let C = e R r − q ) α S β +1 d S C. (52)23hen Equation (51) yields d C d S − ν C = 0 , (53)where ν = ( q − r ) + α ((2 β − q − r ) + 2 λ ) S β α S β ≡ s ( S ) t ( S ) . (54)We note that Equation (52) implies that C is Liouvillian if and only if C isLiouvillian. Thus it suffices to apply Kovacic’s algorithm to Equation (53).Two remarks are now in order regarding Equation (53):1. Kovacic’s algorithm can be applied to Equation (53) if ν is a ratioof relatively prime polynomials s ( S ) and t ( S ) of S . This implies inparticular that 2 β has to be an integer.2. Kovacic’s algorithm starts with assigning orders to the poles of ν , i.e.to the zeros of t ( S ), and the order of the pole is the multiplicity of thezero of t ( S ). It also assigns order to ∞ . The order of ν at ∞ , o( ∞ ),is defined as (Equation (50)) o( ∞ ) = max(0 , o s − d o t ) where d o s and d o t denote the leading degree of s and t . This is related to theorder of ∞ as a zero of ν . The cases to be considered in order to decideif Equation (53) admits Liouvillian solutions depend crucially on theorders associated to the poles of ν and to ∞ .The first remark restricts β to take half − integer values only. The secondremark implies that if different half − integer values of β result in differentorders of the poles of ν , and of ∞ , then these half − integer values of β haveto be considered separately.In the case of Equation (53), ν has one single pole, namely the number0. Let o(0) denote the order of 0. An easy calculation gives o(0 ) o( ∞ ) β = 2,3,... | β | β = 1 β = 0 β = − β = − − | β | The structure of Kovacic’s algorithm is such that the cases 2 β = 2 , , ... and2 β = − , − , ... can be considered together. Thus we are left with five cases24e have to consider st case 2 β = 2,3,...2 nd case 2 β = 13 nd case 2 β = 0 (55) th case 2 β = − th case 2 β = − − The implementation of Kovacic’s algorithm in each one of these casesis by no means trivial and it poses its own problems and challenges. Wepresent here in detail the application of Kovacic’s algorithm in the firstcase and report the results in all the other cases. To give in detail theapplication of Kovacic’s algorithm in all cases would take us too far afield.Detailed application of Kovacic’s algorithm in the other four cases will begiven elsewhere [66]. β = 2,3,... We apply Kovacic’s algorithm to Equation (53)d C d S = ν C . We note that in the Appendix we give only part of Kovacic’s algorithm. Itis the part which is necessary for the application of Kovacic’s algorithm tothe case 2 β = 2 , , ... . For the complete version of Kovacic’s algorithm weuse in this paper see [60]. Input:
From Equation (54) we have ν = ( q − r ) + α ((2 β − q − r ) + 2 λ ) S β α S β ≡ s ( S ) t ( S ) . The partial fraction expansion for ν is the following ν ( S ) = (cid:0) q − rα (cid:1) S β + (2 β − q − r )+2 λα S β . (56) First step: a. From Equation (54) we obtain t ( S ) = S β . Hence , m = 2 + 4 β Γ ′ = { } , Γ = { , ∞} (57)o(0) = 2 + 4 β d o s = 2 β d o t = 2 + 4 β (58)o( ∞ ) = max (0 , β − − β ) = max (0 , − β ) = 0 (59) m + = max ( m, o( ∞ )) = max (2 + 4 β,
0) = 2 + 4 β Γ = {∞} (60)Γ β = { } . (61) Equations (57) give γ = 0 and γ = γ = 0 . (62) Equations (59) and (62) implyL = { } . (63) n = 1 . (64) Second step:2a. ∞ ∈ Γ and consequently E ∞ = h (n) { , , ..., n } = h (1) { , } .Therefore E ∞ = { , } . (65) We have n=1 (Equation (64)), 0 ∈ Γ β ) (Equation (61)),q = 1 + 2 β > , since 2 β = 2 , , ... . Equation (128) gives (cid:2) √ ν (cid:3) = a S β + X i=2 β µ i , S i . (66)From Equations (56) and (88) by using undetermined coefficients we obtaina = (cid:18) q − rα (cid:19) , (67) µ β, = µ β − , = µ β − , = ... = µ , = 0 . (68)Thus there are two possibilities for a , one being the negative of the other,and any of these may be chosen. We choosea = q − r α . (69)26quation (129) yields ν − (cid:2) √ ν (cid:3) = b S β + O (cid:18) S β (cid:19) . (70)Equations (56), (67), (68) and (70) implyb = (2 β − q − r ) + 2 λα . (71)From Equations (69) and (71) we haveE = (cid:26) (cid:18) q + ǫ b a (cid:19) | ǫ = +1 (cid:27) = (cid:26) β + λq − r , − λq − r (cid:27) . (72)From Equation (131) it follows that the function “Sign” S with domain E is defined as followsS (cid:18) β + λq − r (cid:19) = 1 , S (cid:18) − λq − r (cid:19) = − . (73) Third Step:3a.
For each family e = (e c ) c ∈ Γ of elements e c ∈ E c we calculate thedegree d (e) of the corresponding, prospective polynomial P. The sets E c aregiven by Equations (65) and (72). e e ∞ d = 1 − X c ∈ Γ e c Families − λq − r λq − r F11 − λq − r − λq − r F22 β + λq − r − β − λq − r F32 β + λq − r − β − λq − r F4In the last column we enumerated the different families.
If d is a non − negative integer n, the family should be retained,otherwise the family is discarded. This makes λ , in each of the families F1,F2, F3, and F4, to become a function of q, r, β, and n . amilies λ F1 λ = n( q − r ) (74)F2 λ = (n + 1)( q − r ) (75)F3 λ = − (2 β + n − q − r ) (76)F3 λ = − (2 β + n)( q − r ) (77) . For each family retained from step , we form the rationalfunction θ given by Equation (133). Since from Equations (60), (61) and(64) we have respectively Γ = {∞} , Γ β = { } , and n=1, Equation (133)implies θ = e S + S(e ) (cid:2) √ ν (cid:3) , (78)where e c denotes any element of E c , S(e ) are given by Equations (73) and[ √ ν ] is given by Equations (88), (67), and (68). By making use of Equation(78) for each of the retained families we obtain θ Families − λq − r S − q − rα S β F1 (79)1 − λq − r S − q − rα S β F2 (80)2 β + λq − r S + q − rα S β F3 (81)2 β + λq − r S + q − rα S β F4 (82)The functional form of θ is the same both in Families F1 and F2 (Equations(79) and (80)) and in Families F3 and F4 (Equations (81) and (82)). Howeverthe functions θ are different both in Families F1 and F2 and in Families F3and F4 since λ is different in all Families F1, F2, F3, and F4 (Equations(74), (75), (76), and (77)). Fourth step - Output: = − P , (83)P = P ′ + θ P , (84)and P − = 0 = P ′′ + 2 θ P ′ + ( θ + θ ′ − ν )P . (85)Combining Equations (46), (83) and (84) yields ω = P ′ P + θ. (86)For each of the retained families, i.e. for each of the families F1, F2, F3,and F4, we search for a polynomial P of degree d (as defined in step )such that Equation (85) is satisfied. If such a polynomial P is found then ω is given by Equation (86) and the function η = C = e R ω d r = e R ( P ′ P + θ )d r = P e R θ d r (87)is a Liouvillian solution of Equation (53). Then the change of the dependentvariable (52) C = e − R r − q ) α S β +1 d S C gives a Liouvillian solution to the CEV ODE (51).The application of the Fourth step of Kovacic’s algorithm, when 2 β =2 , , ... , is given in detail in Section 9.Now the algorithm can be considered complete when n=1. − integervalues of β For each of the retained families, i.e. for each of the families F1, F2, F3,and F4, we search for a polynomial solution P of degree d, as defined in step , such that Equation (85)P ′′ + 2 θ P ′ + ( θ + θ ′ − ν )P = 0is satisfied. 29he function ν is given by Equation (54) ν = ( q − r ) + α ((2 β − q − r ) + 2 λ ) S β α S β . Given that λ in the four families F1, F2, F3, and F4, is given respectivelyby Equations (74), (75), (76), and (77), the function ν in the four familiesreads ν Families (cid:0) q − rα (cid:1) S β + (2( β +n) − q − r ) α S β F1 (88) (cid:0) q − rα (cid:1) S β + (2( β +n)+1)( q − r ) α S β F2 (89) (cid:0) q − rα (cid:1) S β − (2( β +n) − q − r ) α S β F3 (90) (cid:0) q − rα (cid:1) S β − (2( β +n)+1)( q − r ) α S β F4 (91)and the function θ in the four families reads θ Families − n S − q − rα S β F1 (92) − n S − q − rα S β F2 (93)1 − n S + q − rα S β F3 (94) − n S + q − rα S β F4 (95)where n is a non − negative integer.With the function ν given in the four families F1, F2, F3, and F4, byEquations (88), (89), (90), and (91) respectively, and the function θ given inthe four families F1, F2, F3, and F4, by Equations (92), (93), (94), and (95)respectively, Equation (85) in the four families F1, F2, F3, and F4 reads30 ′′ + 2 − n S − q − rα S β ! P ′ + n(n − S P = 0 F1 (96)P ′′ + 2 − n S − q − rα S β ! P ′ + n(n + 1) S P = 0 F2 (97)P ′′ + 2 − n S + q − rα S β ! P ′ + n(n − S P = 0 F3 (98)P ′′ + 2 − n S + q − rα S β ! P ′ + n(n + 1) S P = 0 F4 (99)We search for polynomial solutions P of degree n, where n is a non − negativeinteger, to the Equations (96), (97), (98), and (99). It is appropriate at thispoint to recall the following definition: The function F ( e, q ; u ) = ∞ X k =0 ( e ) k ( q ) k u k k ! , (100)where the symbol ( w ) k , is the Pochammer’s symbol, and is defined by( w ) k = w ( w + 1) ... ( w + k − , (101)is called confluent hypergeometric function of the first kind or Kummer’sfunction of the first kind. When e = − m, m being a non − negative integer, F ( e, q ; u ) is truncated and it reduces to a polynomial P m ( u ) = F ( − m , q ; u )of degree m F ( − m , q ; u ) = 1 − m q u + m(m − q ( q + 1) u
2! + ... + ( − m m! q ( q + 1) ... ( q + m − u m m! . (102)We easily find that the two − dimensional solution spaces of Equations(96), (97), (98), and (99) are spanned respectively by the following four pairsof functions 1 st pair1 st pair1 st pair f ( S ) = S n F (cid:18) − n2 β , − β ; r − qa βS β (cid:19) (103) f ( S ) = S n − F (cid:18) − n − β , β ; r − qa βS β (cid:19) (104)31e note that f ( S ) truncates when n is a multiple of 2 β and becomes apolynomial of degree n, and f ( S ) truncates when n − β , β = 1 , , , , ..., and it becomes a polynomial of degree n − f ( S ) and f ( S ) span the solution space of Equation (96).2 nd pair2 nd pair2 nd pair f ( S ) = S n F (cid:18) − n2 β , β ; r − qa βS β (cid:19) (105) f ( S ) = S n+1 F (cid:18) − n + 12 β , − β ; r − qa βS β (cid:19) (106)We note that f ( S ) truncates when n is a multiple of 2 β and becomes apolynomial of degree n, and f ( S ) truncates when n + 1 is a multiple of 2 β , β = 1 , , , , ..., and it becomes a polynomial of degree n + 1. f ( S ) and f ( S ) span the solution space of Equation (97).3 nd pair3 nd pair3 nd pair f ( S ) = S n F (cid:18) − n2 β , − β ; q − ra βS β (cid:19) (107) f ( S ) = S n − F (cid:18) − n − β , β ; q − ra βS β (cid:19) (108)We note that f ( S ) truncates when n is a multiple of 2 β and becomes apolynomial of degree n, and f ( S ) truncates when n − β , β = 1 , , , , ..., and it becomes a polynomial of degree n − f ( S ) and f ( S ) span the solution space of Equation (98).4 th pair4 th pair4 th pair f ( S ) = S n F (cid:18) − n2 β , β ; q − ra βS β (cid:19) (109) f ( S ) = S n+1 F (cid:18) − n + 12 β , − β ; q − ra βS β (cid:19) (110)We note that f ( S ) truncates when n is a multiple of 2 β and becomes apolynomial of degree n, and f ( S ) truncates when n + 1 is a multiple of 2 β , β = 1 , , , , ..., and it becomes a polynomial of degree n + 1. f ( S ) and f ( S ) span the solution space of Equation (99).Therefore we do get polynomial solutions to Equations (96), (97), (98),and (99), two for each equation, provided that n, n −
1, or n+1 is a multiple32f 2 β , depending on the family under consideration. In fact in this caseeach f i ( S ), i = 1 , , ..., , represents a family of polynomials. For example f ( S ) truncates to a polynomial of degree n, for each value of n in theset { β, β, β, ... } , β = 1 , , , , ... . Hereafter when we refer to f i ( S ), i = 1 , , ..., , we will mean the associated denumerably infinite family ofpolynomials. Since we look for polynomial solutions of degree n to Equations(96), (97), (98), and (99), we conclude that there remain four families wehave to consider, namely the families f ( S ) , f ( S ) , f ( S ) , and f ( S ) . By combining Equations (26), (29), (52),(74), (75), (76), (77), and (87),we obtain that the four families of polynomial solutions f ( S ) , f ( S ) , f ( S ) , and f ( S ) , to Equations (96), (97), (98), and (99), give rise respectively tothe following four classes of elementary function solutions to the CEV PDE(25) 1 st class1 st class1 st class S , n ( S ) = SF (cid:18) − n2 β , − β ; r − qa βS β (cid:19) e − n( q − r ) t , (111)where n is a multiple of 2 β , β = 1 , , , , ... . The first class arises from thefamily of polynomials f ( S ) given by Equation (103). We easily check thatfor n=0 we also obtain an elementary function solution to the CEV PDE(25). 2 nd class2 nd class2 nd class S , n ( S ) = F (cid:18) − n2 β , β ; r − qa βS β (cid:19) e − (n+1)( q − r ) t , (112)where n is a multiple of 2 β , β = 1 , , , , ... . The second class arises fromthe family of polynomials f ( S ) given by Equation (105). We easily checkthat for n=0 we also obtain an elementary function solution to the CEVPDE (25). 3 nd class3 nd class3 nd class S , n ( S ) = e r − qa βS β SF (cid:18) − n2 β , − β ; q − ra βS β (cid:19) e (2 β +n − q − r ) t , (113)where n is a multiple of 2 β , β = 1 , , , , ... . The third class arises from thefamily of polynomials f ( S ) given by Equation (107). We easily check thatfor n=0 we also obtain an elementary function solution to the CEV PDE(25). 33 th class4 th class4 th class S , n ( S ) = e r − qa βS β F (cid:18) − n2 β , β ; q − ra βS β (cid:19) e (2 β +n)( q − r ) t , (114)where n is a multiple of 2 β , β = 1 , , , , ... . The fourth class arises from thefamily of polynomials f ( S ) given by Equation (109). We easily check thatfor n=0 we also obtain an elementary function solution to the CEV PDE(25). This completes the consideration of the problem when 2 β = 1 , , , ... . There remain four cases to consider 2 β = 1 , β = 0 , β = − , and2 β = − , − , ... (Equation(55)). The implementation of Kovacic’s algorithmin each one of these cases is by no means trivial and it poses its own problemsand challenges [66]. The results in all cases are summarized in the followingTheorem. Theorem 2
The CEV PDE (Equation (25)) ∂V∂t + 12 α S β +2 ∂ V∂S + ( r − q ) S ∂V∂S − ( r − q ) V = 0 ,
1. When β = ..., − , − , − , − , , , , , ... admits the following fourclasses of elementary function solutions S , n ( S ) = S F (cid:18) − n2 β , − β ; r − qa βS β (cid:19) e − n( q − r ) t , S , n ( S ) = F (cid:18) − n2 β , β ; r − qa βS β (cid:19) e − (n+1)( q − r ) t , S , n ( S ) = e r − qa βS β S F (cid:18) − n2 β , − β ; q − ra βS β (cid:19) e (2 β +n − q − r ) t , S , n ( S ) = e r − qa βS β F (cid:18) − n2 β , β ; q − ra βS β (cid:19) e (2 β +n)( q − r ) t , where n= 0, or n is any multiple of β ,2. When β = admits the classes S , n ( S ) and S , n ( S ) , where n= 0, or n is any positive integer, When β = − admits the classes S , n ( S ) and S , n ( S ) , where n= 0, or n is any negative integer, When β = 0 admits the elementary function solution e − λt (cid:18) gS ξ ζ σ + hS ξ − ζ σ (cid:19) , where, ξ = σ − r − q ) , ζ = p ( σ − r − q )) + 8 σ ( r − q + λ ) , and g, h, and λ are arbitrary real numbers (Equations (39), (41)).Since the CEV PDE is linear an immediate Corollary of Theorem 2 isthe following Corollary 1
The CEV PDE (Equation (25)) ∂V∂t + 12 α S β +2 ∂ V∂S + ( r − q ) S ∂V∂S − ( r − q ) V = 0 ,
1. When β = ..., − , − , − , − , , , , , ... admits the following ele-mentary function solution F ( S ) = q X k =1 m k S ,µ k ( S )+ r X k =1 z k S ,ζ k ( S )+ t X k =1 w k S ,ω k ( S )+ j X k =1 d k S ,δ k ( S ) , (115) where q, r, t, j are any positive integers, equal to one or greater thanone, m k , z k , w k and d k are arbitrary real numbers, and µ k , ζ k , ω k and δ k are any multiples of β or 0,2. When β = admits the following elementary function solution F ( S ) = g X k =1 e k S ,ε k ( S ) + h X k =1 y k S ,υ k ( S ) , (116) where g and h are any positive integers, equal to one or greater thanone, e k and y k are arbitrary real numbers, and, ε k and υ k are anypositive integers or 0,3. When β = − admits the following elementary function solution F ( S ) = c X k =1 p k S ,π k ( S ) + l X k =1 r k S ,̺ k ( S ) , (117) where c and l are any positive integers, equal to one or greater thanone, p k and r k are arbitrary real numbers, and, π k and ̺ k are anynegative integers or 0. S , n ( S ), given by Equation (111), is contained implicitly inthe solution C ( S ) given by Equation (31), the class S , n ( S ), given byEquation (112), is contained implicitly in the solution C ( S ) given byEquation (32), whereas the classes S , n ( S ) and S , n ( S ), given respec-tively by Equations (113) and (114), are not contained in the solutions C ( S ) and C ( S ) .
2. When β = 0 the general solution of the CEV ODE is Liouvillian(Equations (39), (41)). Thus it is no surprise that in this case Kovacic’salgorithm gives the two Liouvillian solutions S ξ ζ σ and S ξ − ζ σ whichspan the solution space of the CEV ODE.3. It is well known (Equation (43)) that for every Liouvillian solution f ( S ) to the CEV ODE (Equation (27))12 α S β +2 d C d S + ( r − q ) S d C d S − ( r − q + λ ) C = 0there is a second Liouvillian solution L ( S ) = f ( S ) R e − R r − q ) α S β +1 d S f ( S ) d S to the CEV ODE which also gives rise (Equation (26)) to a solution L ( S ) = f ( S ) R e − R r − q ) α S β +1 d S f ( S ) d Se − λ t to the CEV PDE (Equation(25)) ∂V∂t + 12 α S β +2 ∂ V∂S + ( r − q ) S ∂V∂S − ( r − q ) V = 0 . In all cases the integral R e − R r − q ) α S β +1 d S f ( S ) d S is nonelementary and it canonly be evaluated by using Taylor’s series, or numerically, by usingquadrature rules such as the Simson’s rule or Gaussian rule. Con-sequently some kind of approximation scheme is needed in order toextract information from the solutions L ( S ) for the pricing of the fi-nancial instruments they describe. For this reason we do not includethe solutions L ( S ) to the classes of elementary function solutions tothe CEV model we derive in this paper.
10 Conclusion and Future Development
Analytical tractability of any financial model is an important feature.Existence of a closed − form solution definitely helps in pricing financial in-36truments and calibrating the model to market data. It also helps to verifythe model assumptions, check its asymptotic behavior and explain causality.In fact, in mathematical finance many models were proposed, first based ontheir tractability, and only then by making another argument.If the true asset price process was geometric Brownian motion with con-stant volatility, then the BSM PDE could be used to find out this volatilityby equating the model price of a standard option to its market price (impliedvolatility).Empirically, we find that the implied volatilities computed from marketprices of options with different strike prices are not constant but vary withstrike price. This variation is observed across a wide range of markets andunderlying assets and is known as the implied volatility smile or frown de-pending on its shape. The lognormal assumption with constant volatility ofthe BSM model does not capture this effect. The CEV model is capable ofreproducing the volatility smile observed in the empirical data.Our classes of elementary function solutions to the CEV model allow fastand accurate calculation of prices of various financial instruments under theCEV process. Moreover they will facilitate further the use of the CEV modelas a benchmark for the pricing of various types of financial instruments.In future research we will study the Lie point symmetries of the CEVPDE and we will use the classes of elementary function solutions to theCEV PDE we derived in this paper in order to obtain more solutions to theCEV model in terms of elementary functions. The study of the financialinstruments which they describe will be useful for trading by using CEVdiffusion in all cases. APPENDIX
Notations.
Let L max = { , , , , } and let h be the function defined onL max by h (1) = 1 , h (2) = 4 , h (4) = h (6) = h (12) = 12 . (118) Input:
A rational function ν ( x ) = s ( x ) t ( x ) (Equation(49)) . The polynomials s, t ∈ C [ x ] are supposed to be relatively prime.The differential equation under consideration is z ′′ − νz = 0 (Equation(47)) . irst step: The set L of possible degrees of ω. We are interested in Equation (47) where ν ( x ) is given by (49). . If t ( x ) = 1, set m = 0. Otherwise, factorize t ( x ) t = t t t ...t mm (119)where the t i , i = 1 , , ..., m, are relatively prime two by two and each t i either is equal to one or has simple zeros. LetΓ ′ = { c ∈ C , t ( c ) = 0 } and Γ = Γ ′ ∪ {∞} , (120)where ∪ denotes set − theoretic union. Associate orders with the elements ofΓ: o( c ) = i (121)for all c ∈ Γ ′ ,where i is such that t i ( c ) = 0, ando( ∞ ) = max(0 , o s − d o t ) . Let m + = max( m, o( ∞ )) . (122)For 0 ≤ i ≤ m + let Γ i = { c ∈ Γ | o( c ) = i } . (123) . If m + ≥ γ and γ by γ = | Γ | and γ = γ + |∪ Γ k | , k is odd , and 3 ≤ k ≤ m + (124)respectively, where, if S is a set then | S | denotes the number of elements ofS. . Construct L, a subset of L max , as follows1 ∈ L ⇐⇒ γ = γ , (125)2 ∈ L ⇐⇒ γ ≥ , and , (126)4 , ∈ L ⇐⇒ m + ≤ . (127) . If L = ∅ go to the stage of the algorithm END. Otherwise, let n beequal to the smallest element of L. Second step:
The sets E c associated to the singular points. Construction of the sets E c , c ∈ Γ. . If ∞ ∈ Γ then E ∞ = h (n) { , , ..., n } . b . When n=1, for each c ∈ Γ with q ≥
2, calculate one of the two“square roots” [ √ ν ] c of ν defined as follows:If c ∈ C , (cid:2) √ ν (cid:3) c = a c ( x − c ) q + X i=q − µ i ,c ( x − c ) i , (128)and ν − (cid:2) √ ν (cid:3) c = b c ( x − c ) q+1 + O (cid:18) x − c ) q (cid:19) , (129)where the O symbol has its usual meaning, f ( x ) = O( g ( x )) if | f ( x ) | ≤ M g ( x )for some positive constant M; g ( x ) is assumed to be positive.Let E c = (cid:26) (cid:18) q + ǫ b c a c (cid:19) | ǫ = +1 (cid:27) . (130)Define a function “Sign” S with domain E c as followsS (cid:18) (cid:18) q + ǫ b c a c (cid:19)(cid:19) = (cid:26) ǫ if b c = 01 otherwise . (131)In Equation (128) [ √ ν ] c is the sum of terms involving ( x − c ) − i for 2 ≤ i ≤ qin the Laurent series for √ ν at c . In practice, one would not form the Laurentseries for √ ν , but rather would determine [ √ ν ] c by using undetermined coef-ficients, i.e. by equating ([ √ ν ] c ) = (cid:0) a c ( x − c ) − q + µ q − ,c ( x − c ) − (q − + ... + µ ,c ( x − c ) − (cid:1) with the corresponding part of the Laurent series expan-sion of ν at c . There are two possibilities for [ √ ν ] c , one being the negativeof the other, and any of these may be chosen. In Equation (129), ν denotesits Laurent series expansion at c . This equation defines b c . Third step:
Possible degrees for P and possible values for θ . . For each family e = (e c ) c ∈ Γ of elements e c ∈ E c calculated (e) = n − n h (n) X c ∈ Γ e c (132) . Retain the families e which are such thatd(e) ∈ N where N is the set of non − negative integers andIf none of the families e is retained, go to the stage of the algorithm CONTINUATION. . For each family e retained from step , form the rational function θ = n h (n) X c ∈ Γ ′ e c x − c + δ X c ∈ ∪ Γ q ≥
2S (e c ) (cid:2) √ r (cid:3) c , (133)39here δ is the Kronecker symbol. Fourth Step:
Tentative computation of P.
Search for a polynomial P of degree d (as defined in step ) such thatP n = − P · · · P i − = − P ′ i − θ P i − (n − i)(i + 1) ν P i+1 · · · P − = 0 , (134)where the prime in P ′ ( x ) denotes differentiation with respect to the inde-pendent variable x . Output:
OUTPUT1: If such a polynomial is found and ω is a solution of the irre-ducible algebraic equation n X i=0 P i ( x )(k − i)! ω i = 0 (Equation (46)) , where the rational functions P i ( x ) are defined in (134), then the function η = e R ω is a Liouvillian solution of the equation under consideration z ′′ = νz (Equation(47)) . If no such polynomial is found for any family retained from step , go tothe stage of algorithm CONTINUATION . CONTINUATION : If n is different from the largest element of L then setn equal to the next (in increasing order) element of L and go to
Step 2 .OUTPUT2: Equation z ′′ = νz has no Liouvillian solutions. References [1] F. Black and M. Scholes, The valuation of option contracts and a testof market efficiency Journal of Finance Journal of Political Economy −
654 (1973)[3] R. C. Merton, Theory of Rational Option Pricing,
The Bell Journal ofEconomics and Management Science −
183 (1973)404] R. C. Merton, Continuous Time Finance Blackwell (1990)[5] G. Marsaglia, Evaluating the Normal Distribution
Journal of StatisticalSoftware (2004)[6] J. Cox, Notes on option pricing I: Constant elasticity of variance diffu-sions Working paper, Stanford University (1975) (Reprinted in J. Portf.manage. −
17 (1996))[7] J. Cox and S. Ross, The valuation of options for alternative stochasticprocesses
Journal of Financial Economics −
166 (1976)[8] D. C. Emanuel and J. D. MacBeth, Further results on the constantelasticity of variance call option pricing model
Journal of Financialand Quantitative Analysis −
554 (1982)[9] M. Schroder, Computing the constant elasticity of variance option pric-ing formula
Journal of Finance −
219 (1989)[10] S. G. Kou, A Jump − Diffusion Model for Option Pricing
ManagementScience − − Ho Wang,Asymptotics of Implied Volatility in Local Volatility Models
Mathemat-ical Finance published online (2010)[12] K. C. Yuen, H. Yang and K. L. Chu, Estimation in the Constant Elas-ticity of Variance Model British Actuarial Journal −
292 (2001)[13] R. W. Farebrother, Algorithm AS 231: The Distribution of a Noncen-tral χ Variable with Nonnegative Degrees of Freedom
Applied Statistics −
405 (1987)[14] Harry O. Posten, An Effective Algorithm for the NoncentralChi − Squared Distribution Function
American Statistician − − Central χ Distribution Function
Applied Statistics −
482 (1992)[16] L. Kn¨usel and B. Bablok, Computation of the Noncentral Gamma Dis-tribution,
SIAM Journal of Scientiflc Computing − Computational Statistics and Data Analysis − − Square Distribution forthe Cox − Ingersoll − Ross Process
Computational Economics − Journal of Portfolio Management (1994)[20] A. Neuberger, The log contract and other power contracts, chapter 7in The Handbook of exotic options: Instruments, Analysis and Appil-cations
McGraw − Hill (1996)[21] E. Boussard, Trading and Hedging Volatility, Presentation at Deriva-tives ’97, Risk Conference, Brussels, February 1997[22] J. Kovacic, An algorithm for solving second order linear homogeneousdifferential equations
J. Symb. Comp. −
43 (1986)[23] L. Bachelier, Th´eorie de la speculation − neutral skewness: Evidence from stockoptions Journal of Financial and Quantitative Analysis − Review of Financial Studies −
42 (2000)[26] J. D. MacBeth and L. J. Merville, Tests of the Black − Scholes and Coxcall option valuation models
Journal of Finance −
301 (1980)[27] J. C. Dias and J. P. V. Nunes, Pricing Real Options under the CEVDiffusion
Journal of Futures Markets −
250 (2011)[28] D. Davydov and V. Linetsky, Pricing and Hedging Path − DependentOptions under the CEV Process
Management Science − Finance and Stochastics Operations Research − De-cisions in Economics and Finance −
153 (2006)[32] P. P. Boyle, Y. Tian, and J. Imai, Lookback options under theCEV process: a correction
Journal of Financial and QuantitativeAnalysis Unpublished Appendixes, Notes, Comments, and Correctionshttp://depts.washington.edu/jfqa/ (1999)[33] P. P. Boyle and Y. Tian, Pricing lookback and barrier options underthe CEV process
Journal of Financial and Quantitative Analysis −
264 (1999)[34] P. Wilmott, J. Dewynne, and S. Howison, Option Pricing: Mathemat-ical models and computation Oxford Financial Press, Oxford (1993)[35] W. T. Shaw, Modelling Financial Derivatives with Mathematica Cam-bridge University Press (1998)[36] P. J. Olver, Applications of Lie Groups to Differential Equations,Springer − Verlag (1993)[37] S. Beckers, The Constant Elasticity of Variance Model and Its Impli-cations for Option Pricing
The Journal of Finance −
673 (1980)[38] HANDBOOK OF MATHEMATICAL FUNCTIONS with Formulas,Graphs, and Mathematical Tables Edited by M. Abramowitz and I.A. Stegun Addison − Wesley, U.S.A., and Scientific Computing ServiceLtd., Great Britain (1962)[39] P. P. Boyle, W. Tian and F. Guan, The Riccati Equation in Mathe-matical Finance
J. Symbol. Comput. −
355 (2002)[40] I. Hern´andez, C. Mateos, J. N. Vald´es and A. F. Tenorio, Lie Theory:Applications to problems in Mathematical Finance and Economics
Ap-plied Mathematics and Computation −
452 (2009)[41] E. Galois, M`emoire sur les conditions de r`esolubilit`e des `equations parradicaux (Dated 6 January 1831)
Journal de Math`ematiques pures etappliqu`ees (Ed. J. Liouville) (1846)4342] E. Galois, Des `equations primitives qui sont soubes par radicaux
Jour-nal de Math`ematiques pures et appliqu`ees (Ed. J. Liouville) (1846)[43] E. Galois, Lettre `a Auguste Chevalier (Dated 29 May 1832) Publishedin
Revue encyclop´edique (Sept. 1832)[44] S. Lie, Die Grundlagen f¨ur die Theorie der unendlichen kontinuier-lichen Transformationsgruppen. 1. Abhandlung
Abh. Ges. Wiss. Leipzig (1891) reprinted in Lie, S.: Gesammelte Abhandlungen Leipzig/Oslo300 −
330 (1927)[45] S. Lie, Untersuchungen ¨uber unendliche kontinuierliche Gruppen
Abh.Ges. Wiss. Leipzig (1895) reprinted in Lie, S.: Gesammelte Abhand-lungen Leipzig/Oslo 396 −
493 (1927)[46] E. Picard, Trait´e d’analyse (Deuxieme Ed.) Gauthier − Villars (1908)[47] E. Vessiot, Sur l’int´egration des ´equations diff´erentielles lin´eaires
An-nales Scientifiques de l’ ´Ecole Normale Sup´erieure Encyclop´edie des sci-ences math´ematiques pures et appliqu´ees Gauthier − Villars & Teub-ner, 58 −
170 (1910)[49] N. H. Abel, M´emoire sur une classe particuli`ere d’ ´equations r´esolublesalg´ebriquement
Journ. Crelle (1829)[50] B. J. Cantwell, Introduction to Symmetry Analysis Cambridge Univer-sity Press (2002)[51] W. R. Oudshoorn and van der Put M., Lie symmetries and differentialGalois groups of linear equations Math. Comp. −
361 (2002)[52] N. H. Ibragimov, A bridge between Lie symmetries and Galois groupsPublished in
Differential Equations: Geometry, Symmetries and Inte-grability
Proc. 5th Abel Symposium, Troms¨o, Norway, Hune 17 − −
172 (2009)[53] D. BL ´ A ZQUEZ − SANZ, J. J. MORALES − RUIZ and J. A. WEIL, Dif-ferential Galois Theory and Lie Symmetries
SIGMA
092 (2015)[54] E. R. Kolchin, Algebraic matric groups and the Picard-Vessiot theoryof homogeneous linear ordinary differential equations
Annals of Math-ematics −
42 (1948) 4455] E. R. Kolchin, Differential Algebra and Algebraic Groups AcademicPress (1976)[56] A. Boulanger, Contribution `a l’´etude des ´equations lin´eaireshomog`eenes int`egrables alg`ebriquement
Journal de l’ ´ E cole Polytech-nique, Paris −
122 (1898)[57] M. F. Singer, Liouvillian solutions of nth order homogeneous lineardifferential equations
American Journal of Mathematics − Grundlehrender mathematischen Wissenschaften
Springer (2003)[59] M. F. Singer, Introduction to the Galois Theory of Linear DifferentialEquationsThis is an expanded version of the 10 lectures given as the 2006 Lon-don Mathematical Society Invited Lecture Series at the Heriot − WattUniversity 31 July − − Richaud, A propos de l’algotitme de Ko-vacic
Tech. Rep.
Universite de Paris − Sud Mathematiques Orsay France(1989)[61] A. Duval and M. Loday − Richaud, Kovacic’s Algorithm and its Ap-plication to Some Families of Special Functions
AAECC − J. Symbol. Com-put. (2) 179 −
200 (1996)[63] M. F. Singer and F. Ulmer, Liouvillian and algebraic solutions of secondand third order linear differential equations
J. Symbol. Comput. (1)37 −
73 (1993)[64] M. F. Singer and F. Ulmer, Necessary conditions for Liouvillian solu-tions of (third order) linear differential equations
AAECC (1) 1 − J. SymbolComput. −
609 (1998) 4566] E. Melas Picard − Vessiot −−