Update README.md
Signed-off-by: David Rotermund <54365609+davrot@users.noreply.github.com>
This commit is contained in:
David Rotermund
GitHub
parent
511d67cbae
commit
29a7ab89cb
1 changed files with 42 additions and 36 deletions
|
|
@ -25,11 +25,11 @@ $$\hat{x}(\omega) = {1 \over 2\pi} \int_{-\infty}^{+\infty} x(t) \exp\left( -i\
|
|||
|
||||
Attention: in the literature sometimes a pre-factor of $1/\sqrt{2\pi}$ is used. The factor in our definition is chosen such that it is compatible with the normalization of the Matlab-FFT which is described later on. The Fourier transformation is reversible; the correspondent reverse transformation is
|
||||
|
||||
$$x(t) = \int_{-\infty}^{+\infty} \hat{x}(\omega) \exp\left( i\omega t \right) \, d\omega \,$$.
|
||||
$$x(t) = \int_{-\infty}^{+\infty} \hat{x}(\omega) \exp\left( i\omega t \right) \, d\omega \,$$ .
|
||||
|
||||
If $x(t)$ is a real valued function, the transform $\hat{x}(\omega)$ is usually complex. For the existence of the Fourier transformation it is sufficient that $x(t)$ is absolutely integrable, i.e. $\int_{-\infty}^{+\infty} | x(t) | dt \, < \, \infty$.
|
||||
If $x(t)$ is a real valued function, the transform $\hat{x}(\omega)$ is usually complex. For the existence of the Fourier transformation it is sufficient that $x(t)$ is absolutely integrable, i.e. $\int_{-\infty}^{+\infty} | x(t) | dt \, < \, \infty$ .
|
||||
|
||||
The Fourier transformation of a function is a representation of that function as a superposition of harmonic oscillations with angular frequency $\omega=2\pi f$ and coefficients $\hat{x}(\omega)$. It is thus nicely suited for frequency analysis of temporal signals. The complex coefficients $\hat{x}$ contain information about phase $\phi$ and amplitude $A$ of the oscillation, and yield by using Euler's identity:
|
||||
The Fourier transformation of a function is a representation of that function as a superposition of harmonic oscillations with angular frequency $\omega=2\pi f$ and coefficients $\hat{x}(\omega)$ . It is thus nicely suited for frequency analysis of temporal signals. The complex coefficients $\hat{x}$ contain information about phase $\phi$ and amplitude $A$ of the oscillation, and yield by using Euler's identity:
|
||||
|
||||
$$\hat{x} = \Re(\hat{x}) + i \Im(\hat{x})$$
|
||||
|
||||
|
|
@ -41,11 +41,11 @@ Here, $i$ is the imaginary unit and $\Re$ and $\Im$ denote the real respectively
|
|||
|
||||
If $x(t)$ is periodic, for example in $2\pi$, or if $x(t)$ is only defined in the interval $[0, 2\pi]$, this can be expressed as the Fourier series with coefficients $\hat{x}_k$:
|
||||
|
||||
$$\hat{x}_k = {1 \over 2\pi} \int_{0}^{2\pi} x(t) \exp\left( -ikt \right) \, dt \,$$. (9.1)
|
||||
$$\hat{x}_k = {1 \over 2\pi} \int_{0}^{2\pi} x(t) \exp\left( -ikt \right) \, dt \,$$ . (9.1)
|
||||
|
||||
The reverse transformation is written as an infinite sum:
|
||||
|
||||
$$x(t) = \sum_{k=-\infty}^{+\infty} \hat{x}(t) \exp\left( ikt \right) \,$$.
|
||||
$$x(t) = \sum_{k=-\infty}^{+\infty} \hat{x}(t) \exp\left( ikt \right) \,$$ .
|
||||
|
||||
By rescaling of the t-axis by $t'=2\pi t/T$, a function $a(t)$ defined on an interval $[0, T]$ can always be mapped on a function $x(t')$ in the interval $[0, 2\pi]$, and thus can be transformed by equation (9.1).
|
||||
|
||||
|
|
@ -54,15 +54,15 @@ For further properties of the Fourier transformation you may consult the instruc
|
|||
### Discrete Fourier Transformation
|
||||
In the computer, functions are defined on discrete sampling points in a finite interval. Let us assume a function $a(t)$ is sampled at $N$ equidistant points $t_n=Tn/N$, and the values at these points are $a_n=a(t_n)$. By the transformation $t'=(2\pi/T) t$ we transfer $a(t)$ to a function $x(t')$ in the interval $[0, 2\pi]$ and approximate the integral from equation (9.1) with the midpoint rule (see chapter Integration an Differentiation):
|
||||
|
||||
$$\hat{x}_k = {1 \over 2\pi} \int_{0}^{2\pi} x(t') \exp\left( -ikt' \right) \, dt' \approx {1 \over 2\pi} \sum_{n=0}^{N-1} a_n \exp\left( -ik 2\pi t_n/T \right) \Delta t' \,$$.
|
||||
$$\hat{x}_k = {1 \over 2\pi} \int_{0}^{2\pi} x(t') \exp\left( -ikt' \right) \, dt' \approx {1 \over 2\pi} \sum_{n=0}^{N-1} a_n \exp\left( -ik 2\pi t_n/T \right) \Delta t' \,$$ .
|
||||
|
||||
Here, $\Delta t'$ is given by $(2\pi/T)(T/N)=2\pi/N$, which means that
|
||||
|
||||
$$\hat{x}_k \approx A_k = {1 \over N} \sum_{n=0}^{N-1} a_n \exp\left( -i 2\pi nk/N \right) \,$$. (9.2)
|
||||
$$\hat{x}_k \approx A_k = {1 \over N} \sum_{n=0}^{N-1} a_n \exp\left( -i 2\pi nk/N \right) \,$$ . (9.2)
|
||||
|
||||
This equation describes the discrete Fourier transformation, the implementation of which we will discuss more extensively in the following paragraph. The corresponding reverse transformation is:
|
||||
|
||||
$$a_n = \sum_{k=0}^{N-1} A_k \exp\left( i 2\pi nk/N \right) \,$$.
|
||||
$$a_n = \sum_{k=0}^{N-1} A_k \exp\left( i 2\pi nk/N \right) \,$$ .
|
||||
|
||||
### Fast Fourier Transform (FFT)
|
||||
The discrete Fourier transformation is invoked in Matlab by the command a_trans = fft(a);, where $a$ is a function or time series. The fast Fourier transform is an algorithm that is optimized for a fast execution of the transformation. It is based on the idea to first divide a function of $N$ sampling points into two function of $N/2$ sampling points, and then to compute two Fourier transformations on the partial functions and combining the two results to the sought transformation of the full function. This procedure can be iteratively extended: further partial function then have $N/4$, $N/8$, $\ldots$ sampling points. The reverse transformation, the inverse Fourier transformation, is invoked by the command ifft.
|
||||
|
|
@ -124,7 +124,7 @@ a_sine = 42*sin((0:(n-1))*dt*2*pi*f);
|
|||
a_step = (0:(n-1)) >= n/2;
|
||||
```
|
||||
|
||||
The signals are sampled with a frequency of $1/dt$ for $T=1 s$, i.e. the discretization of the functions has the step width $dt$. The FFT is calculated by the following:
|
||||
The signals are sampled with a frequency of $1/dt$ for $T=1 s$ , i.e. the discretization of the functions has the step width $dt$ . The FFT is calculated by the following:
|
||||
|
||||
```matlab
|
||||
fft_a_sine = 1/n*fft(a_sine);
|
||||
|
|
@ -134,11 +134,11 @@ fft_a_step = 1/n*fft(a_step);
|
|||
|
||||
Now we want to investigate how much energy is contained in the different frequency bands. For this we calculate the so called power spectrum $\Phi(f)$, i.e the square of the sum of oscillation amplitudes for a specific frequency $f$. Here one has to be aware of two possible stumbling blocks:
|
||||
|
||||
* In the first place we have to make sure, how Matlab stores the calculated coefficients $A_k$ in the result vector: The first entry, e.g. a_sine(1) contains $A_0$, i.e. the mean of the function $a(t)$. The $k$-th entry then holds the coefficient $A_{k-1}$, for example a_sine(17) would be $A_{16}$.
|
||||
* Which frequency does the coefficient $A_k$ correspond to? $A_k$ corresponds to the harmonic oscillation, for which exactly $k$ periods fit into the time interval of length $T$. The frequency is thus $f_k = k/T$.
|
||||
* Let us now look at the coefficients $A_k$ and $A_{N-k}$. From equation (9.2) we can see that $\exp(-i2\pi (N-k)n/N)) = \exp(i2\pi kn/N)\exp(-i2\pi n) = \exp(i2\pi kn/N)$. Hence, $A_{N-k}$ is equivalent to the coefficient $A_{-k}$ with negative wave number $-k$ and thus corresponds to the same frequency than $A_k$.
|
||||
* In the first place we have to make sure, how Matlab stores the calculated coefficients $A_k$ in the result vector: The first entry, e.g. a_sine(1) contains $A_0$ , i.e. the mean of the function $a(t)$ . The $k$-th entry then holds the coefficient $A_{k-1}$ , for example a_sine(17) would be $A_{16}$ .
|
||||
* Which frequency does the coefficient $A_k$ correspond to? $A_k$ corresponds to the harmonic oscillation, for which exactly $k$ periods fit into the time interval of length $T$. The frequency is thus $f_k = k/T$ .
|
||||
* Let us now look at the coefficients $A_k$ and $A_{N-k}$ . From equation (9.2) we can see that $\exp(-i2\pi (N-k)n/N)) = \exp(i2\pi kn/N)\exp(-i2\pi n) = \exp(i2\pi kn/N)$ . Hence, $A_{N-k}$ is equivalent to the coefficient $A_{-k}$ with negative wave number $-k$ and thus corresponds to the same frequency than $A_k$ .
|
||||
|
||||
If we want to know the amplitude of the oscillation with frequency $f=5 Hz$, we have to look at the coefficient with the wave number $k = T f_k = 5$ and $k = -5$. The correspondent value $\Phi(f)$ is found in Matlab-notation as:
|
||||
If we want to know the amplitude of the oscillation with frequency $f=5 Hz$ , we have to look at the coefficient with the wave number $k = T f_k = 5$ and $k = -5$ . The correspondent value $\Phi(f)$ is found in Matlab-notation as:
|
||||
|
||||
```matlab
|
||||
phi_5Hz = (abs(fft_a_sine(1+5))+abs(fft_a_sine(n-5)))
|
||||
|
|
@ -156,35 +156,36 @@ ylabel('power');
|
|||
```
|
||||
|
||||
|
||||
|
||||
Figure 9.1.: Power spectrum of a step function.
|
||||
|
||||
## The Convolution Theorem
|
||||
|
||||
The fast Fourier transform is not only useful to determine the frequency content of temporal signals, but also to perform so called 'convolutions'. A convolution of two functions $f(t)$ and $g(t)$ is defined by the following equation:
|
||||
|
||||
$h(t) = \int_{-\infty}^{+\infty} f(t') g(t-t') \, dt'$ (9.3)
|
||||
$$h(t) = \int_{-\infty}^{+\infty} f(t') g(t-t') \, dt'$$ (9.3)
|
||||
|
||||
Here, $g$ is often referred to as convolution kernel. The procedure of a convolution can be visualized as follows: We interpret $f$ as a arbitrary function, that is compared to a template function $g$. To account for possible shifts on the time axis, the template is shifted between $t=+\infty$ and $t=-\infty$ and the comparison is made at every position $t$. If both functions match perfectly, the result $h(t)$ is big. If the shifted $g$ does not match with $f$, the result is small.
|
||||
|
||||
Applications of this equation are numerous in physics, and reach from simple filter operations to timely resolved frequency analysis (examples will be discussed later). First, we want to understand the connection between convolutions and the Fourier transformation. As a short notation of the application of the Fourier transform $F$ to a function $f$ (resp. its reverse transformation $F^{-1}$) we introduce:
|
||||
|
||||
$\hat{f}(k) = F\left.\left.\left[ f(t) \right]\right( k \right)$
|
||||
$$\hat{f}(k) = F\left.\left.\left[ f(t) \right]\right( k \right)$$
|
||||
|
||||
$f(t) = F^{-1}\left.\left.\left[ \hat{f}(k) \right]\right( t \right)$
|
||||
$$f(t) = F^{-1}\left.\left.\left[ \hat{f}(k) \right]\right( t \right)$$
|
||||
|
||||
Now we apply the Fourier transform to both the left and the right side of the Definition (9.3) and gain after short computation,
|
||||
|
||||
$\hat{h}(k) = 2\pi \hat{f}(k) \hat{g}(k) \, ,$
|
||||
$$\hat{h}(k) = 2\pi \hat{f}(k) \hat{g}(k) \, ,$$
|
||||
|
||||
or, in short notation, $\hat{h} = 2\pi F[f] F[g]$. To get the sought result $h(t)$, we apply the inverse Fourier transform to the equation, which results in
|
||||
|
||||
$h = F^{-1}\left[ 2\pi F[f] \cdot F[g] \right] \, .$ (9.4)
|
||||
$$h = F^{-1}\left[ 2\pi F[f] \cdot F[g] \right] \, .$$ (9.4)
|
||||
|
||||
The advantage of equation (9.4) over (9.3) lies in the computation speed of FFT: despite a three-fold transformation, calculating equation (9.4) is faster than computing the integral in (9.3), since the convolution integral corresponds to an element-wise multiplication of the Fourier coefficients in the Fourier space $k$ -- try to verify!
|
||||
|
||||
### Non-periodic Functions
|
||||
|
||||
The discrete Fourier transformation can be applied to any arbitrary function $f$, periodicity is not required. The convolution theorem, on the other hand, does strictly speaking only apply for the Fourier integral. It thus demands either a definition of the function in $[-\infty, +\infty]$, or a periodic function on a finite interval. Only in the latter case, fft / ifft can be directly used, i.e. without 'post-treatment', to evaluate a convolution integral.
|
||||
The discrete Fourier transformation can be applied to any arbitrary function $f$, periodicity is not required. The convolution theorem, on the other hand, does strictly speaking only apply for the Fourier integral. It thus demands either a definition of the function in $[-\infty, +\infty]$ , or a periodic function on a finite interval. Only in the latter case, fft / ifft can be directly used, i.e. without 'post-treatment', to evaluate a convolution integral.
|
||||
|
||||
Still, the FFT can be used with non-periodic functions, provided the convolution kernel $g$ is sufficiently narrow, i.e. its absolute value almost vanishes within $m << n$ sampling points around $t=0$. This is for example given for exponentially decaying functions or Gaussian distributions. Unwanted boundary effects can be suppressed by simply discarding the first and last $m$ elements of the result vector. Another possibility is the continuation of a function above its borders by concatenating a vector of $k$ zeros. The command is fft(f, n+k), where $n$ is the length of the vector $f$. The extension of $f$ by zeros is only sensible in specific situations (see also remarks further down), for example if signal was 0 before and after a period of data acquisition, or when a filter is successively 'inserted' into a beginning signal.
|
||||
|
||||
|
|
@ -246,6 +247,7 @@ To end this section, some further examples for the application of the convolutio
|
|||
Fourier transformations are ill suited for signals, whose spectrum changes over time. This also includes signals which do contain periodic oscillations, but where the frequency is subject to fluctuations. In wavelet transformations, a signals is convoluted with a filter $g$, that only contains few oscillations of a specific period and decays to 0 outside of this range. An example is depicted in the following picture:
|
||||
|
||||
|
||||
|
||||
Figure 9.4.: Wavelet filter $g(t)$ for the frequency analysis of non-stationary signals
|
||||
|
||||
|
||||
|
|
@ -268,7 +270,7 @@ $$C(\tau) = F^{-1}\left[\left.\left. 2\pi F[ f(t) ](k) F[ g(t) ](-k) \right]\rig
|
|||
The following examples come from the theoretical neurophysics:
|
||||
|
||||
#### Reverse Correlation.
|
||||
A simple model for the response properties of neurons in the visual system assumes that the response $r(t)$ of a neuron is constructed by a linear superposition of a stimulus $s(x, y, t)$ at the location $(x, y)$ on the retina, multiplied with a weight function $w(x, y, \tau)$. $\tau$ is the delay between stimulus and neuronal response, and $g[\,\,\,]$ an additional point-non-linearity:
|
||||
A simple model for the response properties of neurons in the visual system assumes that the response $r(t)$ of a neuron is constructed by a linear superposition of a stimulus $s(x, y, t)$ at the location $(x, y)$ on the retina, multiplied with a weight function $w(x, y, \tau)$ . $\tau$ is the delay between stimulus and neuronal response, and $g[\,\,\,]$ an additional point-non-linearity:
|
||||
|
||||
$$r(t) \propto g\left[ \int_x dx \int_y dy \int_\tau d\tau w(x, y, \tau) s(x, y, t-\tau) \right]$$
|
||||
|
||||
|
|
@ -276,7 +278,7 @@ Again, $w$ can be interpreted as a filter. With this insight, the inner integral
|
|||
|
||||
When the stimulus $s(x, y, t)$ and the response $r(t)$ is known, the unknown filter $w$ can be determined under special conditions. This becomes specifically easy, when the stimulus is uncorrelated white noise:
|
||||
|
||||
$w(x, y, \tau) \propto \int r(t) s(x, y, t-\tau) dt$
|
||||
$$w(x, y, \tau) \propto \int r(t) s(x, y, t-\tau) dt$$
|
||||
|
||||
After application of the convolution theorem, we get:
|
||||
|
||||
|
|
@ -287,15 +289,15 @@ Neuronal networks usually have a one- or two-dimensional topology, where a neuro
|
|||
|
||||
$$\tau\dot{A}(x, t) = -A(x, t)+g[ I(x, t) ]$$
|
||||
|
||||
$I(x, t)$ denotes the incoming current of the neuron at position $x$; this is a sum across the weighted activities of all other neurons:
|
||||
$I(x, t)$ denotes the incoming current of the neuron at position $x$ ; this is a sum across the weighted activities of all other neurons:
|
||||
|
||||
$I(x, t) = \int w(x, x') A(x', t) dx'$
|
||||
$$I(x, t) = \int w(x, x') A(x', t) dx'$$
|
||||
|
||||
If the weight are invariant to translation via $w(x, x') = w(x-x')$, the solution of this integral is again a case for the well known convolution theorem.
|
||||
If the weight are invariant to translation via $w(x, x') = w(x-x')$ , the solution of this integral is again a case for the well known convolution theorem.
|
||||
|
||||
## Probabilities and Distribution Functions
|
||||
### Distributions
|
||||
Often the problem arises to estimate an underlying distribution $\rho(x)$ from a sample $\{x_i\}$ of cardinality $N$ ($i=1,\ldots,N$). The simplest possibility is to create a histogram of the values $x_i$. For this, a range of the $x$-axis between $[x_{min}, x_{max}]$ is divided into $M$ intervals $I_j$ of size $\Delta x = (x_{max}-x_{min})/M$. An approximation of $\rho$ is now given by finding the number of elements $x_i$ that fall within a specific interval $I_j=[x_{min}+(j-1)\Delta x, x_{min}+j\Delta x[$:
|
||||
Often the problem arises to estimate an underlying distribution $\rho(x)$ from a sample $\{x_i\}$ of cardinality $N$ ($i=1,\ldots,N$). The simplest possibility is to create a histogram of the values $x_i$. For this, a range of the $x$-axis between $[x_{min}, x_{max}]$ is divided into $M$ intervals $I_j$ of size $\Delta x = (x_{max}-x_{min})/M$. An approximation of $\rho$ is now given by finding the number of elements $x_i$ that fall within a specific interval $I_j=[x_{min}+(j-1)\Delta x, x_{min}+j\Delta x[$ :
|
||||
|
||||
$$h_j = \int_{x_{min}+(j-1)\Delta x}^{x_{min}+j\Delta x} \sum_i^N \delta(x-x_i) dx$$
|
||||
|
||||
|
|
@ -372,6 +374,7 @@ xlabel('x'); ylabel('\rho_{normal}(x)');
|
|||
```
|
||||
|
||||
|
||||
|
||||
Figure 9.5.: Uniform and normal distribution.
|
||||
|
||||
|
||||
|
|
@ -382,12 +385,13 @@ To generate random numbers from arbitrary distributions $\rho(x)$, we will schem
|
|||
Let $\rho_{max}$ be the maximum of the distribution $\rho$. To generate random numbers from $\rho(x)$, we have to make sure that the value $x$ occurs with the relative frequency $\rho(x)/\rho_{max}$. This can be assured by a simple procedure that uses the random number generator for uniformly distributed random numbers: we first draw a candidate $x_1$ for a random number from the area of definition of $\rho(x)$. Then we draw a second number $x_2$ from $[0, 1[$. If this one is smaller than $\rho(x_1)/\rho_{max}$, the first number is accepted, otherwise the procedure is repeated. See also the following picture: only pairs of numbers $(x_1, x_2)$ are accepted, that fall into the area under the renormalized distribution $\rho$.
|
||||
|
||||
|
||||
|
||||
Figure 9.6.: Hit-and-miss method: draw two uniformly distributed random numbers $x_1$ and $x_2$, and accept $x_1$ if both numbers fall into the blue area. By this, the probability distribution marked by the red line.
|
||||
|
||||
|
||||
The hit-and-miss method can be very time consuming. It is easier if the inverse function of the primitive $F$ of $\rho$ can be calculated. The idea is to start out with a uniformly distributed random number $y$ and find a transformation $x = g(y)$ such that the probabilities $\rho_u(y)\,dy$ and $\rho(x)\,dx$ are equal. With $\rho_u(y)=1$ one can derive
|
||||
|
||||
$$y = F(x) = \int_{-\infty}^x \rho(x') \, dx' \,$$,
|
||||
$$y = F(x) = \int_{-\infty}^x \rho(x') \, dx' \,$$ ,
|
||||
|
||||
and thus $x = g(y) = F^{-1}(y)$. Is $F^{-1}$ known, one simply draws the uniformly distributed random numbers $y$ and calculates the sought random numbers $x$ from the inverse function. If $F^{-1}$ can not be found analytically, there are two possibilities:
|
||||
|
||||
|
|
@ -397,7 +401,7 @@ The primitive $F$ can be numerically approximated as the integral of $\rho$ with
|
|||
##### Tabulation.
|
||||
If many random number have to be drawn, the inverse function $F^{-1}(y)$ can first be tabulated. For this one is advised to use equidistant sampling points $y_j$. In order to then draw a random number $x$ from $\rho(x)$, again first a random number $y_1$ is taken from a uniform distribution, and then the neighboring values $y_j$ and $y_{j+1}$ are found, for which $y_j \leq y_1 < y_{j+1}$. With another random number $y_2$ from a uniform distribution in the interval $[0, 1[$ the random number $x$ is:
|
||||
|
||||
$$x = y_2 (F^{-1}(y_{j+1})-F^{-1}(y_j)) + F^{-1}(y_j) \,$$.
|
||||
$$x = y_2 (F^{-1}(y_{j+1})-F^{-1}(y_j)) + F^{-1}(y_j) \,$$ .
|
||||
|
||||
## Regression Analysis
|
||||
Analysis of measured data is an important part in evaluating physical theories and models. Usually, some parameters of the physical models have to be fitted to the data and then the adapted model shall be scrutinized whether it describes the data well. The mandatory techniques of parameter fitting as well as the evaluation of the goodness of a fit will be discussed in the following.
|
||||
|
|
@ -410,11 +414,12 @@ $$Y(x;a_1,a_2) = a_1+a_2 x$$
|
|||
to a data set.
|
||||
|
||||
|
||||
|
||||
Figure 9.7.: Fitting a curve to data points.
|
||||
|
||||
It is the aim to adjust the values of these parameters such that the function $Y(x;a_1,\ldots,a_M)$ represents the date as well as possible. Intuitively we would assume that in case a fit is good, the graph of the function $Y(x;a_1,\ldots,a_M)$ will lie close to the points $(x_i,y_i)$ (see figure 9.7). We can quantify this statement by measuring the deviations between the data points and the function
|
||||
|
||||
$$\Delta_i = Y(x_i;a_1,\ldots,a_M)-y_i \,$$.
|
||||
$$\Delta_i = Y(x_i;a_1,\ldots,a_M)-y_i \,$$ .
|
||||
|
||||
We choose our fitting criterion such that the sum of the squares of the deviations becomes minimal. This means that we have to tune our parameter values $a_1,\ldots,a_M$ such that the "cost function"
|
||||
|
||||
|
|
@ -424,12 +429,12 @@ is minimized. This method is dubbed least squares method. It is not the only pos
|
|||
|
||||
Usually the data points have an estimated error bar (the confidence interval), that we denote as $y_i\pm \sigma_i$. In this case we should change our fitting criterion such that points with bigger error bars have a smaller weight. Following this logic we define
|
||||
|
||||
$$ \chi^2(a_1,\ldots,a_M) = \sum_{i=1}^N\left(\frac{\Delta_i}{\sigma_i}\right)^2 = \sum_{i=1}^N\frac{[Y(x_i;a_1,\ldots,a_M)-y_i]^2}{\sigma_i^2} \,$$.
|
||||
$$ \chi^2(a_1,\ldots,a_M) = \sum_{i=1}^N\left(\frac{\Delta_i}{\sigma_i}\right)^2 = \sum_{i=1}^N\frac{[Y(x_i;a_1,\ldots,a_M)-y_i]^2}{\sigma_i^2} \,$$ .
|
||||
|
||||
### Linear Regression
|
||||
We now want to fit a straight line to the data points,
|
||||
|
||||
$$Y(x;a_1,a_2) = a_1+a_2 x \,$$.
|
||||
|
||||
$$Y(x;a_1,a_2) = a_1+a_2 x \,$$ .
|
||||
|
||||
This type of fitting is called linear regression. The two parameters $a_1$ and $a_2$ have to be chosen such that
|
||||
|
||||
|
|
@ -449,7 +454,7 @@ $$\Sigma_{x^2} = \sum_{i=1}^N\frac{x_i^2}{\sigma_i^2}, \quad \Sigma_{xy} = \sum_
|
|||
|
||||
These sums are calculated directly from the data points, are thus known constants. Hence we can rewrite the previous system of equations as
|
||||
|
||||
$$a_1\Sigma + a_2\Sigma_x-\Sigma_y = 0$$
|
||||
$$a_1\Sigma + a_2\Sigma_x-\Sigma_y = 0$$
|
||||
|
||||
$$a_1\Sigma_x+a_2\Sigma_{x^2}-\Sigma_{xy} = 0$$
|
||||
|
||||
|
|
@ -463,7 +468,7 @@ $$\sigma_{a_j}^2 = \sum_{i=1}^N\left(\frac{\partial a_j}{\partial y_i}\right)^2\
|
|||
|
||||
Insertion of the equations (9.6) yields
|
||||
|
||||
$$\sigma_{a_1} = \sqrt{\frac{\Sigma_{x^2}}{\Sigma\Sigma_{x^2}-(\Sigma_x)^2}} , \quad \sigma_{a_2} = \sqrt{\frac{\Sigma}{\Sigma\Sigma_{x^2}-(\Sigma_x)^2}} \,$$.
|
||||
$$\sigma_{a_1} = \sqrt{\frac{\Sigma_{x^2}}{\Sigma\Sigma_{x^2}-(\Sigma_x)^2}} , \quad \sigma_{a_2} = \sqrt{\frac{\Sigma}{\Sigma\Sigma_{x^2}-(\Sigma_x)^2}} \,$$ .
|
||||
|
||||
As an example we consider:
|
||||
|
||||
|
|
@ -472,6 +477,7 @@ $$\hat{x}_k = {1 \over 2\pi} \int_{0}^{2\pi} x(t') \exp\left( -ikt' \right) \, d
|
|||
Here the fit parameters are $a_1 = 0.1529 \pm 0.2633$ and $a_2 = 1.0939\pm 0.0670$. Note that the error bars $\sigma_{a_j}$ do not depend on the $y_i$. These error bars are thus no quantifier of the goodness of the fit.
|
||||
|
||||
|
||||
|
||||
Figure 9.8.: An example of linear regression. A straight line $Y(x)$ is fitted to the data points $(x_i,y_i)$.
|
||||
|
||||
### Goodness of the Fit
|
||||
|
|
@ -481,11 +487,11 @@ $$|y_i-Y(x_i)| \approx \sigma_i \, .$$
|
|||
|
||||
Insertion into the definition of $\chi^2$ in
|
||||
|
||||
$$\chi^2(a_1,\ldots,a_M) = \sum_{i=1}^N\left(\frac{\Delta_i}{\sigma_i}\right)^2 = \sum_{i=1}^N\frac{[Y(x_i;a_1,\ldots,a_M)-y_i]^2}{\sigma_i^2} \,$$.
|
||||
$$\chi^2(a_1,\ldots,a_M) = \sum_{i=1}^N\left(\frac{\Delta_i}{\sigma_i}\right)^2 = \sum_{i=1}^N\frac{[Y(x_i;a_1,\ldots,a_M)-y_i]^2}{\sigma_i^2} \,$$ .
|
||||
|
||||
yields $\chi^2 \approx N$. The more parameters are used, the better the fit will be. For the case of $N=M$, the fit will be exact. For the case of the straight line, this simply means that to $N=2$ points, we can always exactly fit a straight line ($M=2$). Thus, a good criterion for the goodness of the fit is
|
||||
|
||||
$$\chi^2 \approx N-M \,$$.
|
||||
$$\chi^2 \approx N-M \,$$ .
|
||||
|
||||
As an example we again refer to fig. 9.8. Here, $\chi^2 \approx 4.5$ and $N-M = 8$. The goodness of the fit is thus quite high.
|
||||
|
||||
|
|
@ -504,7 +510,7 @@ $$Y(x;a_1,a_2) = a_1+a_2 x \,$$.
|
|||
|
||||
A second example are power laws of the type
|
||||
|
||||
$$Z(t;\alpha,\beta) = \alpha t^\beta \,$$.
|
||||
$$Z(t;\alpha,\beta) = \alpha t^\beta \,$$ .
|
||||
|
||||
the variable transformation
|
||||
|
||||
|
|
|
|||
Loading…
Reference in a new issue