|
|
|
Middle of the page Bottom of the page Selected Topics on Electrogastrography 4. METHODS TO PROCESS GASTRIC ELECTRICAL ACTIVITY 4.1. INFORMATION EXTRACTED FROM GASTRIC ELECTRICAL ACTIVITY. 4.1.1. MIGRATING MYOELECTRICAL COMPLEX AND CONTRACTIONS. It was pointed out earlier (see Section 1) that spikes indicate the presence of gastric contractions. In experiments in vivo, spikes are best seen in serosal SDB recordings, and this is the most reliable method of assessing contractions from gastric electrical signals (11, 12, 42). As a direct consequence of that a histogram of the Migrating Myoelectrical Complex (MMC) can be built from SDB recordings. Each value in this histogram represents the ratio between the number of ECA waves followed by spikes, and the total number of ECA waves in a certain time interval (usually 2 minutes). Whether contractions can be reliably assessed from LDB and EGG recordings is subject of a separate section. 4.1.2. ELECTRICAL COUPLING BETWEEN DIFFERENT PARTS OF THE STOMACH. Electrical activity in the normal stomach has a firm temporal organization. While the velocity of propagation from the distal corpus to the terminal antrum can differ individually, it has been shown that corresponding ECAs recorded from different areas of a normal stomach maintain a consistent time shift between them. This phenomenon has been called electrical coupling (12, 43). Even in cases of significant frequency irregularities, the pattern of electrical coupling in the stomach may be preserved (43). Normal electrical coupling between different parts of the stomach corresponds to a distal direction of propagation of gastric electrical and contractile activity (41, 43, 57). Electrical coupling can easily be studied invasively with a set of serosal SDB electrode pairs implanted in the areas of interest (42, 43). Although LDB signals cover greater stomach areas, they can also be used to assess electrical coupling (42). A recent study (41) reported that the direction of propagation could be assessed from EGG recordings as well. In another study, analysis of the EGG waveform has been suggested as an indirect method for non-invasive assessment of direction of propagation of GEA (39). However, a comparative frequency study of SDB and EGG signals to objectively prove these assumptions has not been done. 4.1.3. GASTRIC ELECTRICAL FREQUENCY. The frequency of GEA can be assessed by all recording techniques, but with different degrees of reliability. It could be very informative in identifying gastric electrical irregularities, which are believed to be related to irregularities in contractile activity (17, 39). Some investigators assume that tachygastria (gastric electrical signals with frequencies above 4 cpm) and bradygastria (below 2.25 cpm) could be related to certain gastric motility disorders (11, 17). Biphasic serosal SDB signals have very well defined time characteristics and are very reliable in assessing frequency (39, 42). It has been suggested that EGG signals can also show gastric electrical irregularities (9, 11, 39, 42, 43, 45). However, the relationship between EGG frequency and gastric electrical frequency obtained from SDB recordings has not been sufficiently explored. 4.2. METHODS TO ASSESS THE MIGRATING MYOELECTRICAL COMPLEX (MMC) AND CONTRACTIONS. 4.2.1. BUILDING THE MMC FROM SEROSAL SDB RECORDINGS. The biphasic shape of Electrical Control Activity (ECA) in SDB recordings of GEA is very well defined both in terms of amplitude and duration when appropriate amplification bandwidth is used (Fig.4.1). Each ECA wave is surrounded with quiescent intervals in which no spikes or remnants of the plateau of ERA are recorded. Spikes, on the other hand, usually have less power and are packed between consecutive ECA waves. This makes both ECA and spikes easily identifiable and therefore automatic computation of the MMC from serosal SDB recordings is relatively easy. Computer implementation of an algorithm to do that is relatively easy. A sample canine MMC obtained with such algorithm is shown in Fig.4.2.
Figure 4.1. SDB signal divided into three zones: A - the zone of ECA; B - the quiescent zones; C - the zone of spikes.
Figure 4.2. A typical canine MMC obtained digitally. All four phases of the MMC are clearly seen. 4.2.2. QUANTITATIVE ASSESSMENT OF EGG AMPLITUDES. Many authors (9, 11, 12, 14, 45) believe that the appearance of spikes in internal recordings of GEA corresponds to an increased EGG amplitude, i.e. gastric contractions can be assessed by amplitude analysis of EGG. Changes in EGG amplitudes are relatively easy to follow in the time domain, if the electrodes are positioned properly and the amplifier has an adequate bandwidth. Usually the average of the absolute values of all amplitudes is calculated. EGG amplitude is considered increased if it exceeds 25% of this average (42). The dynamics of EGG amplitudes can also be assessed in the frequency domain (10). This will be subject of special attention when frequency analysis is discussed. 4.2.3. OTHER METHODS TO ASSESS CONTRACTIONS AND MMC. The most reliable method to measure contractions is with implanted serosal force transducers (FT; 58, 59). Such recordings show MMC directly without any need for processing. The method is invasive and is usually used on unconscious dogs. In humans, contractions can be assessed indirectly by measuring changes in gastric intraluminal pressure (ILP, 60). This is usually done by means of special manometric tubes introduced orally and perfused with water. Once in the stomach, the intraluminal tube assembly is connected externally to a pneumo-hydraulic system (Arndorfer pump). Thus the pressure transducer records the changes in the stomach pressure only. Although less reliable than FT, ILP measurements are widely used, because they are less harmful for the patient. The Migrating Myoelectric Complex can be recognized in the time domain ILP recordings. Similarly to the changes in EGG amplitudes, the amplitude changes recorded using both ILP and FT methods can be assessed quantitatively by calculating the average of the absolute values of all amplitudes and taking 25% of this average as a threshold level of recognition (42). 4.3. METHODS TO ASSESS ELECTRICAL COUPLING. 4.3.1. VISUAL ANALYSIS OF THE TIME SHIFTS BETWEEN DIFFERENT CHANNELS. This is the simplest way of assessing electrical coupling and is illustrated on Fig.4.3. As long as the lines that connect the corresponding waves remain parallel, the coupling is intact.
Figure 4.3. Electrically coupled SDB signals from proximal (A) and distal (B) antrum. The lines connecting the peaks of ECA remained parallel even when a minor frequency irregularity occurred between the third and the fourth wave. 4.3.2. ANALYSIS OF THE TIME SHIFTS USING CROSSCORRELATION. Correlation between two groups of data means that they change with respect to each other in a structured way. If, for example, signal a always increases and decreases, as signal b increases and decreases for a certain amount of time, then x and y are positively correlated in this time interval (61). It is obvious that the crosscorrelation function of two shifted signals with similar periods is also periodic. It is also well known (see Challis and Kitney, 61) that the phase lag (time shift) between the two signals is coded as the time shift from k=0 to the first local maximum of the crosscorrelation function. Circular crosscorrelation considers the signal y to be symmetrical and continuous outside the index range (0, N-1). The crosscorrelation sequence in that case has length N, while for linear crosscorelation the length of the sequence is 2N-1 and zero-padding is necessary to avoid errors (61). Both methods produce different results, which are illustrated in the following two examples based on simple numerical sequences (65): (a) linear correlation: {1 1 0 0} {1 3 3 1} = {1 4 6 4 1 0 0} (b) circular correlation: {1 1 0 0} {1 3 3 1} = {2 4 6 4}. It is crucial to use linear correlation when processing non-periodic signals, because if circular crosscorrelation is used the signals would artificially be made periodic (61, 65). GEA signals, however, can be considered quasi periodic, so circular crosscorrelation would be acceptable. In this study it is preferred because its sequence length is less. The signals x and y have been analyzed in a given time interval, and only N points of each signal have been taken to calculate the crosscorrelation function between them. Consequently, the obtained function is just an estimate, based on the two given segments of the processes x and y. This affects the obtained function when calculating for greater shifts k (61). Of course, the longer the time intervals (respectively, the greater the number of points N), the better the estimate of the crosscorrelation function. Time shifts measured from the crosscorrelation functions calculated from successive time intervals of the two studied signals can be presented as a function of the recording time in a time shift plot. 4.3.2. ANALYSIS OF THE WAVEFORM OF EGG. Familoni at al. (39) suggested that EGG waveform analysis could be used for indirect assessment of direction of propagation. They pointed out that in normal subjects the descending portion of EGG waves usually dominates, which probably indicates distal direction of propagation. During periods of uncoupling the direction of propagation is disturbed. Therefore uncoupling could be reflected by a change in EGG waveform. However, no quantitative evaluation of this suggestion has been made. The simplest method to study the dynamics of EGG waveforms in a certain time interval is to divide the number of points with smaller amplitudes than their predecessor by the number of points with greater amplitudes than their predecessor in the whole interval. A ratio greater than unity would indicate that the descending arms of the waves dominated in this interval. These ratios could be obtained for successive time intervals and arranged versus time in a gradient plot. 4.4. METHODS TO ASSESS GASTRIC ELECTRICAL FREQUENCY. 4.4.1. QUANTITATIVE VISUAL ANALYSIS. This is simply visual counting the number of waves in a certain time interval. It is quite reliable when dealing with SDB signals, but some lower amplitude or irregular LDB and EGG waves can easily be missed. This reduces the reliability of the method in comparative studies. 4.4.2. FREQUENCY ANALYSIS AND TIME-FREQUENCY PLOTS. While counting SDB waves in the time domain can easily be implemented by computer, a little bit more sophisticated mathematics is required for the analysis of LDB and EGG signals in the frequency domain. Fourier Transforms The pair of Fourier transforms (63) converts signals that are functions of one variable (time or frequency) into functions of the other (frequency or time). The two variables are independent (63). The two transforms are not entirely similar, because of the different sign in the exponent. It is quite clear that transforming the the signal V(t) from the time domain into the frequency domain to obtain F(f), the integration should be done along the whole time interval in which the function V(t) is defined. In general, if the function is continuous the integration is from -infinity to +infinity. However, if the integral is to be solved practically, some real limits are required. Although these limits can be made very large, they can never reach infinity. Therefore, the obtained function F'(f) when the integral is limited is an estimate of the actual function F(f). The phenomenon that comes as a result of this is called "leakage" (64) and will be considered later. On the other hand even if V(t) is a real signal, F(f) is, in general, a complex function which can be expressed as [65]:
where |F(f)| is the magnitude (or amplitude) spectrum of V(t), and f(f) is the phase spectrum. The square of the magnitude spectrum is known as the power spectrum. These three entities are very useful for practical evaluation of F(f). Hartley Transforms Recently (65, 66) a new pair of transforms was suggested as a real alternative to the complex Fourier transforms. The idea for the Hartley transforms has been suggested by Ralph Hartley in 1942, and it was more recently further by Bracewell (65). The two transforms are defined in (65). They constitute time-to-frequency and frequency-to-time transforms. In addition, they are real and the function multiplier in the integral is the same. These two important qualities distinguish the two Hartley transforms making them easier to implement in practice with one and the same procedure. The problem is to define the relationship between the two functions H(f) and F(f) and respectively that between Fourier and Hartley transforms. Relationship between Fourier and Hartley Transforms. Bracewell (65) points out that if H(f) = E(f) + O(f), where E(f) and O(f) are the even and odd parts of H(f), then:
It is then easy to notice that:
Consequently, the Fourier transform is the even part of the Hartley transform minus j times the odd part. The Hartley transform is the real part of the Fourier transform minus the imaginary part. It is obvious that one set of Hartley coefficients can produce one and only set of Fourier coefficients and vice versa. Power (and respectively magnitude) spectra can easily be obtained from the Hartley transform using
Phase spectrum can be found with
Convolution and Correlation using Frequency Transforms. It is well known that frequency transforms can be used for calculating convolution and correlation, features that are defined in time domain but have easy to handle representations in frequency domain (61-65) The image of the crosscorrelation function in the frequency domain is usually called the cross power spectrum [62]. Since limiting the integral of Fourier is nothing else but a time-domain multiplication of the signal with a rectangular window, the leakage effect mentioned before can now be redefined as the effect of convolution in the frequency domain between the desired spectrum and the spectrum of the rectangular limiting window. As a result of this convolution, spectral components would spread away from the correct frequency, resulting in an undesirable modification of the actual spectrum of the signal (62, 64). The impact of the leakage can be reduced (if it is critical to the processing procedures), using limiting windows that are not rectangular (62, 64). Top of the page Bottom of the page Selected Topics on Electrogastrography Discrete Transforms. The transition from an integral to a discrete Fourier or Hartley transforms has been discussed in many studies (61, 62, 63, 64, 65). One of the most important effects of sampling an analog signal is so called "aliasing". It occurs as a direct consequence of not complying with the Nyquist's (or Shannon's) sampling theorem (62, 64), which states that
where Fs is the sampling frequency, Ts - the sampling interval, and Fh is the highest frequency contained in the signal. If the sampling is done ignoring equation [4.6], the high frequency components of the signal will show up on the left side of the spectrum, thus misleading the investigator that some low frequency components are present. The negative impact of this effect can be avoided by appropriate low pass filtering before the sampling, as well as by choosing an adequate sampling frequency. One of the main properties of the discrete frequency transforms is that they repeat themselves in the frequency domain at intervals of N points. This symmetry implies that the multiplications between the functions in frequency domain will result in their circular convolution in time domain, which might not always be desired, especially when the signal is not periodic. If linear convolution is to be done, the signals in the time domain with a length of N samples each, have to be padded with zeros up to the (2N-1)th sample prior to the frequency conversion (see Challis and Kitney, 62). Fast Algorithms. Calculating DFT (or DHT) directly is a tedious work which involves very many multiplications and thus a large amount of computer time. That is why fast algorithms have been suggested first for the DFT (62, 63, 64) and recently for the DHT (65, 66). These algorithms are known as the fast Fourier transform (FFT) and the fast Hartley transform (FHT). There are two basic groups of algorithms, each of which decomposes the sampled signal in a particular way. The first is based on a decimation in frequency, which decomposes the original signal into sections of sequential samples, until a number of two-sample signals are produced (62). For a two-sample signal the transform is very easy:
In order that a digital signal of length N be decomposed into N/2 two-sample DFTs, N should be a power of 2. The second main group of algorithms is based on decimation in time. The samples that form every new group of samples are not sequential, but alternate. The decomposition results again in N/2 pairs of samples, which are then transformed using equation [4.18]. Most FFTs require log2N.N/2 complex multiplications and log2N.N/2 complex additions, while DFT requires N2 complex multiplications and N.(N-1) complex additions. Theoretically FHT requires the same number of operations as the FFT. However, the fact that the two Hartley transforms are identical and real gives the opportunity to cut significantly the real operating time of the transforms. Bracewell (65) points out that if Fourier coefficients are not to be computed, but instead the power spectrum is calculated directly from the Hartley coefficients using
the actual computing time is always less when using the FHT. More details can be found in Chapter 8 of Bracewell's book "The Hartley Transform" (65), where a careful comparison between the actual computing times for FFT and FHT is made. Running Spectrum Analysis. From the previous discussion it is quite clear that the frequency spectrum of a time segment taken from the analyzed signal represents the distribution of frequencies present in this particular segment. If the segment is very big with respect to the period of the signal and this period changes for a small amount of time within this huge time segment, it is quite possible that in frequency domain this change will not be clearly seen, especially if other noise signals are also present. Therefore, a reasonable compromise between the understandable desire to make N as big as possible (to reduce the leakage effect) and the above considerations should be made. It has been widely accepted EGG signals to be analyzed with 512-point FFT(FHT) and the signal to be sampled with a sampling frequency of 2 Hz (71, 72). This indicates that frequency transforms are to be performed on time segments with a length of 256 seconds (4.27 minutes). Frequency resolution in this case is 0.003906 Hz or 0.234 cpm. Obviously, in order to follow the dynamics of the frequency components in the whole recording, spectra from successive 4.27 min. intervals should be obtained. Kingma et al. (67) proposed to arrange these EGG spectra one behind another with a small shift, forming the so called "three dimensional" (magnitude/power-frequency-time) plots. This technique was later used by Van Der Schee et al. (68, 69) and became known as running spectrum analysis. There are two ways of presenting successive spectra in a three dimensional plot: (a) normalize each spectrum separately to 100% in a given frequency window; (b) find the absolute maximum of all spectra and normalize it to 100%. Both methods have advantages and disadvantages. When analyzing the frequency only, the first method is preferrable if huge differences in the signal amplitudes are expected. On the other hand, when the differences in the amplitudes are to be investigated the second option should be available. Instead of taking successive time intervals, it is possible to overlap them, e.g. if the first spectrum is obtaind for the time interval (0-255 sec), the second could be calculated in the interval (128 - 383 sec) instead in (256 - 511 sec). This introduces a 50% overlap between the successive time segments. The overlap allows frequency dynamics of the signal to be followed better, but it also increases the number of spectra in the three dimensional plot. A sample of a 1-hour three dimensional plot of a sine wave with an initial frequency of 0.05 Hz, which changed to 0.07 Hz after the first half hour of recording is shown on Fig.4.4.
Figure 4.4. One-hour three dimensional plot of a sine wave signal which changes its frequency from 0.05 Hz to 0.06 Hz after the first half hour of recording. The frequency window is 0.01-0.16 Hz . Grey Scale and Contrast Plots. Normalizing the absolute maximum of all spectra in the three dimensional plot to a 100% can be equivalent to assigning it a black color. On the contrary, white color could then represent the absolute minimum in the three dimensional plot. All spectral magnitudes in between the two extremes can be represented by the levels of gray color. Thus, one of the dimensions in the three dimensional plot, the magnitude, can be eliminated. In the new two dimensional plot the magnitudes are represented as pixels with different levels of gray from white to black. This plot was first suggested by Van Der Schee et al. (68, 69) and was given the name "gray scale plot". One convenient modification of the gray scale plot for the needs of EGG is the contrast plot (70), in which only two colors are present - black and white. When a given spectral component exceeds a predetermined magnitude level (specified in percents), it is represented in the plot with a black pixel. The pixels remain white if the spectral magnitude is below this level. Time-Frequency Plots. The spectra from Fig.4.5 are individually normalized to 100%. They explicitly show dominant peaks in the given frequency window. In the time-frequency plot (71), the maxima of the dominant peaks from each spectrum are arranged in a plot versus time and connected with lines, i.e. the amplitude information is completely eliminated (Fig. 4.5). It is quite obvious that (a) time-frequency plots are two-dimensional; (b) they are built by points, each of which represents the frequency of the dominant peak in the corresponding spectrum; (c) since the spectra in the three dimensional plot represent successive time intervals, the points that build up the time frequency plot are shifted with this interval; (d) when overlap is introduced, the interval TI between any two successive points in the time-frequency plot is given with:
where OVL is the percentage of overlap used.
Figure 4.5. Time-frequency plot extracted from the three dimensional plot shown on Fig. 4.5. with frequency scale in cycles per minute (cpm). The moment of the change in frequency is clearly seen. 4.4.3. STATISTICAL EVALUATION OF TIME-FREQUENCY PLOTS. For every signal interval t transformed into frequency domain, there is only one dominant frequency peak and respectively only one point in the time-frequency plot. If a signal with duration Ds is digitized with sampling interval Ts, the number of points in the time-frequency plot Ntf is given with:
where N is the number of points of the transform t = N.Ts and OVL is the percentage of overlap between the successive time intervals. The points are later arranged versus time and connected with lines to form the actual time-frequency plot. These points can be subjected to various methods of simple statistical evaluation (61, 72, 73). Mean value, variance and standard deviation were called by Challis and Kitney (61) "basic averages". Mean value Mn is defined as (61):
where Fi are the frequency values of the points in cycles per minute (cpm). The variance Vr can be calculated as:
The standard deviation SD is the square root of the variance and when assessing time-frequency plots is measured in cycles per minute (cpm), similar to the mean value (61, 72). Another way of evaluating the points from the time-frequency plot is to build the so called "probability density function" (61, 73). This function measures the probability (in %) of one unknown process or function to have a given value. In the simplest case the probability density function (pdf) of a sine wave with a period of 20 seconds would be a single 100% line at a frequency of 3 cpm. However, if the period of the sine wave changes with the time, many different frequencies would be present and the pdf would become a curve. When processing digital signals it is convenient to represent the pdf curve as a histogram. Each horisontal step in this histogram is equal to the frequency resolution of the time-frequency plot, while vertically the percentage of the points associated with the given frequency is shown. Similar statistical evaluations can be used to assess time shift and gradient plots. 4.5. DIGITAL FILTERING OF GEA. 4.5.1. FREQUENCY-SAMPLING FILTERS. The theory of frequency-sampling filters is well described in the literature (63, 64). They can be designed directly in the frequency domain. The designer simply specifies the value one in the passband and zero in the stopband without providing any transition band. Although this simple filter has a finite impulse response (FIR), i.e. its phase characteristic, is linear, its response may display significant ripple between frequency sample points. This difficulty can be minimized by providing a certain transition band, i.e. there by including transition points between the last unity value in the passband and the first zero value in the stopband. However, there is not a definite rule what the width of this transition band should be in a given case (64). Frequency sampling filters with a programmable transition band based on the fast Hartley transform will be used in this study. The linearity of the phase characteristic of these filters make them also applicable for time shift studies of EGG. 4.5.2. ADAPTIVE FILTERS. Although adaptive filtering was introduced in the early sixties, only recently has it been applied to canine GEA signals (74) and human (75) GEA signals. There has not been dramatic success in either subject. Figure 4.6 shows the typical block-diagram of an adaptive filter as suggested by Widrow et al. (76). Several simple equations can be written to describe this system (74):
which after squaring gives
Taking the expectation values (the power) at both sides and presuming that s is uncorrelated with n0 and y, it can be shown that
Figure 4.6. Block diagram of an adaptive filter. Recorded signal d is considered a superposition of the primary signal s and the noise n0, which are not correlated.. The reference signal n1 is a noise, correlated with n0. The objective of the filtering process is to extract the primary signal s. Consequently, if the output signal power E(z*z) is minimized the signal s will not be affected, but the noise will be eliminated. The minimization process, known as Widrow-Hoff least mean squares (LMS) algorithm (76) is based on the adaptive adjustment of the coefficients w, with which the reference noise n1 is multiplied to produce the signal y. Unfortunately, when recording GEA signals it is very difficult to find a pure source of noise. That is why the block-diagram from Fig. 4.6 has been adopted for the needs of electrogastrography by Kentie et al. (74) and later utilized by others (75, 77). Instead of using two channels, Kentie et al. (74) suggest using one channel from which the noise channel is derived by a band rejection of the EGG component (Fig. 4.7). This modification, however, is not quite legitimate, because the quality of the noise n1, and especially its relations to the signal s are now dependent on the quality of the band reject filter used. In that sense, converting the band reject filter into a bandpass filter and using it directly to eliminate the noise from signal d would have exactly the same result (10).
Figure 4.7. A modified adaptive filter for EGG. The noise n1 is obtained from the recorded signal d using a band reject filter. The only impact that this modified adaptive algorithm may have on the output signal is the unknown time shift that it would introduce. In general, the time shift introduced by the adaptive filters would differ for different EGG channels and it is not surprising that Chen et al.(41) discovered occasional large time shifts between different EGG channels after using an adaptive algorithm. Top of the page Middle of the page Selected Topics on Electrogastrography
|
|
|