Numerical Calculation
A theme of this book is that Riemann sums are an effective means of analyzing random variability phenomena, enabling a comprehensive theory to be constructed. The traditional theory of probability is based on measurable sets formed by operations involving infinitely many steps. In contrast, Riemann sums have only a finite number of terms. On the face of it, such sums should be relatively easy to calculate. This chapter contains a number of such calculations, using Maple 15 computer software.
Numerical and empirical investigations warrant detailed study in their own right. But such an project is beyond the scope of this book. All that can be given are a few indications and pointers.
In preceding chapters Riemann sums have been used to calculate and analyze measurability of sets and functions, expectation values of random variables, state functions of diffusion systems and quantum mechanics, Feynman diagrams, valuation of share options, and strong and weak stochastic integrals.
In order to produce a robust theory, considerable subtlety has been built into the construction of Riemann sums. By making the Riemann sums conform to construction rules called gauges it has been possible to establish rules and criteria for sophisticated mathematical operations by which many classical results could be deduced, and also many new results which are beyond the reach of traditional theory.
But the essential simplicity of Riemann sums is a constant feature. Furthermore, in the various themes discussed in this book, the expressions encountered in the construction of Riemann sums were the familiar polynomial, exponential, trigonometric, and logarithmic functions. These are not expressions of an exotic or pathological kind which require very delicately constructed partitions.
In one dimension, integrals of polynomial functions can be estimated with Riemann sums of a straightforward kind, without having to resort either to the simple functions of Lebesgue theory, nor even to the δ-gauges of Riemann-complete theory. Similarly, good estimates of the Henstock integrals of the functions investigated in this book can often be obtained with regular partitions—or even binary partitions—of the infinite-dimensional domain of integration.
To sum up, the calculations involve finite Riemann sums. And in many cases the partitions involved are particularly amenable to numerical calculation. So it should not be a surprise that many of the steps involved in numerical calculation of the themes of this book can be illustrated with Maple software.
9.1 Introduction
This chapter presents Maple estimates of some of the calculations encountered in Chapters 7 and 8. With 
= ]0, t] and
= ]0, t] and
consider deterministic calculations f performed on the unpredictable or random elements xs which, for 
, are the potential joint outcomes of an experiment or joint observation denoted by X—that is,
, are the potential joint outcomes of an experiment or joint observation denoted by X—that is,
where, for each s, Xs denotes the measurement or process of determination of the individual datum xs.For any potential joint outcome x, the deterministic calculation f(
) is subject to the random variability of each xs for . The random variability in the final outcome f(), induced by the joint random variability of xs (), is measured by a joint potentiality distribution function FX = 
, defined on the events I = I[N] for 
.
) is subject to the random variability of each xs for . The random variability in the final outcome f(), induced by the joint random variability of xs (), is measured by a joint potentiality distribution function FX = , defined on the events I = I[N] for .
For real-valued f, a potential final outcome f() can be regarded as a potential value y 
R. Then y is the unpredictable outcome of an observation Y. And, provided f is measurable and FX is continuous, the likelihoods of values y are measured by a potentiality distribution function FY defined for cells J I(R), with
) can be regarded as a potential value y R. Then y is the unpredictable outcome of an observation Y. And, provided f is measurable and FX is continuous, the likelihoods of values y are measured by a potentiality distribution function FY defined for cells J I(R), with
This is the elementary form Y ≃ y[R, FY] of the contingent joint-basic observable f(
), as described in Sections 5.2 and 5.3.
), as described in Sections 5.2 and 5.3.
We are primarily concerned with Brownian motion X, and, except where otherwise stated, will take FXto be the Brownian distribution function G.
A key point in these illustrations is the representation of the elementary-form observables Y by means of histograms. Using Maple, we will, for particular calculations f on values in 
, demonstrate with histograms both the possible values y of Y and also the distribution function FY (depending on FX and f) which measure the likelihoods in R.
, demonstrate with histograms both the possible values y of Y and also the distribution function FY (depending on FX and f) which measure the likelihoods in R.
Histograms are a useful way of visualizing elementary observables Y ≃ y[R, FY], because they exhibit the range of potential values y R of the observable, and they indicate the values FY(J) of cells J ⊂ R. Constructed, as they are, on a finite number of cells J partitioning the domain R, histograms are particularly apt for representing observables in the Riemann sum frame work of this book. Histograms capture and display the distinguishing features of elementary observables.1
R of the observable, and they indicate the values FY(J) of cells J ⊂ R. Constructed, as they are, on a finite number of cells J partitioning the domain R, histograms are particularly apt for representing observables in the Riemann sum frame work of this book. Histograms capture and display the distinguishing features of elementary observables.1
Of particular interest are the various kinds of calculation f discussed in Chapter 8. These calculations involve
- selection of partition points N = {..., s, s′,...} in ;
- for any potential outcome , forming terms composed of differences or increments such as
- calculating Riemann sums on domain using these terms.
Such calculation depends not just on the potential joint occurrence 
, but also on elements N ⊂ , often in the form of differences sj − sj−1. This is the joint-contingent view, in contrast to the elementary view, and developing the theory more fully in this direction requires ideas from Section A.l of the Epilogue.
, but also on elements N ⊂ , often in the form of differences sj − sj−1. This is the joint-contingent view, in contrast to the elementary view, and developing the theory more fully in this direction requires ideas from Section A.l of the Epilogue.
But if the Riemann sums converge for each , then the limit does not depend on any particular N ⊂ , and when the limit of the Riemann sums is taken, the outcome can be treated as a potential datum of a strong stochastic integral—a well defined elementary observable. The features of such a calculation are displayed in the histogram of Figure 9.4.
, then the limit does not depend on any particular N ⊂ , and when the limit of the Riemann sums is taken, the outcome can be treated as a potential datum of a strong stochastic integral—a well defined elementary observable. The features of such a calculation are displayed in the histogram of Figure 9.4.
When the Riemann sums do not converge, and in default of a suitable theory involving joint-contingent observables with variable N, the “weak stochastic integral” device is available. To visualize an observable defined by this method, take several successive partitions, each one refining the previous one. Then, for each such partition, perform the Riemann sum calculation using a sample of potential joint outcomes , and draw a histogram of the resulting real numbers. By examining the histograms for successive partitions it is possible get a visual sense of the range of possible values y of the weak stochastic integral in each case, along with the “shape” of the corresponding likelihood values FY(J); and also a sense of how these elements change when the partitions are refined successively.
, and draw a histogram of the resulting real numbers. By examining the histograms for successive partitions it is possible get a visual sense of the range of possible values y of the weak stochastic integral in each case, along with the “shape” of the corresponding likelihood values FY(J); and also a sense of how these elements change when the partitions are refined successively.
Figures 9.6, 9.7, and 9.8 illustrate this for 
, where X is standard Brownian motion. Maple allows us to create samples of thousands of Riemann sums, and the histogram of any such sample constitutes a visualization of the particular elementary observable Y ≃ y[R, FY] whose potential datum is
, where X is standard Brownian motion. Maple allows us to create samples of thousands of Riemann sums, and the histogram of any such sample constitutes a visualization of the particular elementary observable Y ≃ y[R, FY] whose potential datum is
These diagrams demonstrate what is meant by “weak stochastic integral”, in which the Riemann sums do not generally converge for particular outcomes y = f(). The same device is applied to Itô’s formula in Section 9.5.
). The same device is applied to Itô’s formula in Section 9.5.9.2 Random Walk
As an introduction, random walk and Brownian path diagrams can easily be produced with Maple.
The following Maple code simulates a random walk whose steps or increments are normally distributed with mean zero, standard deviation 0.1, and variance 0.12 = 0.01. Provided the increments are independent, Theorem 100 says that the variance of the sum of the increments equals the sum of the variances of the increments. Therefore the unit interval = ]0,1] is partitioned by 100 equal intervals of length 0.01.
= ]0,1] is partitioned by 100 equal intervals of length 0.01.restart:
with(Statistics):
randomize():
BrownianIncrements := Sample(Normal(0, 0.1), 100):
ListBrownianIncrements := [seq(BrownianIncrements[j], j=1..100)]:
ordinate[0] := 0:
for j from 1 to 100 do
ordinate[j] := ordinate[j − 1] + ListBrownianIncrements[j]:
end do:
for j to 100 do
abscissa := [seq(0.01*j, j=1..100)]:
end do:
RandomWalk := [[0, 0], seq([abscissa[j], ordinate[j]], j=1..100)]:
plot(RandomWalk);
Line 1 of Calculation 1 cancels any preceding Maple code. Line 2 activates Maple’s statistics functions, such as
Sample. Line 3 resets the random number generation algorithm of the program. Line 4 produces a random sample of 100 numbers, each of them sampled from a normal distribution with mean zero, standard deviation 0.1, and variance 0.12 = 0.01. Line 5 converts these 100 items into a Maple list structure, indexed from 1 to 100. Line 6 assigns 0 as the ordinate (vertical) of the zeroth item. Lines 7 to 9 assign the sum of the first j Brownian increments as the jth ordinate, for j from 1 to 100. Lines 10 to 12 assign the value 0.01 j as the jth abscissa (horizontal coordinate), for j from 1 to 100. Line 13 assembles 101 abscissa–ordinate pairs into a Maple indexed-list structure. Line 14 plots the 101 points of this list, joining each consecutive pair of points with a straight line segment.
Figure 9.1: Random walk with normal increments.
In Figure 9.1, the points are the important thing. The straight line segments joining up the points have no particular significance, but they are helpful in drawing the eye from one point to the next.
For j = 0 to 100, the jth point of the graph is (Sj, x(Sj)), where s = 0.01j and xs is the sum of the first jrandom normal increments; that is,
x(Sj) = xSj =
ordinate [j].
x1 = 0.49372516204936334.
If the Maple program is re-run, 100 random increments produce a different random walk, and a different datum x1. Here is an example:
x1 = −2.268719956985433.
Calculation 4 calculates a sample of 1000 values of x(1). The histogram confirms that, subject to the limitations of the sampling process (which is not truly random), the data values produced by this program are from a normal distribution with mean zero and standard deviation 1.
The following program has 10,000 increments, each with standard deviation 0.01, so their variance is 0.0001; giving us a partition of domain = ]0,1] consisting of 10,000 intervals or steps. The graph of the sample path conforms to the customary representation of Brownian motion.
= ]0,1] consisting of 10,000 intervals or steps. The graph of the sample path conforms to the customary representation of Brownian motion.restart:
with(Statistics):
randomize():
BrownianIncrements := Sample(Normal(0, 0.01), 10000):
ListBrownianIncrements := [seq(BrownianIncrements[j], j=1..10000)]:
ordinate[0] := 0 :
for j from 1 to 10000 do
ordinate[j] := ordinate[j − 1] + ListBrownianIncrements[j]:
end do:
for j from 1 to 10000 do:
abscissa := [seq(0.0001*j, j=1..10000)]:
end do:
BrownianPath := [[0, 0], seq([abscissa[j], ordinate[j]], j=1..10000)]:
plot(BrownianPath);
Figure 9.2: Sample path of Brownian motion.
Figure 9.2 has the familiar2 jagged path shape. By using a bigger sample it is possible to produce a graph that looks less like the random walk of Figure 9.1, and more like the traditional picture of standard Brownian motion.
The graph consists of points joined by straight line segments, just like Figure 9.1. The points have both mathematical and physical significance. The straight line segments joining the points have neither. The only substantial difference between Figure 9.1 and Figure 9.2 is that the latter has more sample points.
With = ]0,1], each of Figure 9.1 and Figure 9.2 is a sample datum of a Riemann sum estimate of the strong stochastic integral
= ]0,1], each of Figure 9.1 and Figure 9.2 is a sample datum of a Riemann sum estimate of the strong stochastic integral
In each case, the partition of is regular (since each cell has the same length), though not binary. The Riemann sum observable in Figure 9.1 is
is regular (since each cell has the same length), though not binary. The Riemann sum observable in Figure 9.1 is
where tj = j10−2, and 
= ](j − 1)10−2, j10−2].
= ](j − 1)10−2, j10−2].
The observables 
and 
have joint-contingent representation
and have joint-contingent representation
respectively, where G is the standard Brownian distribution function. The representation of each of them in elementary form is particularly easy to establish, because of the “telescopic” cancellation of terms in the Riemann sum (in every Riemann sum, in fact). Cancellation gives
The latter has elementary form. Then, writing y = x1, we have
X1 = Y, Y ≃ y[R,N];
where N is the standard normal distribution function N0,1 with mean 0 and standard deviation 1.
These calculations are confirmed by applying Theorems 159 and 160 to the observable 
. They are also confirmed in Maple Calculation 4, and in the larger sample of 10,000 Riemann sums of 100 terms each.
. They are also confirmed in Maple Calculation 4, and in the larger sample of 10,000 Riemann sums of 100 terms each.9.3 Calculation of Strong Stochastic Integrals
In Calculation 3 Maple’s
Sample command is used to select 100,000 random values from a normal distribution with mean zero and standard deviation 0.1. So the variance is 0.01. In Calculation 4 these are taken in groups of 100 values at a time, giving 1000 samples of size 100. If the 100 sampled values were really instances of 100 independent basic normal observables (which, strictly speaking, cannot necessarily be guaranteed by computer software) then the variance of the sum of each group of 100 values would be 100 × 0.01, =1. The following examines the output of 100,000 standard normal increments.restart :
with(Statistics):
BrownianIncrements := Sample(Normal(0, 0.1), 100000):
dx := [seq(BrownianIncrements[i], i=1..100000)]:
Histogram(dx);
Line 1 of Calculation 3 gets rid of the stored results of any previous Maple calculations. Line 2 activates the Maple Statistics functions. Line 3 generates a random sample of 100,000 values of a normal distribution with mean zero and standard deviation 0.1. Line 4 converts this output into a Maple list whose individual elements are indexed from 1 to 100,000. Line 5 produces a histogram of these values, displayed in Figure 9.3.
The theory of random variation can be thought of as a theory of estimation or approximation. Given a single datum (estimate or approximate value), it is not necessarily obvious what it is an estimate of; nor what is the “true value” to which the datum approximates. But if, as in Calculation 3, a great number of estimates can be produced, it may be possible to deduce the “true value”. Calculation 3 provides 100,000 estimates. To see these numbers, just change the colons to semicolons in the Maple code. The program will then display a list of 100,000 numbers—something it may be preferable to avoid.
What is actually needed is not such a list, but a “sense”, or summary, of these numbers as estimates of a “true value”, and that is provided by the histogram, which indicates that:
Figure 9.3: Histogram for 1000 samples, Riemann sum step 0.01.
- the estimates cluster around the value 0; or
- the likelihood that an estimate will be close to zero is much greater than the likelihood it will not be close to zero.
In other words, the histogram is a summary or picture of the observable or random variable for which the numbers produced by Calculation 3 are data.
Inspection of the histogram suggests that the mean value of the sample is zero. Most of the values lie within the range −0.35 to +0.35, suggesting a standard deviation of 0.1. This is confirmed by applying Maple functions to the sample data generated by Calculation 3 :
Mean(dx);
For one particular sample this gave result −0.0003180143098.
StandardDeviation(dx);
This gave output 0.0995857853503455, close to 0.1. Combined with Figure 9.3, this gives some broad confirmation that, for a large sample (100,000), Maple’s random normal increment generator can be relied on, within reason. While statistical independence of the data generated by the
Sample command cannot be guaranteed, this histogram has the appearance of a normally distributed observable with the appropriate mean and standard deviation, allowing us to proceed with due caution.
Calculation 4 below selects 100 consecutive values at a time from the list of 100,000 values produced by Calculation 3, giving 1000 groups of 100 values. The 100 values in each group are taken to be independent Brownian increments. These increments are added in each of the 1000 groups, so the total of the increments in each group can be regarded as:
- the outcome of a random walk of 100 steps, or
- a Riemann sum estimate
of the stochastic integral 
(or 
of the integrand g1 = dxs on the domain = ]0,t], with t = 1; the Riemann sum consisting of 100 terms.
The Maple code produces a sample of 1000 of the potential values of the observable Y ≃ y[R, FY] where each sample value y is obtained by calculating the Riemann sum 
Thus the joint-contingent and elementary forms are, respectively,
Since the strong stochastic integral of integrand g1 exists, with value equal to 
for every partition of , this gives
for every partition of , this gives
When Maple produces the histogram Figure 9.4 of the list of 1000 values of y, it can be expected to confirm the theoretical prediction of Example 59—that occurrences y are sample data of a normally distributed observable with mean 0 and standard deviation 1.
From the latter point of view, the Maple code in Calculation 4 produces Riemann sum estimates of 1000 “independent” instances of the strong stochastic integral
The elementary form of this reduces simply to Y ≃ y[R, N0,1], a normal random variable with mean zero and standard deviation
The histogram of Figure 9.4 provides some empirical confirmation of this, but perhaps not very convincing in visual terms.
for r from 1 to 1000 do
k := 100*r:
j := k − 99:
x1[r] := add(dx[i], i =j..k):
end do:
X1 := [seq(x1[r], r=1..1000)]:
Histogram(X1);
Provided Calculation 4 is run immediately after Calculation 3, the data from the latter feeds into the former. Otherwise Calculation 3 must be re-run. The output of Calculation 4 is the histogram in Figure 9.4.
To check the parameters of this histogram, the command
Mean(X1) gave a result −0.03180143098 from a Maple sample, while StandardDeviation(X1) gave 1.00075019279754. So the mean of a sample of 1000 Riemann sums came to approximately zero, and their standard deviation was approximately 1; while their histogram is at least vaguely reminiscent of the shape of a normal distribution.
Figure 9.4: Histogram for 1000 sums of 100 normal increments.
The cancellation or “telescoping” argument of Section 8.3 establishes that Brownian variability of each joint-data component xs of the contingent observable 
has no effect on the variability of Xt. In other words, simply replicates Xt. This is suggested (even if only faintly) in the histogram of Figure 9.4.
has no effect on the variability of Xt. In other words, simply replicates Xt. This is suggested (even if only faintly) in the histogram of Figure 9.4.
The fact that exists as a strong stochastic integral, equal to its Riemann sum estimates regardless of which partitions of ]0,1] are chosen, makes the preceding remarks unsurprising. It is obvious that the distribution function FX1 of the elementary observable is 
, with t = 1 in this case. So the histogram above provides partial confirmation of the obvious.
exists as a strong stochastic integral, equal to its Riemann sum estimates regardless of which partitions of ]0,1] are chosen, makes the preceding remarks unsurprising. It is obvious that the distribution function FX1 of the elementary observable is , with t = 1 in this case. So the histogram above provides partial confirmation of the obvious.
The differences in shape between the histograms of Figures 9.4 and 9.3 suggest that data generated by Maple should be treated with caution. Randomness and independence are ideal mathematical concepts, and perfect manifestations of them cannot be expected to appear in practice.
In order to get a better simulation, amend the code in Calculation 3 to produce 1,000,000 Brownian increments instead of 100,000, and then amend Calculation 4 to deliver 10,000 (instead of 1000) Riemann sums of 100 terms each. The latter produces the histogram of the Riemann sum values in Figure 9.5. This histogram more plausibly conforms to the theoretical result; that if 
is standard Brownian motion, then, for t = 1, the elementary form of the strong stochastic integral 
has standard normal distribution function N0,1, with mean 0 and standard deviation 1.
is standard Brownian motion, then, for t = 1, the elementary form of the strong stochastic integral has standard normal distribution function N0,1, with mean 0 and standard deviation 1.
Quite large samples are used here as a practical means of compensating for any lack of randomness and independence in Maple-generated sample data. There is a penalty for this, in that the running time of the programs can be relatively long.
9.4 Calculation of Weak Stochastic Integrals
With = ]0,1] this section investigates numerically the weak stochastic integral 
, for which the Riemann sum estimates do not generally converge whenever the partitions of the domain ]0,1] are successively refined.
= ]0,1] this section investigates numerically the weak stochastic integral , for which the Riemann sum estimates do not generally converge whenever the partitions of the domain ]0,1] are successively refined.
Section 8.4 based the meaning of weak stochastic integrals on rth binary partitions of ]0,1] and the corresponding rth binary Riemann sum contingent observables given by the sum of the squares of the rth binary Brownian increments.
For any given r the rth binary Riemann sum contingent observable has elementary form given by
The potential data-values 
are unbounded above, and have zero as lower bound. This holds for all r. Therefore, unlike , none of the observables is normally distributed. (The potential data-values in a normal distribution are unbounded above and below.)
are unbounded above, and have zero as lower bound. This holds for all r. Therefore, unlike , none of the observables is normally distributed. (The potential data-values in a normal distribution are unbounded above and below.)
Figure 9.5: Histogram for 10,000 sums of 100 normal increments.
Loosely speaking, the idea behind the weak stochastic integral is as follows. As r → ∞ the rth binary Riemann sums may converge for some, but not all, of the potential joint outcomes of the Brownian motion . But when the set of rth binary Riemann sums is formed for each 
, the Maple calculations below show that, for integrand 
, these sets tend to cluster ever more closely around some target value as r increases to infinity. This “target value” is the weak stochastic integral of g2 on .
of the Brownian motion . But when the set of rth binary Riemann sums is formed for each , the Maple calculations below show that, for integrand , these sets tend to cluster ever more closely around some target value as r increases to infinity. This “target value” is the weak stochastic integral of g2 on .
The test for weak convergence of the rth binary Riemann sums of any stochastic integrand g(XS) is as follows.
- Find the “target value”; that is, the observable k() which is to be candidate for the weak stochastic integral of g(XS). (In the case of
, the candidate is the constant 1.) - For each r calculate the difference kr() between rth binary Riemann sum of g(XS) and k():
- Determine whether
→ 0 as r → ∞. If it does, then k() is the weak stochastic integral of 
on .
The function 
is an outer measure—a kind of volume—defined on subsets S of 
. If Sr is the set 
, then continuity of G implies
is an outer measure—a kind of volume—defined on subsets S of . If Sr is the set , then continuity of G implies
provided kr is also continuous; so, under these conditions, 
→ 0 implies that the “volume” 
→ 0. The aim of this section is to produce numbers and diagrams which illustrate this.
→ 0 implies that the “volume” → 0. The aim of this section is to produce numbers and diagrams which illustrate this.
The distribution function of the contingent rth binary Riemann observable 
is the Brownian distribution function G. What about 
, the distribution function of the elementary form 
, with z(r) = 
is the Brownian distribution function G. What about , the distribution function of the elementary form , with z(r) =
Denote the inverse function of 
by Q. Provided g is G-measurable, Theorem 227 then gives
by Q. Provided g is G-measurable, Theorem 227 then gives
If the stochastic integrand g is g1 = X(
) = dXs, then, for every r, is simply the standard normal distribution function N.
) = dXs, then, for every r, is simply the standard normal distribution function N.
Other stochastic integrands g may be dealt with on a case by case basis, as in Chapter 8. But without striving3 to find analytical form for distribution functions , numerical calculations can provide a sense of them.
, numerical calculations can provide a sense of them.
In the calculations below it is more convenient to use decimal rather than rth binary partitions. The Maple code in Calculation 5 is intended to give 1000 Riemann sum estimates, each containing 25 terms, of the weak stochastic integral of 
on ]0,1].
on ]0,1].
Calculation 5 selects 25,000 random Brownian increments, each having standard deviation 0.2 and variance 0.04. With t = 1 this is equivalent to partitioning the domain = ]0, t] into 25 equal steps of length 0.04. In each case Maple calculates the Riemann sum of integrand
= ]0, t] into 25 equal steps of length 0.04. In each case Maple calculates the Riemann sum of integrand
so the sample of 25,000 increments generates a sample of 1000 Riemann sums, each containing 25 terms. (A Riemann sum containing 25 terms is not binary—the closest binary partition has r = 5, giving a 25 = 32 term Riemann sum. Either way, the Riemann sum calculation gives an idea of how the resulting values are distributed.) Maple then constructs the histogram, Figure 9.6, of the 25,000 Riemann sum values.
restart:
with(Statistics):
BrownianIncrements := Sample(Normal(0, 0.2), 25000):
dx := [seq(BrownianIncrements[i], i=1..25000)]:
for r from 1 to 1000 do
k := 25*r:
j := k − 24:
RiemannSum[r] := add(dx[i]2, i=j..k) :
end do:
RiemannSums := [seq(RiemannSum[r], r=1..1000)] :
Histogram(RiemannSums);
This histogram in Figure 9.6 gives a sense of the elementary form Z(r) of
Figure 9.6: Histogram of 1000 25-term Riemann sums 
.
.
which are, respectively,
with r = 5. (Taking r = 5 requires partitioning into 32 parts, while, for convenience, the partitioning in Calculation 5 and Figure 9.6 involves only 25 terms.)
into 32 parts, while, for convenience, the partitioning in Calculation 5 and Figure 9.6 involves only 25 terms.)
The observable Z(r) does not appear to be normally distributed. If the Maple command
Mean(RiemannSums);
is now executed, then, depending on the particular Riemann sum values produced by the random sample or simulation of Calculation 5, a result such as 0.9930280578 is obtained. Since this is close to 1, it may constitute some kind of confirmation that (9.1) converges in some sense to 1 for this particular stochastic integrand.
The values of dx[i], given by the random sample in Calculation 5, are normally distributed Brownian increments. Subject to Maple programming constraints, it is hoped that they are also independent. Therefore the elementary-form distribution function indicated in Figure 9.6 corresponds to the distribution function given by (9.2). In other words, Z(r) “inherits”, in some sense, its distribution function FZr in domain R, from the standard Brownian distribution function G of the standard Brownian joint-basic process in .
indicated in Figure 9.6 corresponds to the distribution function given by (9.2). In other words, Z(r) “inherits”, in some sense, its distribution function FZr in domain R, from the standard Brownian distribution function G of the standard Brownian joint-basic process in .
The following calculations use successively finer partitions of = ]0,1] in order to see the shape of the distribution functions as r increases. For convenience, binary partitions are avoided.
= ]0,1] in order to see the shape of the distribution functions as r increases. For convenience, binary partitions are avoided.restart:
with(Statistics):
BrownianIncrements := Sample(Normal(0, 0.1), 100000):
dx := [seq(BrownianIncrements[i], i=1..100000)]:
for r from 1 to 1000 do
k := 100*r:
j := k − 99:
RiemannSum[r] := add(dx[i]2, i=j..k) :
end do:
RiemannSums := [seq(RiemannSum[r], r=1..1000)]:
Histogram(RiemannSums);
Calculation 6 produces the histogram in Figure 9.7. It indicates a clustering of Riemann sum values around a “target” value 1. The range of Riemann sum values obtained in this case is less than the range demonstrated in Figure 9.6 for Calculation 5, indicating “tighter” clustering as r increases.
Calculation 7 has 1000 Riemann sums 
, each containing 400 squares of Brownian increments dxs, each of which in turn has mean zero with standard deviation 0.05 and variance 0.0025; corresponding to random walks of 400 steps in the domain ]0,1].
, each containing 400 squares of Brownian increments dxs, each of which in turn has mean zero with standard deviation 0.05 and variance 0.0025; corresponding to random walks of 400 steps in the domain ]0,1].
Figure 9.7: Histogram of 1000 100-Riemann sums 
.
.
Figure 9.8: Histogram of 1000 400-term Riemann sums .
.restart:
with(Statistics):
BrownianIncrements := Sample(Normal(0, 0.05), 400000):
dx := [seq(BrownianIncrements[i], i=1..400000)]:
for r from 1 to 1000 do
k := 400*r:
j := k − 399:
RiemannSum[r] := add(dx[i]2, i=j..k) :
end do:
RiemannSums := [seq(RiemannSum[r], r=1..1000)] :
Histogram(RiemannSums);
The output from Calculation 7 is the histogram in Figure 9.8.
The three histograms in Figures 9.6, 9.7, and 9.8 can be seen to provide some empirical support for (9.1) converging to zero as r → ∞. They provide some confirmation that the integrand g = g2 = dX2 has (weak) stochastic integral f() = 1.
) = 1.
To test further for weak convergence of the Riemann sums to 1, change the preceding Maple calculations to calculate
instead of
RiemannSum.restart:
with(Statistics):
BrownianIncrements := Sample(Normal(0, 0.2), 25000):
dx := [seq(BrownianIncrements[i], i=1..25000)]:
for r from 1 to 1000 do
k := 25*r:
j := k − 24:
VarSum[r] := abs(add(dx[i]2, i=j..k ) - 1):
end do:
VarSums := [seq(VarSum[r], r=1..1000)]:
Histogram(VarSums);
Mean(VarSums);restart:
with(Statistics):
BrownianIncrements := Sample(Normal(0, 0.1), 100000):
dx := [seq(BrownianIncrements[i], i=1..100000)]:
for r from 1 to 1000 do
k := 100*r:
j := k − 99:
VarSum[r] := abs(add(dx[i]2, i=j..k − 1):
end do:
VarSums := [seq(VarSum[r], r=1..1000)] :
Histogram(VarSums);
Mean(VarSums);
Figure 9.9: Histogram for 25-term Riemann sum observable with 
= 0.04.
= 0.04.restart:
with(Statistics) :
BrownianIncrements := Sample(Normal(0, 0.05), 400000):
dx : = [seq(BrownianIncrements[i], i=1..400000)] :
for r from 1 to 1000 do
k := 400*r:
j := k − 399:
VarSum[r] := abs(add(dx[i]2, i=j..k) − 1)
end do:
VarSums := [seq(VarSum[r], r=1..1000)]:
Histogram(VarSums);
Mean(VarSums);
Each of these three (r = 1,2,3) Maple calculations returns two pieces of output. These are, in each case, a histogram of the rth sample, and the mean of the rth sample. For r = 1,2,3, the rth sample consists of 1000 data values
Figure 9.10: Histogram for 100-term Riemann sum observable with = 0.01.
= 0.01.
Figure 9.11: Histogram for 400-term Riemann sum observable, = 0.0025.
= 0.0025.
say. (For convenience of calculation, the partitions of for the Riemann sum are not rth binary.) The rth sample can be thought of as consisting of the 1000 data values returned by 1000 observables:
for the Riemann sum are not rth binary.) The rth sample can be thought of as consisting of the 1000 data values returned by 1000 observables:
being the corresponding elementary and contingent forms.
The terms dx[i] in the Riemann sums are selected by Maple as independent random Brownian increments; so the sample of contingent data elements hr() are distributed in accordance with the Brownian distribution function G and the deterministic function hr. In other words, the elementary data-values y(r) = hr() inherit from G and hr a distribution function 
in R. Therefore, using (9.2), for any J I(R) the proportion of sample data-values y in J is approximately
) are distributed in accordance with the Brownian distribution function G and the deterministic function hr. In other words, the elementary data-values y(r) = hr() inherit from G and hr a distribution function in R. Therefore, using (9.2), for any J I(R) the proportion of sample data-values y in J is approximately
and, using Theorems 39 and 79,
Each time Calculations 8, 9, and 10 are executed a particular random sample of Brownian increments dx[i] is generated. For a trio of such samples, the three estimates of 
returned by
returned byMean(VarSums);
were
0.2180014861, 0.1101283834, 0.05825507338,
for partitions of using cells of length
using cells of length
0.04, 0.01, 0.0025,
respectively. Though only three terms of the series are examined, the results provide some intuitive support for the statement that 
→ 0 as r → ∞.
→ 0 as r → ∞.9.5 Calculation of Itô’s Formula
For σ R, = ]0,i], and 
, Itô’s formula states that
R, = ]0,i], and , Itô’s formula states that
the final integral being weak stochastic. This can be written as
The “=” sign in Itô’s formula amounts to a statement about the value of the weak stochastic integral on the right-hand side and is therefore a statement about weak convergence of Riemann sums.
In Section 8.12 the formula was verified for the function
For this particular function f, the objective here is to demonstrate numerically and visually the meaning of the weak convergence of Riemann sum estimates of the observables in Itô’s formula. To simplify the Maple calculations take μ = 0.
As before, the calculations are done by constructing large samples of the potential values of the Riemann sum observables and drawing histograms of these samples. The process is repeated, using finer partitions in each case, and the resulting histograms illustrate successively closer clustering of the sample Riemann sum values around a weak limit.
The intention is to provide some numerical and visual sense of what is involved in Itô’s formula, as expressed in (8.36). That is,
where, with integrands from (9.3),
The following Maple code is used to estimate 
.
.restart;
with(Statistics);
BrownianIncrements := Sample(Normal(0, 10), 10000);
dx := [seq(BrownianIncrements[i], i=1..10000)];
for q from 1 to 1000 do
k := 10*q:
j := k − 9:
s[j − 1] := 0:
xs[j − 1] := 0:
RS[j − 1] := 0:
for i from j to k do
s[i] := 100*i:
xs[i] := add(dx[q], q =j..i):
hl[i] := 10*exp(10*xs[i] − (1/2*100)*s[i]
− 10*exp(10*xs[i − 1] − (1/2*100)*s[i − 1]) :
h2[i] := 10*exp(10*xs[i − 1]
− (1/2*100)*s[i − 1]) * (xs [i] − xs[i − 1]):
RS[i] := RS[i − 1] + abs(hl[i] − h2[i]):
end do:
IT0[q] := RS[k]:
end do:
Ito := [seq(IT0[q], q=1..1000)]:
Histogram(Ito);
Mean(Ito);
This simulation uses standard Brownian increments dx = xs′ − xs where
Sample (Normal (0, 10), 10000) selects a random sample of 10,000 “independent” normal increments. The domain = ]0,1000] is partitioned into 10 equal intervals of length σ2 = 100, with partition points Si = s[i]. The process value x(si) = xs[i] is obtained by adding up i normal increments dx. The values of h1 and h2 at partition point Si are given by hi[i] and h2[i] in accordance with (9.4). Calculation 11 has r = 10 (so the partition is not binary). For any k, 1 ≤ k ≤ 1000, the kth sample Riemann sum value RS[k] of 
is given by
which the Maple code finds by accumulating, in successive cycles, the absolute value of terms (
hi [i]− h2 [i]) in the r-partition:RS[i] := RS[i − 1] + abs(hl[i] − h2[i]).
A sample of 1000 of these Riemann sum estimates of 
is listed in the Maple list denoted by
is listed in the Maple list denoted by Ito.
The histogram (Figure 9.12) of the list
Ito gives an impression of the Riemann sum observable for this instance of Itô’s formula.
As explained in Section 9.4, the 
for a particular integer r. According to the theory,
Mean(Ito) command provides a Riemann sum estimate of for a particular integer r. According to the theory,
as r → ∞. With r = 10 in = ]0,1000], a sample of 1000 Riemann sum estimates of produced for Calculation 11 returned an estimated value
= ]0,1000], a sample of 1000 Riemann sum estimates of produced for Calculation 11 returned an estimated value Mean(Ito) = 76.86002866.
To see what happens when = ]0,1000] is partitioned with smaller subintervals, Calculation 12 selects 1,000,000 Brownian increments with σ = σ2 = 1. A sample of 1000 Riemann sums is formed from partitions of = ]0,1000] consisting of 1000 equal intervals of length 1.
= ]0,1000] is partitioned with smaller subintervals, Calculation 12 selects 1,000,000 Brownian increments with σ = σ2 = 1. A sample of 1000 Riemann sums is formed from partitions of = ]0,1000] consisting of 1000 equal intervals of length 1.restart;
with(Statistics);
BrownianIncrements := Sample(Normal(0, 1), 1000000);
dx := [seq(BrownianIncrements[i], i = 1..1000000)] ;
for q from 1 to 1000 do
k := 1000*q:
j := k − 999:
s[j − 1] := 0:
xs[j − 1] := 0:
RS[j − 1] := 0:
for i from j to k do
s[i] := 100*i:
xs[i] := add(dx[q], q=j..i):
h1[i] := 1*exp(1*xs[i] − (
*1)*s[i]
− 1*exp(1*xs[i − 1] − (*1)*s[i − 1]):
h2[i] := 1*exp(1*xs[i − 1]
− (*1)*s[i − 1])*(xs[i] − xs[i − 1]):
RS[i] := RS[i − 1] + abs(h1[i] − h2[i]):
end do:
IT0[q] := RS[k]:
end do:
Ito := [seq(IT0[q], q-1..1000)] :
Histogram(Ito);
Mean(Ito);
The histogram representation of the elementary form of the Riemann sum observable is in Figure 9.13. Mean(Ito) returned 1.152977886 as estimate of 
.
.
The following two calculations are done on domain = ]0,1]. Calculation 13 has σ = 0.2, so = ]0,1] is partitioned into 1/0.04 = 25 subintervals; Maple simulates 1000 Riemann sum estimates for (9.4), with histogram output in Figure 9.14, and Mean(Ito) estimate 0.2185890853 for VhrG[RT].
= ]0,1]. Calculation 13 has σ = 0.2, so = ]0,1] is partitioned into 1/0.04 = 25 subintervals; Maple simulates 1000 Riemann sum estimates for (9.4), with histogram output in Figure 9.14, and Mean(Ito) estimate 0.2185890853 for VhrG[RT].restart;
with(Statistics);
BrownianIncrements := Sample(Normal(0, 0.2), 25000);
dx := [seq(BrownianIncrements[i], i=1..25000)];
for q from 1 to 1000 do
k := 25*q:
j := k − 24:
s[j − 1] := 0:
xs[j − 1] := 0:
RS[j − 1] := 0:
for i from j to k do
s[i] := 0.04*i:
xs[i] : = add(dx[q], q=j..i) :
h1[i ] := 0.2*exp(0.2*xs[i] − (*0.04)*s[i]
− 0.2*exp(0.2*xs[i − 1] − (*0.04)*s[i − 1]) :
h2[i] := 0.2*exp(0.2*xs[i − 1]
− (*0.04)*s[i − 1])*(xs[i] − xs[ i − 1]) :
RS[i] := RS[i − 1] + abs(h1[i] − h2[i]) :
end do:
ITO[q] := RS[k]:
end do:
Ito := [seq(IT0[q], q=1..1000)]:
Histogram(Ito);
Mean(Ito);
Figure 9.12: Histogram for Itô formula, 10-term approximation, = ]0,1000].
= ]0,1000].
Figure 9.13: Itô formula, 1000-term approximation, =]0,1000].
=]0,1000].
Calculation 14 has σ = 0.05, so = ]0,1] is partitioned into 1/0.0025 = 400 subintervals, and Maple simulates 1000 Riemann sum estimates for (9.4), with histogram output in Figure 9.15, and Mean(Ito)estimate 0.1893459552 for VhrG [RT].
= ]0,1] is partitioned into 1/0.0025 = 400 subintervals, and Maple simulates 1000 Riemann sum estimates for (9.4), with histogram output in Figure 9.15, and Mean(Ito)estimate 0.1893459552 for VhrG [RT].restart;
with(Statistics);
BrownianIncrements := Sample(Normal(0, 0.05), 400000);
dx := [seq(BrownianIncrements[i], i=1..400000)] ;
for q from 1 to 1000 do
k := 400*q:
j := k − 399:
s[j − 1] := 0:
xs[j − 1] := 0:
RS[j − 1] := 0:
for i from j to k do
s[i] := 100*i:
xs[i] := add(dx[q], q=j..i) :
h1[i] := 0.05*exp(0.05*xs[i] − (*0.0025)*s[i]
− 0.05*exp(0.05*xs[i − 1] − (*0.0025)*s[i − 1] ) :
h2[i] := 0.05*exp(0.05*xs[i − 1]
− (*0.0025)*s[i − 1])*(xs[i] − xs[i − 1]) :
RS[i] := RS[i − 1] + abs(h1[i] − h2[i]) :
end do:
IT0[q] := RS[k]:
end do:
Ito := [seq(IT0[q], q=1..1000)]:
Histogram(Ito);
Mean(Ito);
Figure 9.14: Histogram for Itô formula, 25-term approximation, = ]0,1].
= ]0,1].
Figure 9.15: Histogram for Itô formula, 400-term approximation, = ]0,1].
= ]0,1].9.6 Calculating with Binary Partitions of RT
The calculations in preceding sections have been done, essentially, for elementary-form observables Y≃ y[R, FY], a form which lends itself to visualization by histograms. In contrast, a joint-contingent observable with sample space RT has representation
f(XT) ≃ f(xT) [RT, FXT],
and its expected value is
E[f(XT)] = 
f(xT)FxT(I[N]).
f(xT)FxT(I[N]).
This calculation is done by forming partitions and Riemann sums in RT.
Binary (r, q)-partitions of infinite-dimensional domains RT are described in Chapter 3, Section 3.5. With T = ]0, t], 2r is the number of partition points 
= j2−r of T, and 
+ 1 is the number of partition points k2−q of 
for each j. As detailed in (3.12), the number of cells in an (r, q)-partition 
rq of RT is 
, which increases greatly as r and q increase. Since this is also the number of terms in the corresponding Riemann sum, the scale and amount of calculation involved quickly escalate as r and q increase.
= j2−r of T, and + 1 is the number of partition points k2−q of for each j. As detailed in (3.12), the number of cells in an (r, q)-partition rq of RT is , which increases greatly as r and q increase. Since this is also the number of terms in the corresponding Riemann sum, the scale and amount of calculation involved quickly escalate as r and q increase.
Take T = ]0,1] and FXT = G, so the process is standard Brownian motion on ]0,1]. The binary partition points of T = ]0,1] are
In Section 7.9 a binary cell K partitioning RT is denoted by
and the value of G on K is denoted by
where, for c = −,
,
the lower and upper limits of integration (ql[j] and qu[j], respectively) being drawn from the cells Kq|kwhich, for each j, partition 
:
:
, j = 1, 2, 3, 4; that is,
With q = 2, then, in accordance with (9.6), and for each j, the one-dimensional domain 
has partition points
has partition points
For 1 ≤ j ≤ 4, the lower and upper integration limits in (9.5), ql[j] and qu[j], respectively, are drawn from (9.7). Therefore, for each j, the corresponding integral in (9.5) can be any one of
By separating (9.7) into two separate lists,
the following Maple code establishes limits ql[j] and qu[j] in domain R]0, 1] for the iterated integrals in (9.5).
restart :
r := 2:
q := 2:
for j from 1 to 2r do S[j] := j*2−r: end do;
s := [seq(S[j], j = 1 .. 2r)];
kmax := q*2(q + 1):
for k from 1 to kmax do qList[k] := − q + (k − 1)*2−q : end do:
q1 := [−∞, seq(qList [k], k = 1 .. kmax), q];
qu := [seq(qList[k], k = 1 .. kmax), q, ∞];
− 
= 2
Then
G(K) = βα.
This is calculated in the next section.
9.7 Calculation of Observable Process in RT
Section 7.10 shows how to use binary partitions rq to estimate the expectation of a joint-contingent observable
rq to estimate the expectation of a joint-contingent observable
f(XT) ≃ f(xT)[RT, G].
In Theorem 207, a Riemann sum estimate for E[f(XT)] is
where f(rq)(xT) is a step function approximation of f(xT) This section calculates a single term of this Riemann sum. Other terms are calculated similarly.
To simplify the Maple calculation, instead of a binary (r, q)-partition suppose T = ]0,1] is partitioned as
and, at each of the partition points of T, suppose R is partitioned as
] − 10, −0.6], ] − 0.6, −0.3], ] − 0.3, 0], ]0,0.3], ]0.3,0.6], ]0.6, 10]
in order to calculate the integral in (9.9). (With slight re-wording, Theorem 209 ensures convergence of the step functions for any regular partition, not just the binary ones.)
Lower and upper limits of −10 and 10 are used here in place of −∞ and ∞. This makes little difference to the calculated values of the Brownian distribution function G(I[N]), as the variance of the Brownian increments is
The corresponding partition of RT consists of 63 = 216 cells K; for example,
Then, with y0 = 0 and
Maple can calculate the distribution function value G(K) = βα, as follows.
restart:
y[0] := 0:
g := mul(exp(− (3/2)*(y[j] − y[j −1])2), j=1.,3) :
α := evalf10(int(g, y[1]=0.. 0.3, y[2] =−0.3 ..0, y[3]=0.3..0.6)):
β := evalf10((2*Pi)
* 
):
G := β * α;
0.004336025795
Thus G(K) = 0.004336025795. The Brownian likelihood of other cells can be similarly calculated. For any cell K of the regular partition, any one of the three lower limits of integration in Line 4 of Calculation 16 is any one of the six numbers
10 − 0.6 − 0.3 0 0.3 0.6;
and in each case, the corresponding upper limit of integration in Line 4 of Calculation 16 is
−0.6 −0.3 0 0.3 0.6 10.
Once the lower limits of the integration are established, the corresponding upper limits are then determined.
To make a list of the lower limits, it is necessary to list all possible permutations, with repetition, of the numbers 10, −0.6, −0.3, 0, 0.3, 0.6. Maple has a permute command which lists all permutations of nitems, taken r at a time without repetition. The following sequence of commands will generate the listings needed for permutations of six items taken three at a time, with repetition.
restart;
l[1] := − 10: l[2] := − 0.6: l[3] := − 0.3:
l[4] := 0: l[5] := 0.3: l[6] := 0.6:
u[l] := − 0.6: u[2] := − 0.3: u[3] := 0:
u[4] := 0.3: u[5] := 0.6: u[6] := 10:
L := [seq(l[j], j = 1 .. 6)];
[− 10, − 0.6, − 0.3, 0, 0.3, 0.6]
U := [seq(u[j], j = 1 .. 6)];
[− 0.6, − 0.3, 0, 0.3, 0.6, 10]
p := 0:
for rl from 1 to 6 do
for r2 from 1 to 6 do
for r3 from 1 to 6 do
p :- p + 1:
P[p] := [L[r1], L[r2], L[r3]]:
Q[p] := [U[r1], U[r2], U[r3]]:
end do:
end do:
end do:
This generates a list of 216 lower integration limits L, and the list Q of corresponding upper integration limits. For instance, the 100th cell K of the 216 member partition of R]0, 1] is determined by the lists
P[100] = [−0.3, −0.3,0.3], Q[100] = [0,0,0.6],
so
K[100] = ] −0.3, 0] × ] −0.3, 0] × ]0.3, 0.6] × 
.
.
A Riemann sum estimate of the expected value of the observable f(XT) is then given by
the sum being taken over the 216-member partition. Suppose
With the continuous modification of G, f(xT) can be taken to be zero if xT is not continuous. Theorem 209 ensures that E [f(XT)] exists, and that it can be approximated to any degree of accuracy by Riemann sums over regular partitions of RT. For the cell K in Calculation 16, Theorem 214 implies that f(xT) can be taken to be
since x0 = 0, and the values x1 and x2 can be taken to be 0.3 and 0. This reduces to exp(−0.03), or 0.9704455335. Therefore this particular term of the Riemann sum reduces to
exp(−0.03)G(K) = 0.9704455335 × 0.004336025795 = 0.004207876866.
If this calculation is replicated for the other 215 cells partitioning R]0, 1], then the resulting Riemann sum calculation gives an estimated value for E [f(XT)]. The same method can be used to obtain estimates of the diffusion function ψ(ξ, 
), given by the marginal density of expectation
), given by the marginal density of expectation
Eξτ [f(XT)].
In this case = 1. Suppose ξ = 0.5. Then, in Calculation 16, take
= 1. Suppose ξ = 0.5. Then, in Calculation 16, take
y[3] = ξ = 0.5,
and perform the integration on variables y[l] and y[2] only. The result given by Maple in this case is
G(K−)= 0.03999644118,
and the corresponding Riemann sum term is exp(−0.03)G(K−) =
The approximate numerical value of E0.5,1(f(XT)) = ψ(0.5, 1) can be built up from calculations like this one.
9.8 Other Joint-Contingent Observables
As further illustration of numerical calculation in sample space RT, suppose the Brownian observable is
f(XT) = sup{x(t),t ]0, 1] = T}.
]0, 1] = T}.
This function of xT is not continuous as x(t)t T varies continuously in its separate components. Nor is this function bounded. Therefore Theorem 209 does not guarantee convergence of its Riemann sums over regular partitions, and further investigation into convergence is required.
T varies continuously in its separate components. Nor is this function bounded. Therefore Theorem 209 does not guarantee convergence of its Riemann sums over regular partitions, and further investigation into convergence is required.
For similar reasons it is not possible to apply Theorem 214, as it stands, in order—as in (9.10)—to produce a version of f which is more amenable to numerical calculation. But if conditions were such that devices like these could be applied, then, for the cell K of Calculation 16, this particular term of the Riemann sum estimate of E[f(XT] has a value of 0.3 or 0.6—depending on the method of calculation—for the potential data-value f(xT) for that particular cell K. Potential data-values f(xT) for other cells can be obtained similarly.
Much of this book is concerned with potentiality distribution functions FX which can take negative or complex values, such as the Feynman distribution function G
, with 
= 
. To perform computations such as those in Calculation 16 and (9.11), the integrand can be separated into real and imaginary parts. For the cell K of Calculation 16, with T = ]0, 1], the real part 
(G (K)) is α1 + α2, where
, with = . To perform computations such as those in Calculation 16 and (9.11), the integrand can be separated into real and imaginary parts. For the cell K of Calculation 16, with T = ]0, 1], the real part (G (K)) is α1 + α2, where
This is obtained by separating
into real and imaginary parts, with n = 3 and tj − tj−1 = 
. These iterated integrals can be calculated as in Calculation 16, with the appropriate trigonometric function replacing the exponential function. As in (9.11), a term of the Riemann sum estimate of E [
(XT)] can be calculated; and again the potential datum-value (xT) for each such term can be separated into real and imaginary parts.
. These iterated integrals can be calculated as in Calculation 16, with the appropriate trigonometric function replacing the exponential function. As in (9.11), a term of the Riemann sum estimate of E [(XT)] can be calculated; and again the potential datum-value (xT) for each such term can be separated into real and imaginary parts.
For any given ξ R and R+, Riemann sum estimates of the quantum mechanical state function ψ(ξ, ) can similarly be obtained, by removing the integration on ξ = yn as in Calculation 16. Likewise, Riemann sum estimates of the values of terms ψ, p(ξ, ) of the series expansion of ψ(ξ, ) on which Feynman diagrams are based—see Section 7.22. Theorem 222 then provides convergence conditions for these Riemann sums over regular partitions such as those in Calculation 16.
R and R+, Riemann sum estimates of the quantum mechanical state function ψ(ξ, ) can similarly be obtained, by removing the integration on ξ = yn as in Calculation 16. Likewise, Riemann sum estimates of the values of terms ψ, p(ξ, ) of the series expansion of ψ(ξ, ) on which Feynman diagrams are based—see Section 7.22. Theorem 222 then provides convergence conditions for these Riemann sums over regular partitions such as those in Calculation 16.9.9 Empirical Data
At the start of this book Table 1.1 gives simple numerical data from which an estimate of the distribution function for the observable is deduced. No claim is made regarding the accuracy or reliability of the estimate.
Is it possible to estimate a joint distribution function from sample data xT of a joint observable or random process XT? And why should the attempt be made? This section seeks to provide support for empirical approaches.
In contrast, the basic assumption of Section 8.16 is (8.48). This states that asset price processes ZTsatisfy ZT ≃ zT 
so the increments ln Zt – ln Zt′ are assumed to be
so the increments ln Zt – ln Zt′ are assumed to be- independent, and
- normally distributed.
In Ross [198], the end-of-day price of a barrel of crude oil, as traded on the New York Mercantile Exchange, is listed for each trading day from 3 January 1995 to 19 November 1997. A chi-square statistical test of the independence hypothesis is applied to increments ln Zt − ln Zt′, and the independence hypothesis fails.
In general, independence is intuitively implausible. Realistically it is hard to conceive of any pair of real-world events for which there is no possible connecting chain of causal events, however remote, linking the two together.
Does this mean that the financial theory presented in Section 8.16 is worthless? Independence is a useful mathematical abstraction. The fact that, in geometry, the mathematically perfect circular shape rarely if ever actually occurs in physical reality does not detract from the crucial role of the mathematical circle in understanding the world. But it is unwise, in the “real” world, to gamble on occurrence of a mathematically perfect abstraction.
The assumption that the increments of the logarithms of asset prices are normally distributed has also been challenged. It is argued that large changes in asset prices are much more common than the theory of Section 8.16 predicts.
This section presents series of asset prices from4 Thomson Reuters Datastream [224] in order to assess these issues further. Datastream provides data in Excel format. To import an Excel file FileName. xls into Maple the following Maple commands should be executed:
with(ExcelTools):
Import(FileLocation/FileName.xls):
Figure 9.16 is a Maple graph displaying a series of daily end-of-day share prices of Glanbia, a food processing company. Prices are in pounds sterling. The series consists of 5219 terms, running from 8 March 1991 to 9 March 2011, beginning and ending as follows:
0.8, 0.817, 0.836, ···, 3.782, 3.699, 3.692.
Figure 9.16: Twenty-year graph of Glanbia daily share prices.
Take the logarithm of each of the daily share prices, and calculate the daily increments by subtraction. This gives a list of 5218 log-increments, with graph Figure 9.17 and histogram Figure 9.18. Maple returns a mean value μ = 0.0002930839150 for the 5218 increments of the logarithms of the share prices, with standard deviation σ = 0.0234953419046857.
If the Glanbia share prices satisfy the assumptions of Section 8.16, Figure 9.18 should approximate to a normal distribution with mean zero and standard deviation 0.0235.
Does Figure 9.18 look like a histogram of a normally distributed random variable? That is what is assumed in Section 8.16.
A chi-square statistical test can be applied to the data to check whether the Glanbia data satisfy the assumptions of Section 8.16. The null hypothesis is that the sample of 5218 data items of Figure 9.18are observations of a normally distributed random variable with mean μ and standard deviation σ. Then the expected frequencies e are obtained, in accordance with the null hypothesis, by calculating
Figure 9.17: Graph of log-increments of Glanbia share prices.
Figure 9.18: Histogram of log-increments of Glanbia share prices.
A simple Maple calculation gives the actual, or observed, frequencies o of the Glanbia data. The observed and expected frequencies are given in Table 9.1.
The tabulation confirms the impression given by Figure 9.16, that extreme changes in share prices are more common than what is predicted under the lognormality assumption.
The chi-square value is
With seven degrees of freedom, the critical values of χ2 for significance levels 0.05, 0.01, and 0.001 are 14.07, 18.48, and 24.32, respectively. Therefore the null hypothesis should be rejected. The evidence of the Glanbia data contradicts the hypotheses of Section 8.16. The Glanbia data demonstrate the “black swan” phenomenon—events which are rare in theory are common in practice.
Figures 9.19 to 9.24 refer to daily series of share price data, each series of twenty years’ duration, selected miscellaneously from Thomson Reuters Datastream [224].
Figures 9.21 and 9.24 show significant numbers of extreme outlying values. Superficial inspection of this miscellaneous selection of share price data indicates that the hypotheses of Section 8.16 cannot be assumed to be valid. In particular, it cannot be assumed, without verification, that a geometric Brownian distribution function 
is applicable to series of share prices. And even if some modified version of is used, it is still necessary to verify, test, or otherwise ensure that the distribution function fits the data.
is applicable to series of share prices. And even if some modified version of is used, it is still necessary to verify, test, or otherwise ensure that the distribution function fits the data.
Figure 9.19: Graph of Rio Tinto share prices.
Figure 9.20: Graph of log-increments of Rio Tinto share prices.
Figure 9.21: Histogram of log-increments of Rio Tinto share prices.
Figure 9.22: Graph of GlaxoSmithKline share prices.
Figure 9.23: Graph of log-increments of GlaxoSmithKline share prices.
Figure 9.24: Histogram of log-increments of GlaxoSmithKline prices.
9.10 Empirical Distributions
Even if the finance theory of Section 8.16 cannot be relied on, it may still be possible to analyze the random variability properties of empirical data.
Section 8.16 is predicated on a particular analytic distribution function . Various modifications of this function are sometimes invoked to deal with the mismatch between the theory and the empirical data such as those in Section 9.9. But it is difficult to imagine a reason why the likelihood distribution—if it exists—of any series of empirical data must have some analytic representation. There seems to be no a priori justification for assuming that there is an “algebraic formula” for such likelihoods, nor that any algebraic formula will give accurate values of the joint likelihoods.
. Various modifications of this function are sometimes invoked to deal with the mismatch between the theory and the empirical data such as those in Section 9.9. But it is difficult to imagine a reason why the likelihood distribution—if it exists—of any series of empirical data must have some analytic representation. There seems to be no a priori justification for assuming that there is an “algebraic formula” for such likelihoods, nor that any algebraic formula will give accurate values of the joint likelihoods.
Does this mean that it is impossible to carry out accurate analysis of the random variability of such data?
If there is no underlying likelihood 
operating on the joint data series, then the series is not an observable, and the methods of this book have nothing to offer in the form of guidance. The entrails of a goat would be just as good.
operating on the joint data series, then the series is not an observable, and the methods of this book have nothing to offer in the form of guidance. The entrails of a goat would be just as good.
But if there is an underlying joint likelihood , then it can be estimated numerically, and the likelihood estimates can then be used, as in Sections 9.7 and 9.8, to estimate expected values, marginal expectation densities, and the like—just as we sought to use the distribution function . The accuracy of the results will then depend on the accuracy of the estimates of .
, then it can be estimated numerically, and the likelihood estimates can then be used, as in Sections 9.7 and 9.8, to estimate expected values, marginal expectation densities, and the like—just as we sought to use the distribution function . The accuracy of the results will then depend on the accuracy of the estimates of .
In Brownian motion there is a theoretical foundation—such as Theorem 217—which gives some guarantee of the ultimate accuracy of the results of estimates based on regular and binary partitions of RT. Such estimates are good when the processes involved are Brownian or geometric Brownian, and when the contingent observable integrand function f satisfies the conditions of Theorem 217. But what guarantee of accuracy is there when these assumptions are dropped?
Brownian and geometric Brownian processes (Xt)t T are defined for continuous t, and the data-values xt are continuous and unbounded. Therefore results such as Theorem 217 are needed as a basis for discrete-time, discrete-value estimates from regular or binary partitions of RT.
T are defined for continuous t, and the data-values xt are continuous and unbounded. Therefore results such as Theorem 217 are needed as a basis for discrete-time, discrete-value estimates from regular or binary partitions of RT.
But the empirical price processes (Xt) of Section 9.9 do not have continuous t. As recorded in Thomson Reuters Datastream [224], t is discrete, and each datum value xt is constant from one end-of-day observation to the next. And even if the constantly changing moment-to-moment prices are taken into account, these too are momentarily constant from one “time-tick” to the next—unlike the values taken by a Brownian or geometric Brownian process.
Also, the prices xt are discrete, not continuous, since these are units of money and are not infinitely divisible. The value of xt is bounded below by zero, and while theoretically there is no upper bound, in practice some very large upper cutoff value can be applied. So in any finite interval of time there is only a finite number of time-ticks t, and there is, in effect, only a finite number of possible values xt.
In other words, while the insights and perspectives of advanced mathematics can be helpful, investigation of the random variability of share prices is inherently simpler than Brownian motion and quantum mechanics. In fact, results like Theorem 217 are irrelevant if, when analyzing price processes, the false assumption of continuous time/continuous values is abandoned, and the more realistic discrete view is adopted.
Thus, if numerical estimates of joint likelihood are used, the accuracy of any final calculation depends on the accuracy of the estimates. In other words, accuracy of results depends on the amount of effort, insight, and skill used in estimating joint likelihoods.
are used, the accuracy of any final calculation depends on the accuracy of the estimates. In other words, accuracy of results depends on the amount of effort, insight, and skill used in estimating joint likelihoods.9.11 Calculation of Empirical Distribution
Though detailed investigation of such estimates of joint likelihood is beyond the scope of this book, some rudimentary calculation is presented here for the Glanbia series from Thomson Reuters Datastream [224].
The series starts on Friday 8 March 1991, with Glanbia share price standing at 0.8 at the end of that day, and concludes on Wednesday 9 March 2011, with share price 3.692. Let ′, ) denote the start and finish date, respectively, of this series, and write T = ]′, ].
′, ) denote the start and finish date, respectively, of this series, and write T = ]′, ].
Assuming the existence of joint likelihood potentiality for the Glanbia prices, the series xT from Thomson Reuters Datastream [224] is a joint datum of basic observable
for the Glanbia prices, the series xT from Thomson Reuters Datastream [224] is a joint datum of basic observable
Now suppose Wednesday 9 March 2011 is the present, and suppose a contract in Glanbia shares matures at a point four months in the future (in fact, sixteen weeks later), at end of Wednesday 29 June: in other words, after sixteen five-day trading weeks.
Assuming continued existence of joint likelihood potentiality for the Glanbia prices, let t represent the “future” date Wednesday 29 June 2011 (end-of-day), so = ], t] represents the eighty-day trading period 9 March 2011 to 29 June 2011. Then the Glanbia joint observable for the period is
for the Glanbia prices, let t represent the “future” date Wednesday 29 June 2011 (end-of-day), so = ], t] represents the eighty-day trading period 9 March 2011 to 29 June 2011. Then the Glanbia joint observable for the period is
where, for I[N] I (R), the distribution function values are
I (R), the distribution function values are
To perform some analysis of the joint random variability of the share prices in the forthcoming sixteen weeks, Thursday 10 March 2011 to Wednesday 29 June 2011, estimates of the joint distribution function FX, = FX, are needed.
, are needed.
Suppose is partitioned by times s1, s2, and s3, each time increment being twenty days or four five-day trading weeks. Then s0 = is the present (end of 9 March 2011), s4 = t is 29 June, s1 represents end of Wednesday 6 April 2011; and so on.
is partitioned by times s1, s2, and s3, each time increment being twenty days or four five-day trading weeks. Then s0 = is the present (end of 9 March 2011), s4 = t is 29 June, s1 represents end of Wednesday 6 April 2011; and so on.
Suppose the daily end-of-day prices are partitioned by intervals of length one-tenth of a pound. Then the domain R = R], t] is partitioned by cells
= R], t] is partitioned by cells
where each Ij has the form ]d, d + 0.1]. If marginal densities at time t = s4 are required, then a typical partitioning cell for domain R− = R],t[ is
− = R],t[ is
K− = I1 × I2 × I3 × R\{s1, s2, s3},
\{s1, s2, s3},
where x(t) is now some given, fixed positive number. “Present time” is s = ; and present price is xr = x() = 3.692. Using intervals of one-tenth of a pound, the components of K could be
; and present price is xr = x() = 3.692. Using intervals of one-tenth of a pound, the components of K could be
The next task is to find an estimate of the joint distribution function value FX(K), using the historic joint data-values of the series xT from Thomson Reuters Datastream [224].
One way of making such an estimate is, by examining the series values from 8 March 1991 to 9 March 2011, to find the proportion of occasions when, at four-week intervals, the prices had increments or transitions corresponding to those in the partitioning cell K.
Subtracting 80 from 5219, the number of such tests or comparisons that can be made in the component elements of xT is 5139. Thus, for each of 5139 days, count the number of occasions that, after a further twenty days, the comparison price increases by no more than 0.1, followed by a decrease of no more than 0.1 from the comparison price after 40 days from the comparison date; and so on.
The Maple code in Calculation 18 below evaluates this proportion.
Fgb is the estimated value of FX(K). In this case, fgb = 25, so Maple finds 25 occurrences of joint increments or decrements of the kind appearing in K, and returns
FX(K) = 0.004864759681
as the value of Fgb. Other cells K in the partition can similarly be calculated. A Maple listing of the cells in the partition can be generated by permutation-with-repetition, as in Calculation 17. Once the distribution function estimates FX(K) are found for the cells of the partition, then expected value of a contingent observable f(X) can be estimated; also marginal density of the expectation.
) can be estimated; also marginal density of the expectation.
To obtain a risk-neutral model it is fairly simple to adjust the distribution function estimates in order to simulate, in expectation estimates, a specified growth rate such as the riskless rate of interest ρ.
Calculation 18 is not necessarily suitable5 for every series, and there are many ways in which estimates can be improved. But whatever method is used, care must be taken to ensure that the distribution function values satisfy consistency. The counting procedure in Calculation 18 ensures consistency in this case.
To see how consistency is ensured in this example, suppose dimension s1 is partitioned, with
where the cells 
are disjoint. The relative frequency counting argument ensures that FX (I2 × I3 × I4 × R]0, 1]\{s2, s3, s4}) is given by
are disjoint. The relative frequency counting argument ensures that FX (I2 × I3 × I4 × R]0, 1]\{s2, s3, s4}) is given by
Conditioning is built into Calculation 18, so dubious assumptions of the independence of log-increments kind are averted.
What about other potentially bogus assumptions? It is possible the assumptions (what are they?) in Calculation 18 are invalid, in which case some other estimation method must be tried. But assumptions that might perhaps be valid in some circumstances are preferable to assumptions that are demonstrably false.
for j from 1 to 5139 do
if gb[j + 20] > gb[j] and gb[j + 20] ≤ gb[j] + 0.1
and gb[j + 40] > gb[j] - 0.1 and gb[j + 40] ≤ gb[j]
and gb[j + 60] > gb[j] and gb[j + 60] ≤ gb[j] + 0.1
and gb[j + 80] > gb[j] + 0.1 and gb[j + 80] ≤ gb[j] + 0.2
then fgb := fgb + 1:
end if:
end do:
Fgb := evalf10(fgb/5139);
__________
1 Unfortunately, the corresponding joint-contingent observables f(X) ≃ f(x)[RT, FX] do not have such a convenient and familiar visual representation in multiple dimensions. Figure 1.3 shows the limitations, even when there are only two random variables involved.
3As previously mentioned, the extension of measurability described in Section A.l of the Epilogue may provide a simpler and more efficient way of dealing with this.
5Calculation 18 is, essentially, the relative frequency counting procedure of Table 1.1 with which this book commences.
A.1 Measurability
Without invoking a theory of measure, Chapter 4 demonstrates that Riemann sums, constructed from point/cell/dimension-set association and regulated by gauges, deliver a theory of integration which has the characteristics, properties, and theorems required for the purposes of real analysis, and in particular for the purposes of analysis of random variability.
Sections 4.7 and 4.8 of Chapter 4 introduced variation as a form of outer measure of subsets S of RT. In Section 4.9 this concept was extended from sets S in RT to sets S × 
in RT × 
, and the extended concept of outer measure was used in Theorem 52. This in turn was neded when establishing the integrability of functions that arise in quantum mechanics—in Theorem 219, for instance.
in RT × , and the extended concept of outer measure was used in Theorem 52. This in turn was neded when establishing the integrability of functions that arise in quantum mechanics—in Theorem 219, for instance.
The Henstock or Burkill-complete integral of a function h of elements x, N, I[N], structured by a multifunctional relation of association between these elements, is defined by forming Riemann sums of the function values, the terms of the Riemann sums being selected by “gauge” rules. Thus the concept of measurability is not required in order to define the integral. However, in proving certain properties of the integral, or in calculating the integrals of particular functions, the analytical device of step functions can be very useful.
Measurability is related to the construction of step functions and limits of step functions. These are formed from measurable sets, and they are measurable functions. Existence of the integrals of particular functions can often be established by means of step functions, using measurable sets and measurable functions. Accordingly, this section considers measurability of sets and functions in RT, in a further extension of the meaning given to the measurability concept in Section 4.9.
Suppose h(x, N, I[N]) is a real- or complex-valued function, and suppose A is a subset of 
defined as follows. Given
defined as follows. Given
; so
Some of the elements (x, N, I[M]) of A may be associated, with M = N and x(N) a vertex of I(N) in RN; where I(t) is a proper subset of R for t 
N, and I(t) = R for t T\N. But it is not assumed that the elements (x, N, I[M]) of A are associated.
N, and I(t) = R for t T\N. But it is not assumed that the elements (x, N, I[M]) of A are associated.
Measurability has been defined and used in preceding chapters. The following definition extends the meaning of measurability.
Definition 58 Suppose a function h of associated elements (x, N, I[N]) is given. Then a set A is h-measurable if the function
1A(x, N, I[N])h(x, N, I [N])
is integrable on RT.
Example 70 To relate this definition to preceding chapters, take h to be a distribution function F defined on I(RT), and let A be a subset of RT. Then A is F-measurable if and only if 1A(x)F(I) is integrable on RT. In terms of Definition 58, take
so 
and 
. Then 1A(x)F(I) is integrable if and only if 1A(x, N, I[N])F(I) is integrable. 
and . Then 1A(x)F(I) is integrable if and only if 1A(x, N, I[N])F(I) is integrable.
Example 71 If none of the triples (x, N, I[N]) A = 
is associated, then, for any function h(x, N, I[N]), the set A is not h-measurable. This is because the function 1A(x)h(x, N, I[N]) is not integrable if A has no associated elements.
A = is associated, then, for any function h(x, N, I[N]), the set A is not h-measurable. This is because the function 1A(x)h(x, N, I[N]) is not integrable if A has no associated elements.
If 1A(x, N, I[N])h(x, N, I[N]) is integrable on RT, write
If h(x, N, I[N]) is integrable on RT, write
The following results show that the function Ph has some probability-like properties.
Ph[\A] = Ph [RT] – Ph [A].
By assumption, both functions on the right are integrable on RT, so the result follows by Theorem 13.
Ph [A1 ⋃ A2] = Ph [A1] + Ph [A2].
Proof. Write A = A1 ⋃ A2. We have
By Theorem 13,
giving the result.
This result can be extended to any finite union of disjoint sets. For the union of a countably infinite family of disjoint sets the following holds.
Theorem 249 Suppose h and |h| are integrable on RT. If A1, A2, A3, . . . are disjoint h-measurable subsets of 
, then 
is h-measurable, and
, then is h-measurable, and
Proof. Let 
and let 
. Then, for each (x, N, I[N]), as k → ∞,
and let . Then, for each (x, N, I[N]), as k → ∞,
|1A (x, N, I[N])h(x, N, I[N]) – 1Bk (x, N, I[N])h(x, N, I[N])| < |h(x, N, I[N])|
with
1Bk (x, N, I[N])h(x, N, I[N]) < |h(x, N, I[N])|
for each k. By hypothesis, |h(x, N, I[N])| is integrable on RT, so Theorem 61 (dominated convergence theorem) gives the result.
as in
The countable additivity of the function Ph on disjoint measurable sets is, ultimately, a consequence of the properties of Riemann sums restricted by gauges.
Next, we define h-measurability of a function f.
Definition 59 Given a real-valued function f(x, N, I[N]) and a real- or complex-valued function h(x, N, I[N]), we say that f is h-measurable if, for each a R, each of the sets
R, each of the sets
is h-measurable.
If f is complex-valued then its h-measurability is defined by taking the real and imaginary parts of f.
Theorem 250 Suppose h(x, N, I[N]) is real-valued, with h and |h| integrable on RT. Suppose f(x, N, I[N]) is real-valued, with fh, f|h|, and |fh|, integrable on RT. Then f is h-measurable.
Proof. Suppose υ R. Write
R. Write
k(x, N, I[N]) = |f(x, N, I[N])h(x, N, I[N])| + |υh(x, N, I[N])|.
The function k is integrable because, by hypothesis, each of the two functions of which it is composed is integrable. We have
By Theorem 59, the functions
are integrable for every υ R. Thus, for u < υ, the function p(x, N, I[N]) given by
R. Thus, for u < υ, the function p(x, N, I[N]) given by
max {min {f(x, N, I[N])h(x, N, I[N]), υh(x, N, I[N])}, uh(x, N, I[N])} is integrable. We have
whenever, respectively,
Let S denote
{(x, N, I[N]) : f(x, N, I[N]) ≥ υ}.
Write
q(x, N, I[N]) = p(x, N, I[N]) – uh(x, N, I[N]).
For u → υ–,
monotonically; so for any ε > 0 there exists uε so that υ > u > uε implies
is h-measurable for every υ. The other conditions required for h-measurability of f are proved similarly.
Theorem 251 Suppose h and |h| are integrable on RT. Suppose A1, A2 are h-measurable subsets of 
. Then A1 ∩ A2 and A1 ⋃ A2 are h-measurable subsets of .
. Then A1 ∩ A2 and A1 ⋃ A2 are h-measurable subsets of .
Proof. The functions 1AJ (x, N, I[N])h(x, N, I[N]) are integrable on RT for j = 1, 2. Therefore, with
f(x, N, I[N]) = 1A1(x, N, I[N]) + 1A2(x, N, I[N]),
the function f(x, N, I[N])h(x, N, I[N]) is integrable. Since this function satisfies the conditions of Theorem 250, the set
{(x, N, I[N]) : f(x, N, I[N]) ≥ 2}
is h-measurable. In other words A1 ∩ A2 is h-measurable. Then
1A1 ⋃ A1 (x, N, I[N]) = f(x, N, I[N]) – 1A1 ∩ A2 (x, N, I[N]),
and this is h-integrable, giving the second result.
Theorem 252 Suppose h and |h| are integrable on RT. If, for each a R, any one of the sets in Definition 59 is h-measurable, then the other sets are also h-measurable; and the function f(x, N, I[N]) is h-measurable.
R, any one of the sets in Definition 59 is h-measurable, then the other sets are also h-measurable; and the function f(x, N, I[N]) is h-measurable.
Proof. Denote the sets in Definition 59 by S1, S2, S3, S4, respectively. Assume S1 is h-measurable. Then
Since h and 1S1 h are integrable on RT, so is 1S2 h. As this holds for each a R, S2 is h-measurable. Given a R, for j = 1, 2, 3, . . . , let
R, S2 is h-measurable. Given a R, for j = 1, 2, 3, . . . , let
and, for k = 1, 2, 3, . . . , let 
. As k → ∞,
. As k → ∞,
1Bk(x, N, I[N]) → 1S3(x, N, I[N]).
Let
gk(x, N, I[N]) = 1Bk(x, N, I[N])h(x, N, I[N]).
For each k, gk is integrable on RT, with
|gk(x, N, I[N])| ≤ |h(x, N, I[N])|
Given ε > 0, there exists kε so that k ≥ kε implies
As |h| is integrable, Theorem 61 (dominated convergence theorem) implies that
is integrable on RT, so S3 is h-measurable; and this holds for each a R. Finally,
R. Finally,
giving h-measurability of S4. Thus h-measurability of S1 implies h-measurability of Sj for j = 2, 3, 4. The other cases are proved similarly.
Theorem 253 Suppose h(x, N, I[N]) is real- or complex-valued, with h and |h| integrable on RT, and suppose f(x, N, I[N]) is h-measurable. Given α R, the functions f(x, N, I[N]) + α and αf(x, N, I[N]) are h-measurable.
R, the functions f(x, N, I[N]) + α and αf(x, N, I[N]) are h-measurable.
Proof. For each a R, f(x, N, I[N]) < a implies f(x, N, I[N]) + α < a + α and
R, f(x, N, I[N]) < a implies f(x, N, I[N]) + α < a + α and
αf(x, N, I[N]) ≤ αa or αf(x, N, I[N]) ≥ αa
depending on whether α and a are, respectively, positive, negative, or zero. The result then follows from Theorem 252.
Theorem 254 Suppose h(x, N, I[N]) is real- or complex-valued, with h and |h| integrable on RT, and suppose f(x, N, I[N]) is h-measurable. Then the functions (f(x, N, I[N]))2 and |f(x, N, I[N])| are h-measurable. Suppose, in addition, that f ≠ 0. Then 
is h-measurable.
is h-measurable.
Proof. If a < 0 then 
is the set
is the set
{(x, N, I[N]) : |f(x, N, I[N])| > a}, = R,
which is h-measurable since h is integrable. Suppose a ≥ 0. Denote the sets
by A1, A2, A3, A4, respectively. Each of these sets is h-measurable by hypothesis. Then
With f ≠ 0, the set
so is h-measurable.
is h-measurable.
Theorem 255 Suppose h and |h| are integrable on RT and suppose each of the functions f(x, N, I[N]) and g(x, N, I[N]) is h-measurable. Then the set
A = {(x, N, I[N]) : f(x, N, I[N]) < g(x, N, I[N])}
is h-measurable.
Proof. Let {a1, a2, a3, . . .} be an enumeration of the rational numbers in R, and let
Then 
and, as k → ∞,
and, as k → ∞,
with |1Bk (x, N, I[N])h(x, N, I[N])| ≤ |h(x, N, I[N])| for each k, where both 1Bk (x, N, I[N])h(x, N, I[N]) and |h| are integrable. Given ε > 0 we can find kε so that, for k ≥ kε,
|1A(x, N, I[N])h(x, N, I[N]) – 1Bk(x, N, I[N])h(x, N, I[N])| < ε|h(x, N, I[N])|.
Therefore by Theorem 61 (dominated convergence theorem),
1A(x, N, I[N])h(x, N, I[N])
is integrable on RT, giving the result.
Theorem 256 Suppose h and |h| are integrable on RT, and suppose the functions f and g are h-measurable. Then the functions f – g, f + g, and fg are h-measurable. If g ≠ 0 then 
is h-measurable.
is h-measurable.
Proof. For any a R, the function a + g(x, N, I[N]) is h-measurable. Then, by Theorem 255, the set
R, the function a + g(x, N, I[N]) is h-measurable. Then, by Theorem 255, the set
is h-measurable, so the function f – g is h-measurable. Since f + g = f – (–g), the function f + g is h-measurable. Since
the preceding results imply that the function fg is h-measurable. For the final part use Theorem 254.
Theorem 257 Suppose h and |h| are integrable on RT, and suppose the real-valued functions fj are h-measurable for j = 1, 2, 3, . . . . If, for each (x, N, I[N]),
fj(x, N, I[N]) → f(x, N, I[N])
as j → ∞, then f(x, N, I[N]) is h-measurable.
Proof. Given ε > 0 there exists jε so that, for j > jε,
|f(x, N, I[N]) – fj(x, N, I[N])| < ε.
For any a R write A = {(x, N, I[N]): f(x, N, I[N]) > a}. Then, for each (x, N, I[N]),
R write A = {(x, N, I[N]): f(x, N, I[N]) > a}. Then, for each (x, N, I[N]),
|1Afjh – 1Afh| < ε|h(x, N, I[N])|;
the function 1A(x, N, I[N]) fj (x, N, I[N])h(x, N, I[N]) is integrable for each j; and, for each (x, N, I[N]),
1Afjh → 1Afh
as j → ∞. Theorem 61 then implies that the function
1A (x, N, I[N]) f (x, N, I[N]) h (x, N, I[N])
is integrable. Since this holds for each a, f is h-measurable.
The following is a partial converse to Theorem 252.
Theorem 258 Suppose k(x, N, I[N]) ≥ 0 is h-integrable on RT. If real-valued f(x, N, I[N]) is h-measurable with f ≥ k then f is h-integrable on RT.
f(x, N, I[N]) + k(x, N, I[N]) ≤ 2k(x, N, I[N]),
so, if we replace f by f + k, we can take f(x, N, I[N]) ≥ 0. For each positive integer m write
and define
fm(x, N, I[N]) = aj for aj ≤ f(x, N, I[N]) < aj + 1, j = 0, 1, . . . , m2m,
with fm(x, N, I[N]) = 0 if f(x, N, I[N]) ≥ m2–m. Then
fm(x, N, I[N]) ≤ f(x, N, I[N]) < fm + 1(x, N, I[N])
and, for each (x, N, I[N]), fm (x, N, I[N]) is monotone increasing as m increases. By measurability of f, fm is h-integrable, and, writing
Aj = {(x, N, I[N]) : aj ≤ f(x, N, I[N]) < aj + 1},
we have
Then fm → f monotonically as m → ∞, and, given ε > 0, for each (x, N, I[N]) there exists mε so that m > mε implies
|f(x, N, I[N])h(x, N, I[N]) – fm(x, N, I[N])h(x, N, I[N])| < εk(x, N, I[N]).
Provided f ≥ 0 Theorem 57 then implies fh is integrable on RT. If f does not satisfy f ≥ 0 we can apply the preceding argument to the function f + k. Then the h-integrability of k and f + k implies that (f + k) – k is h-integrable.
Instead of assuming |f| ≤ k the preceding result remains valid if it is assumed that
- f ≥ 0 and
is bounded as m → ∞.
Theorem 259 Suppose f(x, N, I[N]) is real-valued and f is h-measurable, and suppose g is a continuous real-valued function defined on domain R. Then the composite function g o f is h-measurable.
Proof. If J I(R), then, by continuity of g, g–1(J) I(R). Write g–1(J) = K. By h-measurability of f, the set
I(R), then, by continuity of g, g–1(J) I(R). Write g–1(J) = K. By h-measurability of f, the set
is h-measurable. Thus, for each J I(R), the set
I(R), the set
is h-measurable.
A.2. Historical Note
The background to most of the themes in this book is amply documented in books and articles. The purpose of this note is to draw together some of the newer themes.
The neo-Riemannian approach to the theory of integration associated with R. Henstock, J. Kurzweil, and E.J. McShane originated in the 1950s. The concepts of generalized Riemann (or -complete) integration appeared in 1955 [85], where they provided a solution1 to a problem of integrals converging or diverging in the presence of convergence factors, a problem first posed by Henstock to his Ph.D. supervisor Paul Dienes in 1944. Using his new techniques of Riemann-, Stieltjes- and Burkill-complete integration, Henstock solved a number of problems in papers published in the 1950’s and 1960’s, culminating in 1963 [93], Theory of Integration, which gives a comprehensive account of an early version2 of the Burkill-complete integral.
Addressing the 1962 International Congress of Mathematicians in Stockholm, Henstock declared: “Lebesgue is dead!”. Nonetheless, Henri Lebesgue (1875–1941) set the standard [137, 138] for integration in the twentieth century. But Henstock’s flamboyant announcement that a new, deeper, and more effective theory of integration had arrived was already being given additional substance by J. Kurzweil, who independently discovered Riemann-complete integration, and whose 1957 paper [129], Generalised ordinary differential equations and continuous dependence on a parameter, was creating a new paradigm for the theory of differential equations.
However, this was not the end of the story, just the beginning. The Riemann-complete integral provides a resolution of the problem of integrating derivatives, a significant but relatively small area of the subject — one that had already been resolved [82] by Denjoy and Perron, using other methods.
Integrability of derivatives is part of the broader subject of integrability of functions, strands of which can be found [101, 102] in the work of authors such as C-J. de la Valleé Poussin, W.H. Young, J.C. Burkill3, L.C. Young, A.J. Ward, and S. Saks. Prom the point of view of this book a particularly important strand is the role of integration in the theory of approximation.
About 1962, a physics researcher introduced Henstock to Feynman integrals. These constructs had been used successfully in quantum electrodynamics and quantum mechanics since 1948. But no adequate mathematical method had been found to establish whether or not they “converge”; that is, whether or not they actually exist in a specifically mathematical sense.
Feynman “integrals” worked in practice but not in theory, so to speak. The Lebesgue integral method could not resolve this issue. Could the new developments in integration theory do it?
Independently, E.J. McShane [11] engaged with this problem. His search for a rigorous mathematical foundation for quantum field theory led him to sustained and productive involvement in the theory of integration. Seminal works by McShane include Order-Preserving Maps and Integration Processes in 1953 [153], Integrals devised for special purposes in 1963 [155], and Stochastic Calculus and Stochastic Models in 1974 [158] which has influenced much subsequent research in this area.
With such issues in mind, over the following decade Henstock embarked on construction of a more complete and comprehensive theory of integrability. Page 221 of Henstock [94] has a sketch of the historical background to some of these developments. His paper [101] (The Lebesgue syndrome) is also useful.
In 1968 Henstock published Linear Analysis [94], containing an initial exposition of the division system (or integration basis) theory described in Chapter 4 of the present book, and designated here as the Henstock integral. The Riemann-complete integral is a modification of the nineteenth century Riemann integral with more delicate selection of Riemann partitions. The Henstock integral is based on
- classes of objects to be selected for partitions, and
- classes of rules for determining selections.
This system ties the integral to partitions (or divisions), while setting it free from any particular method of integration such as Riemann’s or Lebesgue’s. The basic idea can be described simply as the Riemann sums method in integrability.
Linear Analysis [94] also included an outline (Exercises 43.12 and 43.13, pages 223−224) of integration in infinite-dimensional spaces, partly inspired by a 1934 paper [116] by B. Jessen.
The Henstock integral, or division system, was presented more fully [96] in 1969. Henstock’s 1973 paper [97] showed how the system can be used to define the Feynman integral as a manifestation of an extended version of Brownian motion, thus bringing this construct into the realm of mathematical integration and probability—which is the spirit in which it was originally presented by Feynman in 1948 [64].
Henstock [97] is a historic landmark in the theory of integration and probability. But definition is one thing, existence is another. It was not shown that the Feynman integrals exist. And if they do happen to exist, means of performing more than a finite number of operations on them were not provided. The question of taking limits under the integral sign for highly oscillating integrands such as these was broached in this paper, but not resolved there. Henstock concluded [97] with the statement: “This shows how difficult it will be to introduce Lebesgue-type limit theorems into Feynman integration theory. ... The present paper is only the beginning of the theory.”
A General Theory of Integration in Function Spaces (Muldowney [162], 1987) extended the division system idea, providing a broader structure of associated elements of integration with the capacity to deal with a large class of integrands in infinite-dimensional spaces. This is, essentially, the extension of the Henstock integral used in the present book. It was established in Muldowney [162] that the Feynman integrals of step function integrands exist. But to go further requires taking limits under the integral sign.
Criteria for taking limits under the integral sign—actually a departure from the “Lebesgue-type limit theorems” mentioned in the conclusion of Henstock [97]—were constructed by Henstock for the division system theory in his 1991 book The General Theory of Integration [105]. These criteria (similar to Theorems 62, 64, and 65 in this book) were used by Muldowney in 2000 [165] to establish properties of the Feynman integral.
The partitioning theorem (Theorem 4) used in this book was proved by Hen-stock, Muldowney, and Skvortsov [107]. This result simplifies and strengthens the definition of the Henstock integral in infinite-dimensional spaces.
The following assessment of the Feynman integral problem was given by Gian-Carlo Rota in his 1997 book Indiscrete Thoughts [199]:
The Feynman path integral is the mathematicians’ pons asinorum. Attempts to put it on a sound footing have generated more mathematics than any subject in physics since the hydrogen atom. To no avail. The mystery remains, and it will stay with us for a long time.
The Feynman integral, one of the most useful ideas in physics, stands as a challenge to mathematicians. While formally similar to Brownian motion, and while admitting some of the same manipulations as the ones that were made rigorous long ago for Brownian motion, it has withstood all attempts at rigor. Behind the Feynman integral there lurks an even more enticing (and even less rigorous) concept: that of an amplitude which is meant to be the quantum-mechanical analog of probability ...
Rota went on to say: A concept similar to that of a sample space should be brought into existence for amplitudes and quantum mechanics should be developed starting from this concept.
The brief historical outline in this section is a summary of developments in the agenda described by Rota—a work-in-progress but no longer a mystery.
__________
2Henstock’s Theory of Integration uses the term “Riemann-complete” to designate integrals which, in this book, are distinguished from each other as Riemann-complete, Stieltjes-complete and Burkill-complete.
3In [103] Henstock cites the work of W.H. Young, as developed in lectures by J.C. Burkill, as a source of the inspiration for his proof by Riemann sums of Fubini’ s theorem. Henstock attended a course of lectures by J.C. Burkill in 1942. His Ph.D. thesis, Interval Functions and Their Integrals (1948), is an investigation and extension of J.C. Burkill’s idea of the integral.
No hay comentarios:
Publicar un comentario