text:
# DIRECT Distances to Nearby Galaxies Using Detached Eclipsing Binaries and Cepheids. VII. Additional Variables in the Field M33A Discovered with Image Subtraction1footnote 11footnote 1Based on observations obtained with the 2.1m telescope at the Kitt Peak National Observatory.

## 1 Introduction

Starting in 1996 we undertook a long term project, DIRECT (i.e. “direct distances”), to obtain the distances to two important galaxies in the cosmological distance ladder, M31 and M33. These “direct” distances will be obtained by determining the distance of Cepheids using the Baade-Wesselink method and by measuring the absolute distance to detached eclipsing binaries (DEBs). While the cosmological distance scale has been the subject of numerous recent observation campaigns, especially those enabled by the Hubble Space Telescope (HST) and massive variability studies of the Magellanic clouds, M33 has not been re-surveyed since the photographic survey of Kinman, Mould & Wood (1987).

M31 and M33 are the stepping stones to most of our current effort to understand the evolving universe at large scales. First, they are essential to the calibration of the extragalactic distance scale. Second, they constrain population synthesis models for early galaxy formation and evolution and provide the stellar luminosity calibration. There is one simple requirement for all this—accurate distances. These distances are now known to no better than 10-15%, as there are discrepancies of $`0.20.3\mathrm{mag}`$ between various distance indicators (e.g. Huterer, Sasselov & Schechter 1995; Holland 1998; Stanek & Garnavich 1998).

DEBs have the potential to establish distances to M31 and M33 with an unprecedented accuracy of better than 5% and possibly to better than 1%. Detached eclipsing binaries (for reviews see Andersen 1991; Paczyński 1997) offer a single step distance determination to nearby galaxies and may therefore provide an accurate zero point calibration of various distance indicators – a major step towards very accurate determination of the Hubble constant, presently an important but daunting problem for astrophysicists. DEBs have been recently used to obtain accurate distance estimate to the Large Magellanic Cloud (e.g. Guinan et al. 1998; Udalski et al. 1998).

The detached eclipsing binaries have yet to be used as distance indicators to M31 and M33. According to Hilditch (1996) there was only one eclipsing binary of any kind known in M33 (Hubble 1926). The recent availability of large-format CCD detectors and inexpensive CPUs has made it possible to organize a massive search for periodic variables, which will produce a handful of good DEB candidates. These can then be spectroscopically followed-up with the powerful new 6.5-10 meter telescopes.

The study of Cepheids in M33 has a venerable history (Hubble 1926). Freedman, Wilson & Madore (1991) obtained multi-band CCD photometry of some of the Cepheids discovered in photographic surveys, to build a period-luminosity relations in M33. However, the sparse photometry and the small sample (11 Cepheids) do not provide a good basis for obtaining direct Baade-Wesselink distances (see, e.g., Krockenberger, Sasselov & Noyes 1997) to Cepheids—the need for new digital photometry has been long overdue.

As the first step of the DIRECT project we have searched for DEBs and new Cepheids in the M31 and M33 galaxies. We have analyzed five $`11\mathrm{}\times 11\mathrm{}`$ fields in M31, A-D and F (Kaluzny et al. 1998, 1999; Mochejska et al. 1999; Stanek et al. 1998, 1999; hereafter Papers I, IV, V, II, III). A total of 410 variables, mostly new, were found: 48 eclipsing binaries, 206 Cepheids and 156 other periodic, possible long-period or non-periodic variables. We have also analyzed two fields in M33, A and B (Macri et al. 2001a; hereafter Paper VI) and found 544 variables: 47 eclipsing binaries, 251 Cepheids and 246 other variables.

As a second step, we started followup observations of selected DEBs in both the M31 and M33 galaxies with bigger telescopes in order to construct more precise and well sampled light curves for them.

In this paper, seventh in the series, we present a catalog of variable stars found in the same field as the detached eclipsing binary D33J013346.2+304439.9 using followup observations collected at the Kitt Peak National Observatory 2.1m telescope. The paper is organized as follows: Section 2 provides a description of the observations. The data reduction procedure is outlined in Section 3. The catalog of variable stars is presented in Section 4, followed by its brief discussion in Section 5. Section 6 deals with the first overtone Cepheid candidates. The concluding remarks are stated in Section 7.

## 2 Observations

The data discussed in this paper was obtained at the Kitt Peak National Observatory222Kitt Peak National Observatory is a division of NOAO, which are operated by the Association of Universities for Research in Astronomy, Inc. under cooperative agreement with the National Science Foundation. 2.1m telescope equipped with a Tektronix $`2048\times 2048`$ CCD (T2KA camera) having a pixel scale $`0.305\mathrm{}/pixel`$ during two separate runs, from September 29th to October 5th, 1999 and from November 1st to 7th, 1999. The primary observing targets were three detached eclipsing binaries, one in each of the fields M33A, M33B and M31A, discovered previously as part of the DIRECT project (Papers VI and II). For field M33A we collected $`158\times 600s`$ exposures in the $`V`$ filter and $`62\times 600s`$ in the $`B`$ filter.333The complete list of exposures for this field and related data files are available through anonymous ftp on cfa-ftp.harvard.edu, in pub/kstanek/DIRECT directory. Please retrieve the README file for instructions. Additional information on the DIRECT project is available through the WWW at http://cfa-www.harvard.edu/~kstanek/DIRECT/. The exposure times varied slightly to compensate for the changes of seeing conditions. The typical seeing was $`1.\mathrm{}4`$. The field was observed through airmass ranging from 1 to 1.9, with the average at 1.2. The completeness of our data starts to drop rapidly at about 21.2 mag in $`V`$ and 21.8 mag in $`B`$, judging from the magnitude distributions of the variable stars (Fig. 1).

## 3 Data reduction

### 3.1 Correction for the non-linearity of the T2KA camera

During the second night of the first observing run we compared the photometry from two images of NGC 7789 with differing exposure times and found that the camera had a non-linear response. To lessen the impact of non-linearity we reduced the gain setting from 5.1 to 3.6 $`e^{}`$/ADU.

Several sets of calibration frames were collected at both gain settings for the purpose of correcting the data for the non-linearity. Each set consisted of a series of dome flats with intensities from $`500`$ ADU to $`60000`$ ADU. After subtracting the readout bias under IRAF444IRAF is distributed by the National Optical Astronomy Observatories, which are operated by the Association of Universities for Research in Astronomy, Inc., under cooperative agreement with the NSF. the mode intensity $`I_i`$ (in ADU) was computed for each calibration frame.

We adopted the correction function for the non-linearity in the same form as in the IRAF task IRLINCOR:

$$I_e=I_i(c_1+c_2\frac{I_i}{32767}+c_3(\frac{I_i}{32767})^2)$$
(1)
where $`I_i`$ is the mode observed intensity on the calibration frame, $`I_e`$ the intensity corrected to the expected level and the coefficients $`c_1,c_2,c_3`$ are to be determined. The exposure time, $`t`$ corrected for the shutter speed (0.1 s), was used to calculate the expected intensity $`I_e`$ from the measured intensity level $`I_i`$:

$`I_e=I_i{\displaystyle \frac{t}{t_0}}`$
By definition the detector response was assumed to be linear at the intensity level of the image with the shortest exposure time $`t_0`$. The coefficients $`c_1`$, $`c_2`$ and $`c_3`$ were determined from a least-squares fit of Eq. 1 to the measured and expected intensity levels $`I_i`$ and $`I_e`$. We list them in Tab. 1 for the two gain settings used.

Figure 2 shows the difference between the expected intensity $`I_e`$ and the measured intensity $`I_i`$ as a function of the measured intensity for gain 3.6$`e^{}`$/ADU. Each point corresponds to a calibration exposure. Exposures from the same set are marked with like symbols. The points exhibit substantial scatter between the sets and in some cases within the set (ie. open stars and open triangles). Several fits were made to different combinations of sets. The derived coefficients were then tested on pairs of frames of the same object with exposure times differing by a factor of at least 10. The fit corresponding to the best coefficients is shown with a dotted line on Fig. 2. The calibration frames used in the fit are indicated by filled symbols.

A comparison between the photometry derived from a 100s and a 10s exposure of the NGC 7789 open cluster is presented in Fig.3. The left panel shows a comparison for the uncorrected exposures and the right panel for the corrected ones. A significant improvement is seen between the uncorrected and corrected frames. A close inspection of the right panel of the figure will reveal that the non-linearity has not been removed completely.

The standard preliminary processing of the science data was performed with the routines in the IRAF CCDPROC package. The non-linearity was corrected with the IRLINCOR task using the derived coefficients.

### 3.2 Nonlinearity of the T1KA CCD

During an observing run in October 1998 at the KPNO 2.1 meter telescope we have collected data using the T1KA $`1024\times 1024`$ CCD. We have compared a series of images of an open cluster, Be 32, which were obtained one after another in photometric conditions (Fig. 4). The images vary only in exposure times (starting at 6 seconds and going up to 180 seconds). We have found that the magnitudes (obtained via point spread function fitting + aperture correction) do not differ by the expected factor of 2.5 log (time), but rather show increasing deviations from the expected scaling as the difference in exposure time gets larger. Furthermore, the relations show strong non-linearity, which increases with the difference in exposure time. This fact was not known to us during the time of observations and no calibration frames were collected to correct the data. We provide this information to warn the astronomical community of the problem we have encountered with this camera.

### 3.3 Photometry and Selection of Variables

The photometry for the variable stars was extracted using the ISIS image subtraction package (Alard & Lupton 1998, Alard 2000a) from the $`V`$ and $`B`$-band data.

The ISIS reduction procedure consists of several steps. Initially all of the frames are transformed to a common coordinate grid. Next a reference image is created by stacking several frames with the best seeing. For each frame the reference image is convoluted with a kernel to match its point spread function (PSF) and then subtracted. On the subtracted images the constant stars will cancel out and only the signal from variable stars should remain. A median image is constructed of all the subtracted images and the variable stars are identified as bright peaks on it. Finally, profile photometry is extracted from the subtracted images.

All of the computations were performed with the frames internally subdivided into four sections (the ISIS parameters sub\_x and sub\_y set to 2). The differential brightness variations of the background were fit with a second degree polynomial (deg\_bg=2). A convolution kernel varying quadratically with position was chosen (deg\_spatial=2).

An image of particularly good quality was selected as the template frame for the stellar positions. The remaining images were re-mapped to the template frame coordinate system using a third degree polynomial transform. During this step an initial rejection of cosmic rays was also performed. A setting of 1.0 for the cosmic ray threshold (COSMIC\_THRESH) was used.

A reference image was then constructed from 25 images with the best seeing. The minimum, maximum and average seeing of these 25 images were $`0.95\mathrm{}`$, $`1.23\mathrm{}`$ and $`1.11\mathrm{}`$ in $`V`$ and $`0.93\mathrm{}`$, $`1.28\mathrm{}`$ and $`1.13\mathrm{}`$ in $`B`$. The constituent images were transformed to the same seeing and background level as the best image (the template in our case) and stacked by taking a median in each pixel to obtain the reference image virtually free of cosmic rays.

Image subtraction was then applied to all the frames. For each frame the reference image was convoluted with a kernel to match it as closely as possible and then the frame was subtracted from it. As the flux of the non-variable stars on both images should be almost identical, such objects will disappear from the subtracted image. The only remaining signal will come from variable stars.

To obtain a better signal-to-noise ratio a median of the subtracted images was computed for the purpose of detecting variable stars. The variables were then identified as bright peaks on this image. In practice the situation is slightly complicated by the presence of stars saturated on the reference image, as their profiles were not matched exactly to those on frames with inferior seeing. An easy way to distinguish genuine variables from saturated stars is by their profile on the subtracted image: it should be smooth for the variables and usually contains a flat top or a dip at the center for the saturated stars. The variables were identified visually on the median subtracted image and their positions recorded with the IMEXAMINE task under IRAF. The final positions were derived with the daophot FIND subroutine (Stetson 1987). Variables were identified in both bands. The lists were then merged into a master variable list.

The light curves of the variables (in ADU) were extracted from the subtracted images. Preliminary periods were found with the ISIS czerny subroutine, based on the analysis of variance method, as described by Schwarzenberg-Czerny (1989), and later refined using the same algorithm. Stars with chaotic and/or noisy light curves were rejected from the sample upon visual inspection.

In the next step the light curves were converted into magnitudes. PSF-fitting photometry was extracted with the daophot/allstar package from the template image following the approach described in Paper I. Some of the variable stars were not detected on the template, in most cases because of their faintness or crowding. As light curves of these stars could not be converted to magnitudes, they were left expressed in counts.

The template instrumental magnitude $`m_{tpl}`$ of each variable was converted into counts $`c_{tpl}`$, using the Allstar zeropoint of 25 mag. The light curve was then converted point by point to magnitudes $`m_i`$, by computing the total flux $`c_i`$ for the given epoch as the sum of the counts on the template $`c_{tpl}`$ and the counts on the subtracted template image $`\mathrm{\Delta }c_{tpl}=c_{ref}c_{tpl}`$ decreased by the counts corresponding to the subtracted image $`\mathrm{\Delta }c_i=c_{ref}c_i`$.

$$c_i=c_{tpl}+\mathrm{\Delta }c_{tpl}\mathrm{\Delta }c_i$$
(2)
The flux was converted to instrumental magnitudes, using the same zero point as above. To convert the photometric error $`\sigma _{c,i}`$ from counts to $`\sigma _{m,i}`$ magnitudes, we used the following relation:

$$\sigma _{m,i}=2.5\mathrm{log}(\frac{c_i}{c_i+\sigma _{c,i}})$$
(3)
### 3.4 Photometric Calibration and Astrometry

On the nights of October 3/4 and November 3/4 1999 we observed a total of 75 stars in 13 Landolt (1992) fields. We attempted to derive a calibration from these exposures, corrected for non-linearity. The o-c residuals showed a clear dependence on magnitude, especially in the $`V`$-band, where they exhibited a very similar slope on both nights. Several new fits to Eq. 1 with different combinations of data sets were tried. Those that gave some improvement in the calibration, produced noticeably inferior results when the photometry on two images of very different exposure times was compared. We had decided to stay with the initial values for the correction coefficients (Tab. 1) and calibrate our photometry using a different method. The coefficients for the color terms were derived from the comparison of our NGC 6791 photometry with the data from the KPNO 0.9m telescope (Kaluzny & Udalski 1992). The offsets were determined relative to the DIRECT catalog of stellar objects in M33 (Macri et al. 2001b). The following transformations were adopted:

$`v=V5.491+0.039(BV)`$
$`bv=0.178+0.927(BV)`$
$`b=B5.3120.034(BV)`$
The instrumental $`V`$ and $`B`$-band light curves of the variables were transformed to the standard system by adding the appropriate offsets, as the coefficients next to the color terms are very small.

Equatorial coordinates were determined for the $`V`$ template star list, expanded with the variables with no $`V`$-band photometry. The transformation from rectangular to equatorial coordinates was derived using 818 transformation stars with $`V<19.5`$ from the DIRECT catalog of stellar objects in M33 (Macri et al. 2001b). The average difference between the catalog and the computed coordinates for the transformation stars was $`0.\mathrm{}07`$ in right ascension and $`0.\mathrm{}06`$ in declination. As a check we have also obtained an independent transformation to the USNO-A2 catalog (Monet et al. 1996) using 56 transformation stars ($`\mathrm{\Delta }\alpha =0.\mathrm{}26`$, $`\mathrm{\Delta }\delta =0.\mathrm{}23`$). The average differences between the computed coordinates and the M33 catalog for 1334 stars with $`V<19.5`$ were $`\mathrm{\Delta }\alpha =0.\mathrm{}71`$ and $`\mathrm{\Delta }\delta =0.\mathrm{}60`$.

## 4 Catalog of Variables

### 4.1 Classification

The variables we are most interested in are Cepheids and eclipsing binaries (EBs). We therefore searched our sample of variable stars primarily for these two classes of variables. The variable stars were preliminarily classified as eclipsing, Cepheid or miscellaneous by visual inspection, based on the shape of their light curves. The variables for which neither $`V`$ nor $`B`$-band magnitude could be determined (with only flux light curves or having a period in excess of 14 days) were not reclassified further.

In order to obtain as clean a sample of Cepheids as possible, we have inspected the location of the Cepheid variable candidates on a $`V/BV`$ CMD. All of the Cepheid candidates having highly discrepant colors were reclassified as other periodic variables. The candidates with light curves in only one of the bands were checked on the period-luminosity relation for that band. Extreme outliers were also moved to the other periodic variable category.

The EBs with light curves expressed in magnitudes for at least one band, were fitted with a model described in Papers I and II. Within our assumption the light curve of an EB is determined by nine parameters: the period, the zero point of the phase, the eccentricity, the longitude of periastron, the radii of the two stars relative to the binary separation, the inclination angle, the fraction of light coming from the bigger star and the uneclipsed magnitude. Eclipsing binary candidates for which a satisfactory fit was not achieved, were reclassified as other periodic variables.

In the following sections 4.2-4.5 we present the parameters and light curves for the 434 identified variable stars.555The $`BV`$ photometry and $`V`$ finding charts for all variables are available from the authors via the anonymous ftp from the Harvard-Smithsonian Center for Astrophysics and can be also accessed through the World Wide Web. All tables are sorted by increasing period, unless otherwise noted. We adopt a naming convention after Macri et al. (2001a) based on the J2000.0 equatorial coordinates, in the format: D33J$`hhmmss.s`$+$`ddmmss.s`$. The first three fields ($`hhmmss.s`$) correspond to right ascension expressed in hours, the last three ($`ddmmss.s`$) to declination, expressed in degrees, separated by the declination sign. As none of the newly discovered variables are present in previous variable star catalogs, we refer the reader to Tables 5 and 6 in Paper VI for cross-identifications.

### 4.2 Eclipsing Binaries

We have found a total of 63 eclipsing binaries in field M33A. In Table 2 we present the parameters for the 53 EBs with a magnitude light curve in at least one band. For each variable we list its name, period P, magnitudes $`V_{max}`$ and $`B_{max}`$ of the system outside of the eclipse, and the radii of the binary components $`R_1,R_2`$ in the units of the orbital separation. We also give the inclination angle of the binary orbit to the line of sight $`i`$ and the eccentricity of the orbit $`e`$. The reader should bear in mind that the values of $`V_{max},B_{max},R_1,R_2,i`$ and $`e`$ are derived with a straightforward model of the eclipsing system, so they should be treated only as reasonable estimates of the “true” value. Table 3 lists 10 EBs with flux light curves only. For each variable we give its name and period P. Figure 5 presents the phased light curves of 12 sample EBs (see also Table 8).

### 4.3 Cepheids

A total of 305 Cepheid variables were found in field M33A. In Table 4 we present 242 Cepheids with a magnitude light curve in at least one band. For each variable, we list its name, period P, flux-weighted average magnitudes $`V`$ and $`B`$ and the $`V`$ and $`B`$-band amplitudes $`A_V`$ and $`A_B`$. Due to the short time base of our observations, reliable flux-weighted magnitudes could only be determined for variables with periods shorter than 14 days. We have extracted light curves for 36 Cepheids with longer periods, all of them identified previously in Paper VI and made them available via anonymous ftp. In Table 5 we list 27 Cepheids with flux light curves only. For each variable we list its name and period P. Figure 6 presents the phased light curves of 24 sample Cepheids (see also Table 9).

### 4.4 Other Periodic Variables

In Table 6 we present the parameters of 13 periodic variables. For each variable we list its name, period P, the magnitudes $`V`$ and $`B`$ and the $`V`$ and $`B`$-band amplitudes $`A_V`$ and $`A_B`$. The $`V`$ and $`B`$ columns list the magnitudes outside of the eclipses $`V_{max}`$ and $`B_{max}`$ for the eclipsing variables and flux-weighted average magnitudes $`V`$ and $`B`$ for the other variables. We have also found three other variables with periods of the order of 50 days (identified previously in Paper VI), which we do not list in the table. Their light curves are available via anonymous ftp. In Table 10 we list the light curves of all periodic variables.

### 4.5 Miscellaneous Variables

In Table 7 we present the parameters of 50 miscellaneous variables. For each variable we list its name, the average magnitudes $`\overline{V}`$ and $`\overline{B}`$ and the $`V`$ and $`B`$-band amplitudes $`A_V`$ and $`A_B`$. Variable D33J013342.3+304024.7 is probably an EB, with incomplete phase coverage. Several of those variables are most likely periodic, with periods in excess of 14 days. In Table 11 we list the light curves of all miscellaneous variables.

### 4.6 Comparison with the catalog in Paper VI

In Figure 7 we compare our $`V`$ and/or $`B`$ photometry for 81 Cepheids (circles) and 18 EBs (squares) with the values listed in the DIRECT catalog of variables in M33 (Paper VI). As reference we plot in the background similar comparisons for all the stars in the field (dots). The variables, with some exceptions, have $`\mathrm{\Delta }B`$ and $`\mathrm{\Delta }V`$ distributions roughly similar to the rest of the stars, although they do show a marked tendency to be fainter in our catalog. This trend is especially prominent in the $`V`$-band for the EBs, but is absent in their $`B`$-band comparison. We have examined the images of three of the more discrepant EBs and found that in our data we have resolved more nearby companions to these variables. These companions are brighter in $`V`$ than in $`B`$, relative to the EBs. This is not surprising, as the EBs are among of the bluest stars in the field. There are also a few faint Cepheids, which show unusually large differences in photometry, of the order of 0.6-0.8 mag. We have inspected their light curves in both catalogs. The $`V`$-band light curves of one such star, type II Cepheid D33J013356.5+304632.1, are shown in Fig. 8. The left panel shows our light curve, extracted with ISIS from 2.1m telescope data. The right panel shows the daophot light curve from data taken with a 1.3m telescope (Paper VI). It seems that these discrepancies are caused in large part by the fact that fixed position photometry is prone to identify and fit a profile at the supplied position even when the star is below the detection threshold, resulting in the false measurement of a fainter magnitude.

In Figure 9 we plot a histogram of the surface density of Cepheids (N/kpc2) as a function of their period, assuming a distance of 840 kpc to M33 (Freedman, Wilson & Madore 1991). The solid line represents all Cepheids from our catalog with $`P<14`$ days, the dashed line – the Cepheids from Paper VI. The areas covered in this search and in Paper VI are 108 and 222 arcmin2. As our observations were carried out during two one-week runs spaced one month apart, we lack the baseline to detect long period Cepheids. Since our data were collected with an instrument 2.6 times larger in area and with better seeing, we have a higher detection rate for short period Cepheids. The Paper VI catalog, due to the much longer baseline of observations, contains more long period Cepheids.

## 5 Discussion

In Fig. 10 we plot the positions of the variable stars on the $`V/BV`$ CMD. The EBs are denoted by open squares and Cepheids by filled circles in the left panel, the periodic variables by filled squares and miscellaneous by open circles in the right panel.

All of the EBs occupy the upper main sequence, with the exception of D33J013359.5+304037.6 ($`BV=0.37`$). This variable most likely suffers from larger than average reddening, as it is located at the edge of what appears to be a dust lane.

Most Cepheids occupy the region between $`0.3<BV<1.0`$, with several outliers stretching from $`BV=0.10`$ to $`BV=1.24`$. The discrepant colors of some of the Cepheid variables most likely originate from the phenomenon of blending with nearby bright stars of very different colors and/or reddening. The phenomenon of blending occurs when the Cepheid possesses one or more close companions which cannot be separated at the resolving power of the instrument used. A further discussion of blending and its properties can be found in Mochejska et al. (2000; 2001) and Stanek & Udalski (1999).

All but one of the periodic variables with $`BV`$ photometry are bluer than $`BV=0.35`$. The five brightest ones exhibit Cepheid-like light curves of very low amplitude of the order of a few hundredths of a magnitude. The low amplitudes of their variability suggest that these could be Cepheids blended with blue main sequence stars. The remaining four periodic variables are likely EBs. Three of them have blue colors, consistent with this type of variable. The fourth one, D33J013346.1+304658.9, with a period of 6.4640 days seems to be a contact EB, despite its unusual color ($`BV=0.511`$). A more complete light curve would be required to resolve the true nature of this object.

The miscellaneous variables are spread throughout the CMD. Some are located on the upper main sequence, and a few seem to fall on the upper red giant branch. Several seem to occupy the same region as Cepheids. It is possible that some of these variables are Cepheids with periods in excess of 14 days.

In Figure 11 we present the $`B`$ and $`V`$-band P-L diagrams for the Cepheid variables, drawn as open circles proportional in size to their amplitudes. As expected, the amplitudes in the $`B`$-band are on average larger than in $`V`$. A clear relation between the period and magnitude is discernible. The several Cepheids too faint for their periods are probably suffering from greater than average reddening. Some of the faintest Cepheids exhibit quite large amplitudes. This is not a physical effect and is most likely caused by a similar phenomenon to the one discussed in Subsection 4.6 regarding the comparison with Paper VI photometry. If the magnitude used to convert the light curve from flux to magnitudes is measured too faint, the resulting amplitude will be too large.

For periods shorter than 8 days the relations widen upwards considerably, with the brighter Cepheids for a given period having smaller amplitudes. One possibility is that these Cepheids are pulsating in the first overtone. On the other hand these could be fundamental mode Cepheids which are heavily affected by blending: adding a constant flux would tend to increase the brightness of a Cepheid and diminish its amplitude. We will examine these stars in more depth in Section 6.

Figure 12 shows the location of the EBs and Cepheids within the field M33A. The Cepheids are plotted as open circles, proportional in size to their period and the EBs as open squares. To trace the spiral pattern of the galaxy we plot in the background all the stars with $`V<21.5`$ mag as filled dots of size proportional to their magnitude. These two types of variables appear somewhat more plentiful within the spiral arms.

## 6 First Overtone Cepheids

As we have noted in the previous Section, on the $`B`$ and $`V`$-band P-L diagrams in Fig. 11 there are Cepheids which seem too bright for their periods. In addition they possess smaller pulsational amplitudes of variability compared to the normal Cepheids. Their positions on the P-L diagrams would lead us to expect that these should be first overtone (FO) pulsators (as in Fig. 2 of Udalski et al. 1999; hereafter U99). The situation is, however, complicated by the existence of blending. As a result of blending the Cepheid should appear brighter because of the added constant flux and its amplitude, measured in magnitudes, should decrease (Mochejska et al. 2000, 2001; Stanek & Udalski 1999).

In order to try to determine whether these Cepheids are first overtone pulsators, we have checked whether they possess other properties expected of such stars. In addition to being brighter and having a smaller amplitude, Cepheids pulsating in the first overtone should have more symmetric (sinusoidal) light curves and be on average bluer than fundamental mode Cepheids (FU). The most powerful technique for discriminating them from FU Cepheids are the Fourier parameters of their light curves (Antonello & Aikawa 1995; Beaulieu et al. 1995).

To select a sample of FO Cepheid candidates we have made a division on the P-L diagrams roughly parallel to the P-L relation, between the bright low amplitude and fainter high amplitude Cepheids (dotted lines on Fig. 13). In our final sample we included the Cepheids which were above these lines in both of the P-L diagrams. There are 44 such Cepheids in our catalog (filled circles). FO Cepheids should have periods ranging from 1.7 to 6 days. We decided not to make a cutoff at higher periods, although we did not regard it likely that Cepheids with $`P>6`$ days would turn out to be FO pulsators.

On the $`V/BV`$ CMD (Fig. 14) we plot the positions of the FU Cepheids (open circles) and FO candidates (filled circles). As expected of FO Cepheids, our candidates are bluer than the FU ones. We note that on the CMD presented by U99 for the Large Magellanic Cloud (LMC) Cepheids (Fig. 4 therein) there is also some overlap between the loci of those two types of Cepheids.

The relations between the light curve Fourier parameters $`R_{21}=A_2/A_1`$, $`\varphi _{21}=\varphi _22\varphi _1`$ and the period for the FO Cepheid candidates (filled circles) and FU mode pulsators (open circles) are shown in Fig. 15. We have compared them with Fig. 3 in U99 for LMC Cepheids. We notice on the $`R_{21}/\mathrm{log}P`$ diagram for LMC that within our range of periods the FO Cepheids progress upwards in $`R_{21}`$ with increasing period, forming the second branch of the V-shaped pattern and then merge with the FU Cepheid sequence at higher periods. The situation is very similar on the $`\varphi _{21}/\mathrm{log}P`$ diagram in U99, where on the one hand the FU Cepheids are confined to a narrower sequence, but on the other there is more overlap between them and the FO pulsators.

Indeed we do observe in Fig. 15 the FO Cepheid candidates to occupy roughly the predicted positions despite our larger scatter in the Fourier parameters than in U99. From the $`\varphi _{21}/\mathrm{log}P`$ diagram in the right panel we have selected eight FO Cepheids which are most separated from the FU Cepheid sequence. We examined their positions on the $`R_{21}`$ diagram to find that they are also rather well separated from the FU mode Cepheids, especially the seven with shorter periods. These eight bona fide FO Cepheid candidates are denoted by encircled dots on Figs 13-15.

We have checked the brightest of these Cepheids, D33J013359.5+304017.0, which is unusually bright even for an overtone Cepheid, in archival HST data. It has two nearby companions, one of which appears bright enough to affect its photometry. To try to estimate the influence of blending on the Fourier parameters $`R_{21}`$ and $`\varphi _{21}`$ we have generated two new light curves for each of the eight best FO Cepheid candidates, one with an $`S_V=\frac{f_blend}{f_ceph}=0.5`$ blend added and one with the same blend subtracted. The average differences between the original and recomputed parameters were 9% for $`R_{21}`$ and 6% for $`\varphi _{21}`$. These differences are too small to change the positions of the Cepheids considerably on the $`R_{21}`$ and $`\varphi _{21}`$ diagrams. This very weak dependence of the Fourier parameters on blending strengthens our conclusion that these Cepheids are indeed first overtone pulsators.

As for the rest of the candidates, we believe many of them may also be FO Cepheids. More accurate photometry obtained with a larger instrument would be necessary to obtain better Fourier parameters of their light curves and resolve this issue.

In Figure 16 we show the light curves of the eight FO Cepheids. They appear more symmetrical in shape and do not exhibit the fast rise and slow decline typical for FU Cepheids (see Fig. 6). This fact lends further credibility to the notion that these Cepheids pulsate in the first overtone.

## 7 Conclusions

Our search for variable stars in the data from the followup observations of the detached eclipsing binary D33J013346.2+304439.9 in field M33A collected at the 2.1m KPNO telescope resulted in the discovery of 63 eclipsing binaries, 305 Cepheids, and 66 other periodic, possible long period or non-periodic variables. Out of the total 434 variables 280 are not listed in our first M33 catalog (Paper VI). Due to the short time base of our observations, we were limited to variables with periods not exceeding 14 days. Thanks to the use of a larger aperture instrument and a novel method of image reduction – the ISIS image subtraction package, we were more efficient at finding the fainter and lower amplitude variables than in our previous study of this field, especially for short period Cepheids ($`P<8d`$). We have also found a population of Cepheids which are most likely pulsating in the first overtone and for eight of them we present strong arguments in favor of this interpretation.

The method of image subtraction has two main advantages over the classical profile fitting method. It is more efficient in discovering variables: we have discovered 355 periodic variables using ISIS and only 212 with Dophot. Additionally in crowded fields image subtraction can lead to large improvements in the photometric accuracy (Alard 2000b).

During our first observing run we have found that the T2KA CCD camera produced a non-linear response to signal. We have obtained a set of calibration frames to determine a correction for the non-linearity to be applied to the science data. Although the comparison between images of differing exposure times (Fig. 3) and with the Macri et al. (2001b) stellar object catalog (Fig. 7) do not show marked non-linearity, some traces of it still remain. The presence of a residual non-linearity prevented us from obtaining an independent calibration. During the analysis of the data collected at KPNO in October 1998 with the T1KA camera we have found it to be also non-linear. We would like to caution the astronomical community of the problems we have encountered with the T1KA and T2KA cameras.

We thank the TAC of the KPNO for the generous allocation of the observing time. We would like to thank Lucas Macri for his help and for the analysis of T1KA non-linearity, Grzegorz Pojmański for $`lc`$ \- the light curve analysis utility and Wojtek Pych for his software. We thank the referee, Douglas Welch, for a prompt and useful report. BJM and JK were supported by the Polish KBN grant 2P03D003.17. BJM was also supported by the Polish KBN grant 2P03D025.19 and JK by the NSF grant AST-9819787. Support for KZS was provided by NASA through Hubble Fellowship grant HF-01124.01-A from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS5-26555. DDS acknowledges support from the Alfred P. Sloan Foundation and from NSF grant No. AST-9970812.


format: text

text:
# Radio Emission from the X-ray transient XTE J1550-564

## 1 Introduction

The soft X-ray transient XTE J1550$``$564 was discovered by the All-Sky Monitor (ASM) on board the Rossi X-ray Timing Explorer (RXTE) on MJD 51063 (1998 Sept 7; MJD=JD$``$240000.5) with an intensity of $``$0.07 Crab in the 2–12 keV range Smith (1998) and by the Burst and Transient Source Experiment (BATSE) on the Compton Gamma-Ray Observatory in the 20–100 keV range at a flux level of $`1.26\times 10^8\mathrm{erg}\mathrm{s}^1\mathrm{cm}^2`$ Wilson et al. (1998). The RXTE/ASM intensity increased steadily over the next few days, reaching $``$1.7 Crab on MJD 51071 Remillard et al. (1998). XTE J1550$``$564 flared to 6.8 Crab on MJD 51075–51076, making it the brightest X-ray nova observed with RXTE to date Remillard et al. (1998). An optical counterpart with a magnitude of V=16 was identified on MJD 51065 (Orosz, Bailyn & Jain 1998), and an orbital period of $``$1.54 days was attributed to the system Jain et al. (2000). A radio counterpart was discovered at the optical position with a flux density of 10$`\pm `$2.5 mJy on MJD 51065 Campbell-Wilson et al. (1998). Evolving QPOs Cui et al. (1999) and hard lags Wijnands et al. (1999) typical of black hole binaries were observed from XTE J1550$``$564.

In this paper we introduce briefly the radio emission from XTE J1550$``$564 at the time of the 1998 outburst. A more detailed analysis will be presented in Hannikainen et al. (2001, in preparation).

## 2 observations

XTE J1550$``$564 was observed at multiple radio frequencies during the September 1998 outburst. The Molonglo Observatory Synthesis Telescope (MOST) made twelve observations of the source at 843 MHz between MJD 51065 and 51092, while the Australia Telescope Compact Array (ATCA) made observations at 1.4, 2.3, 4.8 and 8.6 GHz between MJD 51073 and 51085. In addition, Very Long Baseline Interferometry (VLBI) images were obtained with the Australian Long Baseline Array (LBA) at 2.29 GHz. The MOST 843 MHz and ATCA 4.8 and 8.6 GHz lightcurves are shown in the upper panel of Figure 1 along with the epochs of the VLBI observations.

### 2.1 MOST & ATCA

Following the initial detection of XTE J1550$``$564 on MJD 51065, the MOST monitored the source at 843 MHz over the next 27 days (Figure 1). Between MJD 51065 and 51073, the source flux density remained between 10 and 30 mJy. After MJD 51073 the flux density began to increase, reaching 168 mJy on MJD 51076, and peaking on MJD 51078 with a flux density of 375 mJy. The flux density then declined to 120 mJy over the next three days, with a continuing decline to 14 mJy on MJD 51092, after which the monitoring ceased.

The ATCA started observing XTE J1550$``$564 on MJD 51073 while the flux density was still low – $`1012`$ mJy at 4.8 and 8.6 GHz (Figure 1). After the rise to 168 mJy detected with the MOST, the ATCA resumed observing the source and monitored it over the next ten days. The 4.8 and 8.6 GHz flux densities appear to peak approximately 0.5 days before the MOST peak, with flux densities of 283 and 234 mJy respectively and linear polarization $``$ 6%. The ATCA flux densities then followed the same decline as the MOST, reaching a level of around 20 mJy at 4.8 GHz on MJD 51086.5 when monitoring was ceased. The linearly polarized flux also evolved during the outburst (Hannikainen et al., in preparation), reminiscent of the behavior of GRO J1655$``$40 during the 1994 ejection episodes Hannikainen et al. 2000a .

### 2.2 VLBI

The VLBI observations were undertaken on MJD 51080.5 and MJD 51081.5 using the LBA Hannikainen et al. 2000b . On the first day, the 34-m DSS45, the 26-m Hobart and the 22-m Mopra antennas were used, while on the second day the 70-m DSS43 replaced the DSS45, and the six 22-m dishes of the ATCA were also used.

Figure 2 shows the 2.29 GHz VLBI images obtained on MJD 51080.5 (top) and MJD 51081.5 (bottom).

The image from MJD 51081.5 shows three distinct components, and the structure has clearly evolved compared to the previous day. Phase-referenced observations were not performed on either day so the zero coordinate in the images is arbitrary. However, changes in the spectral index, $`\alpha `$ ($`S_\nu \nu ^\alpha `$), are consistent with the appearance of a third component (Figure 1). Although the overall 2-point spectral indices, both the 843 MHz–4.8 GHz and the 4.8–8.6 GHz, stay relatively flat throughout the observing period (much flatter, say, compared to GRO J1655$``$40 which had $`\alpha `$ of $`0.5`$ to $`0.9`$ during the 1994 ejection episodes; Hannikainen et al. 2000a) the 4.8–8.6 GHz spectral index does flatten from $`0.45`$ to $`0.1`$ just before the VLBI observations. Thus, the central component in Figure 2 (bottom) is tentatively identified with a new optically thick outburst. If this interpretation is correct, then the apparent separation velocity of the outermost components between MJD 51080.5 and 51081.5 is $`v>2c`$ Hannikainen et al. 2000b , based on a distance estimate of 3.5–5 kpc derived from H i observations.

## 3 Multiwavelength Lightcurves

Figure 3 shows the multiwavelength behavior of XTE J1550$``$564. The significance of the sequence of lightcurve features in each waveband is discussed in Wu et al. (2000) and an in-depth analysis of the multiwavelength data is in Hannikainen et al. (2001, in preparation), so here we shall only summarize the more salient features. The outburst sequence is characterized by an impulsive rise in hard X-rays accompanied by a gradual brightening in the soft X-rays, followed by a period of exponential-like decay in the hard X-ray luminosity, and then a giant flare in all wavebands. The hard and soft X-rays peaked simultaneously, while the radio flare maxima occurred about 1.3–1.8 days later. As the radio emission continued to decline to a low level after the flare, the X-ray intensities settled into a plateau, lingered for another 35 days, and then faded away.

In order to explain the first impulsive BATSE hard X-ray outburst, which does not conform to the standard disk-instability scenario, a model involving a magnetic secondary star is proposed Wu et al. (2000). In this model, coronal magnetic activity lifts matter above the stellar surface of the secondary (which is presumed to be a G/K subgiant; see also Bildsten & Rutledge 1999), and if this happens near the L1 point of the Roche lobe, this matter will fall easily into the Roche lobe of the black hole, as gravity is weak near the L1 point. Provided the magnetic stress/tension is stronger than the gas pressure terms, matter can be trapped there and accumulate over time. Eventually, this magnetic dam may burst leading to field-line breaking and a matter avalanche. The initial angular momentum of the matter with respect to the black hole is insignificant, as the matter is practically in free-fall from zero velocity. This flow will therefore be quasi-spherical, and deviation occurs only when approaching the black hole. An accretion shock forms when the infalling material encounters the centrifugal barrier of the hole, converting a portion of the kinetic energy to thermal energy, giving rise to radiation. These photons are then upscattered by the surging infalling matter to hard X-ray energies, which is observed as the first hard X-ray outburst. As the quasi-spherical flow subsides, the hard X-ray luminosity declines, and at the same time, residual matter with substantial angular momentum begins to condense in the equatorial plane to form an accretion disk – this is seen as the gradual rise in the soft X-rays (see Igumenschev, Illarionov & Abramowicz 1999). Irradiation of the secondary as a consequence of the first hard X-ray outburst will result in an increase of mass transfer, and provoke instabilities, seen in the large fluctuations in the soft X-ray luminosity between MJD 51070 and 51075. As a result, the accretion disk becomes unstable and portions of it collapse, leading to mass ejection, as evidenced by the radio flare and the discrete components observed with the VLBI.

## 4 Summary

The XTE J1550$``$564 hard and soft X-ray flare in 1998 September was accompanied by radio emission. The radio outburst, which reached 375 mJy at 843 MHz, lagged the X-ray flare by 1.8 days. There was spectral evolution during the outburst, and polarization was detected at the 6% level. The VLBI images showed evolving structure confirming that there was an outflow from the system. The inferred superluminal separation velocity of the two outermost components, $`v>2c`$, strengthens the argument for XTE J1550$``$564 being a black-hole binary system similar to GRO J1655$``$40 and GRS 1915+105.

## Acknowledgments

We thank Michael McCollough and the BATSE team for the BATSE data, and the referee Chris Shrader for helpful suggestions. DCH acknowledges the support of a PPARC postdoctoral research grant to the University of Southampton, and financial support from the Academy of Finland, and thanks the Astrophysics Department, University of Sydney for hospitality during her visits. KW thanks Prof. P. Charles for funding his visits to the University of Southampton. MOST is operated by the University of Sydney and supported by grants from the Australian Research Council. The Australia Telescope Compact Array is funded by the Commonwealth of Australia for operation as a National Facility managed by CSIRO. This research has made use of data obtained from the High Energy Astrophysics Science Archive Research Center (HEASARC) provided by NASA’s Goddard Space Flight Center.


text:
# 1 Introduction

## 1 Introduction

Non-forward parton distributions  have recently been proposed as tools to study hard exclusive reactions: deeply virtual Compton scattering and exclusive meson production  in $`\gamma ^{}p`$ reactions at large $`Q^2`$ (the incident photon virtuality), large $`s`$ (centre of mass energy) and small $`t`$ (squared momentum transfer). It has further been proposed that deeply virtual Compton scattering (DVCS) could provide information about quark and gluonic orbital angular-momentum in the proton . It is now well established that the spin structure of the proton is sensitive to the physics of the axial anomaly and dynamical $`U_A(1)`$ symmetry breaking in QCD – for reviews see . In this paper we discuss the role of these phenomena in non-forward parton distributions and in deeply virtual Compton scattering. The important details are the gauge dependence of non-forward matrix elements of the anomalous Chern Simons current $`K_+`$ – even in perturbation theory – and the role of zero modes in dynamical $`U_A(1)`$ symmetry breaking.

Perturbative QCD factorisation applies to DVCS. In perturbative QCD the amplitudes for deeply virtual Compton scattering and hard exclusive meson production at small $`t`$ may be written as the integral over the product of non-forward parton distributions and coefficient functions \[1-5,11\]. The non-forward parton distributions describe the momentum and spin distribution of partons in the target proton and the coefficient functions describe the infra-red safe hard $`\gamma ^{}`$–parton interaction, which may be calculated in perturbative QCD. Working in light-cone gauge $`A_+=0`$ the spin independent and spin dependent non-forward parton distributions are written as the Fourier transform of light-cone correlation functions: the non-forward matrix elements of “$`\overline{q}(0)\gamma _+q(y)`$” and “$`\overline{q}(0)\gamma _+\gamma _5q(y)`$” respectively at $`y_+=y_{}=0`$.

The unpolarised DVCS cross-section receives contributions from both the spin-dependent and spin-independent parts of the DVCS amplitude $`𝒜_{\mathrm{DVCS}}`$ . Experiments with polarised beams and targets may help to resolve the separate spin-dependent and spin-independent contributions to $`𝒜_{\mathrm{DVCS}}`$ . Here we assume that both contributions can be extracted from future experimental data and investigate the role of the axial anomaly in the interpretation of the data.

In QCD the bare light-cone correlation functions are ultra-violet divergent requiring renormalisation . This renormalisation induces the anomalous dimensions ($`Q^2`$ dependence of the distributions). It also means that the spin dependent (non-forward) parton distribution is sensitive to the axial anomaly . It is important to check whether any non-perturbative physics induced by the anomaly is consistent with the original perturbative QCD factorisation.

The main results of this paper are the following. First, the flavour-singlet part of $`𝒜_{\mathrm{DVCS}}^{\mathrm{spin}}`$ is sensitive to the QCD axial anomaly and does not have a simple canonical interpretation. In full QCD (beyond perturbation theory) whether factorisation provides a complete description of this contribution is sensitive to the dynamics of axial U(1) symmetry breaking – the famous $`U_A(1)`$ problem of the $`\eta ^{}`$ mass – which has the potential to induce zero-mode contributions to the flavour-singlet axial-charge $`g_A^{(0)}`$ and its non-forward generalisation. Second, Ji  has shown that information about quark and gluon total angular-momentum in the proton ($`J_z^q`$ and $`J_z^g`$) may be extracted from the forward limit of the spin-independent part of $`𝒜_{\mathrm{DVCS}}`$. Through the axial anomaly in QCD, some fraction of the intrinsic spin content of the nucleon and of the constituent quark is carried by its quark and gluon partons and some fraction is carried by gluon topology through a zero-mode . The zero-mode term $`𝒞`$ is a model independent physical quantity which can, in principle, be measured through $`\nu p`$ elastic scattering . The axial anomaly and zero-modes contribute separately to $`S_z^q`$ and $`L_z^q`$ but cancels in the sum. Possible zero-modes and scheme dependence mean that values of $`L_z`$ extracted from future DVCS data have to be labelled with respect to the degrees of freedom and the perturbative QCD scheme dependence used to define $`S_z^q`$.

The plan of the paper is as follows. To motivate our discussion of spin in non-forward parton distributions we begin in Section 2 with a brief overview of our present understanding of the proton spin problem from polarised deep inelastic scattering. In Section 3 we give a brief overview of non-forward parton distributions – for excellent reviews see . Section 4 discusses the axial anomaly and Section 5 its role in the spin dependent part of $`𝒜_{\mathrm{DVCS}}`$. In Section 6 we discuss the role of the anomaly and gluon topology in the interpretation of information about $`L_z`$ from DVCS. Finally, in Section 7 we make our conclusions.

## 2 The spin structure of the proton

To motivate our discussion of spin effects in DVCS we first briefly review what is known about the quark and gluonic intrinsic spin structure of the proton from the interpretation of polarised deep inelastic scattering data. We start with the simple sum-rule for the spin $`+\frac{1}{2}`$ proton

$$\frac{1}{2}=\frac{1}{2}\mathrm{\Sigma }+L_z+S_z^g.$$
(1)
Here, $`\frac{1}{2}\mathrm{\Sigma }`$ and $`S_z^g`$ are the quark and gluonic intrinsic spin contributions to the nucleon’s spin and $`L_z`$ is the orbital contribution. One would like to understand the spin decomposition, Eq.(1), both in terms of the fundamental QCD quarks and gluons and also in terms of the constituent quark quasi-particles of low-energy QCD. In relativistic constituent quark models $`\mathrm{\Sigma }`$ is given by the flavour-singlet axial-charge $`g_A^{(0)}`$. The value of $`g_A^{(0)}`$ extracted from polarised deep inelastic scattering experiments is $`g_A^{(0)}|_{\mathrm{pDIS}}=0.2`$ – 0.35 , roughly half of the value predicted by the constituent quark models. In QCD the interpretation of the individual quantities in Eq.(1) is quite subtle because of the axial anomaly , issues of gauge invariance  and dynamical $`U_A(1)`$ symmetry breaking . In QCD the axial anomaly induces various gluonic contributions to $`g_A^{(0)}`$. Working in light-cone gauge $`A_+=0`$ one finds

$$g_A^{(0)}=\left(\underset{q}{}\mathrm{\Delta }q3\frac{\alpha _s}{2\pi }\mathrm{\Delta }g\right)_{\mathrm{partons}}+𝒞$$
(2)
Here $`\frac{1}{2}\mathrm{\Delta }q`$ and $`\mathrm{\Delta }g`$ are the amount of spin carried by quark and gluon partons in the polarised proton and $`𝒞`$ measures the gluon-topological contribution to $`g_A^{(0)}`$ .

In Eq.(2) $`\mathrm{\Delta }q_{\mathrm{partons}}`$ is associated with the hard photon scattering on quark and antiquarks with low transverse momentum squared, $`k_T^2`$ of the order of typical gluon virtualities in the proton, and $`\mathrm{\Delta }g_{\mathrm{partons}}`$ is associated with the hard photon scattering on quarks and antiquarks carrying $`k_T^2Q^2`$. Jet studies and semi-inclusive measurements of the nucleon’s spin-dependent $`g_1`$ structure function may be used to determine $`\mathrm{\Delta }g`$ and $`\mathrm{\Delta }q`$ for each flavour . The topology term $`𝒞`$ is associated with dynamical axial U(1) symmetry breaking and has support only at Bjorken $`x=0`$. It is missed by polarised deep inelastic scattering experiments but could, in principle, be measured in elastic $`\nu p`$ scattering . An example how to obtain a finite value for $`𝒞`$ is provided by Crewther’s theory of quark-instanton interactions . There, any instanton induced suppression of $`g_A^{(0)}|_{\mathrm{pDIS}}`$ (the axial charge carried by partons with finite momentum fraction $`x>0`$) is compensated by a shift of axial charge to the zero-mode so that the total axial-charge $`g_A^{(0)}`$ including $`𝒞`$ is conserved.

The quark and gluonic parton decomposition of $`g_A^{(0)}`$ in Eq.(2) is factorisation scheme dependent. Eq.(2) describes the singlet axial charge in the $`k_T^2`$ cut-off parton model and in the AB and JET schemes . In the $`\overline{\mathrm{MS}}`$ scheme the polarised gluon contribution is absorbed into $`\mathrm{\Delta }q`$ so that $`\mathrm{\Delta }q_{\overline{\mathrm{MS}}}=\left(\mathrm{\Delta }q\frac{\alpha _s}{2\pi }\mathrm{\Delta }g\right)_{\mathrm{partons}}`$ .

The QCD $`g_A^{(0)}`$ is measured by the proton forward matrix element of the flavour singlet axial-vector current.

$$2ms_\mu g_A^{(0)}=P,S|J_{\mu 5}^{GI}|P,S_c$$
(3)
where $`|P,S`$ denotes a proton state with momentum $`P`$ and spin $`S`$ and

$$J_{\mu 5}^{GI}=\left[\overline{u}\gamma _\mu \gamma _5u+\overline{d}\gamma _\mu \gamma _5d+\overline{s}\gamma _\mu \gamma _5s\right]_{GI}$$
(4)
is the gauge-invariantly renormalised singlet axial-vector operator. In each of the “partons”, AB and JET schemes $`\mathrm{\Delta }q`$ and $`\mathrm{\Delta }g`$ are measured by the forward matrix elements of the partially conserved axial-vector current and gluonic Chern-Simons current contributions to $`J_{+5}^{GI}`$ between parton states in the light-cone gauge $`A_+=0`$ – see Sections 4 and 5 below. In light-cone gauge the forward matrix elements of the plus components of these gauge dependent currents are invariant under residual gauge degrees of freedom. However, the non-forward matrix elements of these currents are not invariant, even in perturbative QCD. This means that some of the schemes commonly used to analyse polarised deep inelastic data are not gauge invariant when “skewed” or extrapolated away from the forward direction.

In summary, the axial anomaly brings in gauge invariance issues and the possibility of zero-momentum contributions to the individual spin contributions in Eq.(1). How does this physics enter the theory of deeply virtual Compton scattering ?

## 3 Non-forward parton distributions

Non-forward parton distributions are defined as the Fourier transforms of light-cone correlation functions. Consider an incident proton with mass $`M`$ and momentum $`P_\mu `$ which gets given momentum $`\mathrm{\Delta }_\mu =(P^{}P)_\mu `$ and emerges into the final state with momentum $`P_\mu ^{}`$. We use $`\overline{P}_\mu =\frac{1}{2}(P^{}+P)_\mu `$ to denote the average nucleon momentum. In this paper we follow the notation of Radyushkin . We let $`x`$ denote the fraction of light-cone momentum $`k_+=xP_+`$ carried by a parton in the incident proton and $`\zeta =\mathrm{\Delta }_+/P_+`$ denote the amount of light-cone momentum transfered to the proton.

In perturbative QCD the spin averaged cross-section for deeply virtual Compton scattering receives contributions from both spin independent and spin dependent non-forward parton distributions. The formula suggested by perturbative QCD for the deeply virtual Compton amplitude $`M^{IJ}=𝒜_{\mathrm{DVCS}}`$ is:

$`M^{IJ}(\stackrel{}{q}_{},\stackrel{}{\mathrm{\Delta }}_{},\zeta )=`$ (5)
$`\mathrm{e}_{\mathrm{q}}^{}{}_{}{}^{2}{\displaystyle \frac{1}{2\overline{P}^+}}{\displaystyle _{\zeta 1}^{+1}}dx[F_q(x,\zeta ,t,\mu ^2)C_q(x,\zeta ,Q^2,\mu ^2)+F_g(x,\zeta ,t,\mu ^2)C_g(x,\zeta ,Q^2,\mu ^2)`$
$`+\stackrel{~}{F}_q(x,\zeta ,t,\mu ^2)\stackrel{~}{C}_q(x,\zeta ,Q^2,\mu ^2)+\stackrel{~}{F}_g(x,\zeta ,t,\mu ^2)\stackrel{~}{C}_g(x,\zeta ,Q^2,\mu ^2)]+𝒪({\displaystyle \frac{1}{Q}})`$
Here $`(I,J)`$ refer to the polarisation vectors ($``$ or $``$) of the initial and final state photons – here both taken as $``$; $`F_q`$ and $`F_g`$ are the spin-independent non-forward quark and gluonic parton distributions and $`\stackrel{~}{F}_q`$ and $`\stackrel{~}{F}_g`$ are the spin-dependent parton distributions. These distributions are defined in Eqs.(7-10) below; $`\mu ^2`$ denotes the renormalisation scale. In Eq.(5)

$`C_q^{}`$ $`=`$ $`C_q^{}=\left({\displaystyle \frac{1}{xiϵ}}+{\displaystyle \frac{1}{x\zeta +iϵ}}\right)+O(\alpha _s)`$ (6)
$`C_q^{}`$ $`=`$ $`C_q^{}=\left({\displaystyle \frac{1}{xiϵ}}{\displaystyle \frac{1}{x\zeta +iϵ}}\right)+O(\alpha _s)`$
are the non-forward quark coefficient functions with $`C_q^{}=C_q^{}=\stackrel{~}{C}_q^{}=\stackrel{~}{C}_q^{}=0`$. The gluonic coefficients start at order $`\alpha _s`$, viz. $`C_g,\stackrel{~}{C}_gO(\alpha _s)`$. In DVCS $`\zeta `$ plays the role of the Bjorken variable $`\zeta =x_{\mathrm{Bj}}`$. Radyushkin  has derived conditions on the non-forward distributions which must be satisfied in order that QCD factorisation holds beyond perturbation theory. One of these conditions is that the non-forward parton distributions do not contain any singular behaviour at $`x=0`$ or $`x=\zeta `$. We shall return to this point in Section 5 below.

The spin independent non-forward quark and gluon distributions are:

$`F_q(x,\zeta ,t)`$ $`=`$ $`{\displaystyle \frac{dy_{}}{8\pi }e^{ixP_+y_{}/2}P^{}|\overline{q}(0)\gamma _+𝒫e^{i_y_{}^0𝑑y_{}A_+}q(y)|P}|_{y_+=y_{}=0}`$
$`=`$ $`{\displaystyle \frac{1}{2\overline{P}_+}}\overline{U}(P^{})\left[H_q(x,\zeta ,t)\gamma _++E_q(x,\zeta ,t){\displaystyle \frac{i\sigma _{+\alpha }(\mathrm{\Delta }^\alpha )}{2M}}\right]U(P)`$
and

$`F_g(x,\zeta ,t)`$ $`=`$ $`{\displaystyle \frac{1}{xP_+}}{\displaystyle \frac{dy_{}}{8\pi }e^{ixP_+y_{}/2}P^{}|G_{+\alpha }(0)𝒫e^{i_y_{}^0𝑑y_{}A_+}G_+^\alpha (y)|P}|_{y_+=y_{}=0}`$
$`=`$ $`{\displaystyle \frac{1}{2\overline{P}_+^2}}\overline{U}(P^{})\left[H_g(x,\zeta ,t)\gamma _++E_g(x,\zeta ,t){\displaystyle \frac{i\sigma _{+\alpha }(\mathrm{\Delta }^\alpha )}{2M}}\right]U(P)`$
respectively. In the forward limit $`\zeta =0`$ one recovers the spin independent quark and gluon distributions measured in deep inelastic scattering: $`(q\pm \overline{q})(x)=\frac{1}{2}(H_q(x,0,0)H_q(x,0,0))`$ and $`g(x)=\frac{1}{2}(H_g(x,0,0)H_g(x,0,0))`$ respectively. The spin dependent distributions are:

$`\stackrel{~}{F}_q(x,\zeta ,t)`$ $`=`$ $`{\displaystyle \frac{dy_{}}{8\pi }e^{ixP_+y_{}/2}P^{}|\overline{q}(0)\gamma _+\gamma _5𝒫e^{i_y_{}^0𝑑y_{}A_+}q(y)|P}|_{y_+=y_{}=0}`$
$`=`$ $`{\displaystyle \frac{1}{2\overline{P}_+}}\overline{U}(P^{})\left[\stackrel{~}{H}_q(x,\zeta ,t)\gamma _+\gamma _5+\stackrel{~}{E}_q(x,\zeta ,t){\displaystyle \frac{1}{2M}}\gamma _5(\mathrm{\Delta }_+)\right]U(P)`$
and

$`\stackrel{~}{F}_g(x,\zeta ,t)`$ $`=`$ $`{\displaystyle \frac{i}{xP_+}}{\displaystyle \frac{dy_{}}{8\pi }e^{ixP_+y_{}/2}P^{}|G_{+\alpha }(0)𝒫e^{i_y_{}^0𝑑y_{}A_+}\stackrel{~}{G}_+^\alpha (y_{})|P}|_{y_+=y_{}=0}`$
$`=`$ $`{\displaystyle \frac{1}{2\overline{P}_+^2}}\overline{U}(P^{})\left[\stackrel{~}{H}_g(x,\zeta ,t)\gamma _+\gamma _5+\stackrel{~}{E}_g(x,\zeta ,t){\displaystyle \frac{1}{2M}}\gamma _5(\mathrm{\Delta }_+)\right]U(P)`$
In the forward limit $`\zeta =0`$ one recovers the spin dependent quark and gluon distributions measured in deep inelastic scattering: $`\mathrm{\Delta }(q\pm \overline{q})(x)=\frac{1}{2}(\stackrel{~}{H}_q(x,0,0)\pm \stackrel{~}{H}_q(x,0,0))`$ and $`\mathrm{\Delta }g(x)=\frac{1}{2}(\stackrel{~}{H}_g(x,0,0)+\stackrel{~}{H}_g(x,0,0))`$ respectively. The isotriplet combination $`(\stackrel{~}{E}_u\stackrel{~}{E}_d)`$ contains the pion pole and the flavour-singlet combination $`(\stackrel{~}{E}_u+\stackrel{~}{E}_d+\stackrel{~}{E}_s)`$ is sensitive to the axial U(1) problem  through the $`\eta ^{}`$ pole in the pseudoscalar form-factor.

We now focus on the spin-dependent non-forward quark distribution $`\stackrel{~}{F}_q`$ to discuss the construction and renormalisation of these distributions. The parton model and the light-cone correlation functions are normally formulated in the light-cone gauge $`A_+=0`$. Here the path-ordered exponential becomes a trivial unity factor and is dropped from the formalism:

$$P^{}|\overline{q}(0)\gamma _+\gamma _5𝒫e^{i_y_{}^0𝑑y_{}A_+}q(y)|PP^{}|\overline{q}(0)\gamma _+\gamma _5q(y)|P$$
(11)
In semi-classical QCD, before we come to discuss renormalisation and anomalies, the gluonic degrees of freedom have dropped out along with the path-ordered exponential – hence the terminology “quark distribution”. Before we consider subtleties associated with anomaly theory – see Section 4 below – the non-forward parton distributions have the following simple interpretation . Expanding out the quark and gluon field operators in $`q(y)`$ and $`G_{\mu \nu }(y)`$ one finds the following. For $`\zeta <x<1`$ we have the situation where one removes a quark carrying light-cone momentum $`k_+=xP_+`$ and transverse momentum $`\stackrel{}{k}_{}`$ from the initial state proton and re-inserts it into the final state proton with the same chirality but with light-cone momentum fraction $`x\zeta `$ and transverse momentum $`\stackrel{}{k}_{}\stackrel{}{\mathrm{\Delta }}_{}`$. For $`\zeta 1<x<0`$ one finds the situation for removing an antiquark with momentum fraction $`\zeta x`$ and re-inserting it with momentum fraction $`x`$. For $`0<x<\zeta `$ the photons scatter off a virtual quark-antiquark pair in the initial proton wavefunction: the quark of the pair has light-cone momentum fraction $`x`$ and transverse momentum $`\stackrel{}{k}_{}`$, whereas the antiquark has light-cone momentum fraction $`\zeta x`$ and transverse momentum $`\stackrel{}{\mathrm{\Delta }}_{}\stackrel{}{k}_{}`$. The third region $`0<x<\zeta `$ is not present in deep inelastic scattering where $`\zeta =0`$. The points $`x=0`$ and $`x=\zeta `$ correspond to zero-momentum modes. The flavour-singlet spin-dependent parton distributions at these two points are sensitive to the role of zero-modes in dynamical $`U_A(1)`$ symmetry breaking , which have the potential to generate new non-perturbative contributions to the spin-dependent part of $`𝒜_{DVCS}`$ beyond those appearing in the factorisation formula (5). The $`J=0`$ Regge fixed pole  in the spin-independent part of $`𝒜_{DVCS}`$ which generates a constant real term in $`𝒜_{\mathrm{DVCS}}`$ is manifest in Eqs.(5,6) through the $`\zeta 0`$ limit of $`C_q`$.

The $`x`$-moments, $`_{\zeta 1}^{+1}𝑑xx^n`$, of the non-forward parton distributions are evaluated as follows. First one writes $`x^n`$ as a derivative (in $`y_{}`$) acting on $`e^{ixP_+y_{}/2}`$. Integrating by parts (with respect to $`x`$) over the exponential yields a Dirac delta-function $`\delta (y_{})`$. The $`y_{}`$ integral then projects out the non-forward matrix elements of local operators. One finds

$`{\displaystyle _{\zeta 1}^1}𝑑xx^n\stackrel{~}{F}_q(x,\zeta ,t)`$ (12)
$`={\displaystyle \frac{1}{2}}\left({\displaystyle \frac{2}{P_+}}\right)^nP^{}|\overline{q}(0)\gamma _+\gamma _5(i_+)^nq(0)|P|_{y_+=y_{}=0}`$
$`={\displaystyle \frac{1}{2\overline{P}_+}}\overline{U}(P^{})\left[{\displaystyle _{\zeta 1}^{+1}}𝑑xx^n\stackrel{~}{H}_q(x,\zeta ,t)\gamma _+\gamma _5+{\displaystyle _{\zeta 1}^{+1}}𝑑xx^n\stackrel{~}{E}_q(x,\zeta ,t){\displaystyle \frac{1}{2M}}\gamma _5(\mathrm{\Delta }_+)\right]U(P)`$
A caveat is due here: renormalisation means that we have to be careful not to simply set $`y_{}=0`$ in the point-split operator to obtain the local operator. In QCD the bare light-cone correlation functions are ultra-violet divergent requiring renormalisation . This renormalisation means that spin independent (non-forward) parton distribution $`\stackrel{~}{F}_q`$ is sensitive to the axial anomaly.

We now discuss the effect of the anomaly in non-forward parton distributions.

## 4 Renormalisation and anomalies

Going beyond the semi-classical approximation to QCD we have to take into account that the renormalised composite operator

$$\left[\overline{q}\gamma _\mu \gamma _5q\right]_{\mu ^2}^R(0)\overline{q}(0)\gamma _\mu \gamma _5\mathrm{multiplied}\mathrm{by}q(0)$$
(13)
This point is especially important when we evaluate the integral over $`y_{}`$ of the point-split matrix element with $`\delta (y_{})`$. Evaluating the moments of $`\stackrel{~}{F}_q`$ is non-trivial because the point-split operator is highly singular in the limit that the point splitting is taken to zero. In full QCD it is necessary to work with renormalised composite operators instead of their semi-classical prototypes. To see this explicitly consider Schwinger’s derivation  of the axial anomaly using point split regularisation. For $`(z^{}z^{\prime \prime })0`$, one finds that the vacuum to vacuum matrix element of the point-split operator in an external gluon field is:

$$\mathrm{vac}|\left[\overline{q}(z^{})\gamma _+\gamma _5q(z^{\prime \prime })\right]|\mathrm{vac}\frac{ig}{8\pi ^2}\stackrel{~}{G}_{+\nu }\frac{(z^{}z^{\prime \prime })^\nu }{(z^{}z^{\prime \prime })^2}$$
(14)
where the gluonic term arises from pinching a gluonic insertion between $`z^{}`$ and $`z^{\prime \prime }`$. Going to the light-cone ($`z_T=0`$, $`z_+0`$) the factor

$$\frac{(z^{}z^{\prime \prime })_+}{(z^{}z^{\prime \prime })^2}\frac{1}{(z^{}z^{\prime \prime })_{}}$$
(15)
which diverges when $`z^{}z^{\prime \prime }`$. One clearly has to be careful and not ensure that the theory and its interpretation does not assume equality of both sides in Eq.(13).

The problem is resolved  if, working in light-cone gauge, we define

$$P^{},S^{}|\left[\overline{q}(0)\gamma _+\gamma _5q(y_{})\right]_{\mu ^2}^R|P,S_c\underset{n}{}\frac{(iy_{})^n}{n!}P^{},S^{}|\left[\overline{q}\gamma _+\gamma _5(iD_+)^nq\right]_{\mu ^2}^R(0)|P,S_c$$
(16)
and

$$P^{},S^{}|\mathrm{Tr}\left[G_{+\nu }(0)\stackrel{~}{G}_+^\nu (y_{})\right]_{\mu ^2}^R|P,S_c\underset{n}{}\frac{(iy_{})^n}{n!}P^{},S^{}|\mathrm{Tr}\left[G_{+\nu }(iD_+)^n\stackrel{~}{G}_+^\nu \right]_{\mu ^2}^R(0)|P,S_c$$
(17)
That is, we treat the non-local operators in Eqs.(7-10) as a series expansion in terms of renormalised composite local operators in the operator product expansion. The superscript $`R`$ denotes the renormalisation prescription and $`\mu ^2`$ denotes the renormalisation scale. The composite operators $`\left[\overline{q}\gamma _+\gamma _5(iD_+)^{2n}q\right](0)`$ and $`\mathrm{Tr}\left[G_{+\nu }(iD_+)^{2n}\stackrel{~}{G}_+^\nu \right](0)`$ mix under renormalisation. The axial vector and higher-spin axial tensor operators are sensitive to the axial anomaly.

### 4.1 The axial anomaly

The gauge-invariantly renormalised flavour singlet axial-vector current

$$J_{\mu 5}^{GI}=\left[\overline{u}\gamma _\mu \gamma _5u+\overline{d}\gamma _\mu \gamma _5d+\overline{s}\gamma _\mu \gamma _5s\right]_{\mu ^2}^{GI}$$
(18)
satisfies the anomalous divergence equation

$$^\mu J_{\mu 5}^{GI}=2f^\mu K_\mu +\underset{i=1}{\overset{f}{}}2im_i\overline{q}_i\gamma _5q_i$$
(19)
Here

$$K_\mu =\frac{g^2}{16\pi ^2}ϵ_{\mu \nu \rho \sigma }\left[A_a^\nu \left(^\rho A_a^\sigma \frac{1}{3}gf_{abc}A_b^\rho A_c^\sigma \right)\right]$$
(20)
is a renormalised version of the gluonic Chern-Simons current and the number of light flavours $`f`$ is $`3`$. Eq.(20) allows us to write

$$J_{\mu 5}^{GI}=J_{\mu 5}^{\mathrm{con}}+2fK_\mu $$
(21)
where $`J_{\mu 5}^{\mathrm{con}}`$ and $`K_\mu `$ satisfy the divergence equations

$$^\mu J_{\mu 5}^{\mathrm{con}}=\underset{i=1}{\overset{f}{}}2im_i\overline{q}_i\gamma _5q_i$$
(22)
and

$$^\mu K_\mu =\frac{g^2}{8\pi ^2}G_{\mu \nu }\stackrel{~}{G}^{\mu \nu }.$$
(23)
Here $`\frac{g^2}{8\pi ^2}G_{\mu \nu }\stackrel{~}{G}^{\mu \nu }`$ is the topological charge density. The partially conserved current is scale invariant and the scale dependence of $`J_{\mu 5}^{GI}`$ is carried entirely by $`K_\mu `$. When we make a gauge transformation $`U`$ the gluon field transforms as

$$A_\mu UA_\mu U^1+\frac{i}{g}(_\mu U)U^1$$
(24)
and the operator $`K_\mu `$ transforms as

$$K_\mu K_\mu +i\frac{g}{16\pi ^2}ϵ_{\mu \nu \alpha \beta }^\nu \left(U^{}^\alpha UA^\beta \right)+\frac{1}{96\pi ^2}ϵ_{\mu \nu \alpha \beta }\left[(U^{}^\nu U)(U^{}^\alpha U)(U^{}^\beta U)\right].$$
(25)
Gauge transformations shuffle a scale invariant operator quantity between the two operators $`J_{\mu 5}^{\mathrm{con}}`$ and $`K_\mu `$ whilst keeping $`J_{\mu 5}^{GI}`$ invariant.

The non-abelian three-gluon part of $`K_+`$ vanishes in $`A_+=0`$ gauge and the forward matrix elements of $`K_+`$ are invariant under residual gauge degrees of freedom in this gauge. Furthermore the forward matrix elements of $`K_+`$ measure the amount of spin carried by gluonic partons in the target . This leads ultimately to the “partons”, AB and JET scheme decompositions of $`g_A^{(0)}`$ in Section 2. As soon as we go away from the forward direction matrix elements of $`K_+`$ will pick up a gauge-dependent contribution, which must decouple from any physical observable.

### 4.2 The axial anomaly and higher moments

The anomaly is also present in the $`C=+1`$ higher spin axial tensors . In general, for a given choice of renormalisation prescription $`R`$, the renormalised axial tensor operator differs from the gauge invariant operator by a multiple of a gauge-dependent, gluonic counterterm $`K_{\mu \mu _1\mathrm{}\mu _{2n}}`$, viz.

$`\left[\overline{q}\gamma _\mu \gamma _5iD_{\mu _1}\mathrm{}iD_{\mu _{2n}}q\right]_{Q^2}^R(0)`$ $`=`$ (26)
$`\left[\overline{q}\gamma _\mu \gamma _5iD_{\mu _1}\mathrm{}iD_{\mu _{2n}}q\right]_{Q^2}^{GI}(0)`$ $`+`$ $`\lambda _{R,n}^{(K)}\left[K_{\mu \mu _1\mathrm{}\mu _{2n}}\right]_{Q^2}(0)+\lambda _{R,n}^{(G)}\left[G_{\mu \alpha }iD_{\mu _1}\mathrm{}\stackrel{~}{G}_{\mu _{2n}}^\alpha \right]_{Q^2}^{(GI)}(0)`$
Shifting between possible renormalisation schemes, we pick up contributions from counterterms involving the gauge-invariant gluonic operators $`G_{\mu \alpha }iD_{\mu _1}\mathrm{}iD_{\mu _{2n1}}\stackrel{~}{G}_{\mu _{2n}}^\alpha `$ and also the gauge-dependent $`K_{\mu \mu _1\mathrm{}\mu _{2n}}`$. The coefficients $`\lambda _{R,n}^{(K)}`$ and $`\lambda _{R,n}^{(G)}`$ are fixed by the choice of renormalisation prescription. In $`A_+=0`$ gauge the two-gluon part of $`K_{\mu \mu _1\mathrm{}\mu _{2n}}`$ reads

$$K_{+\mathrm{}+(2n+1)}=\frac{\alpha _s}{\pi }\lambda _n^{(K)}ϵ_{+\lambda \alpha \beta }A^\alpha ^\lambda (i_+)^{2n}A^\beta $$
(27)
It is not easy to derive the $`K_{\mu \mu _1\mathrm{}\mu _{2n}}`$ beyond the two gluon term. Unlike the topological charge density $`\frac{\alpha _s}{4\pi }G_{\mu \nu }\stackrel{~}{G}^{\mu \nu }`$ the higher spin operators $`\frac{\alpha _s}{4\pi }G_{\mu \nu }iD_{\mu _1}\mathrm{}iD_{\mu _{2n}}\stackrel{~}{G}^{\mu \nu }`$ are not topological invariants for $`n1`$. It follows that they are not total derivatives and, therefore, one cannot use a divergence equation alone to fix the non-abelian part of $`K_{\mu \mu _1\mathrm{}\mu _{2n}}`$. There is no equation $`^{\mu _j}𝒮𝒦_{\mu \mu _1\mathrm{}\mu _{2n}}=\frac{\alpha _s}{4\pi }G_{\mu \nu }iD_{\mu _1}\mathrm{}iD_{\mu _{2n}}\stackrel{~}{G}^{\mu \nu }`$ for $`n1`$.

## 5 The axial anomaly in $`\stackrel{~}{F}_q`$

The first observation to make is that the non-forward matrix elements of $`K_+`$ (and $`K_{+\mathrm{}+(2\mathrm{n}+1)}`$) are gauge dependent in $`A_+=0`$ gauge, even in perturbation theory. This means that these operators must decouple from any factorisation scheme or “generalised operator product expansion”  for hard exclusive reactions like deeply virtual Compton scattering. This is in contrast to the situation in deep inelastic scattering where the forward matrix elements of $`K_+`$ are invariant under residual gauge degrees of freedom in $`A_+=0`$ gauge.

Since there is no gauge-invariant local gluonic operator with quantum numbers $`J^{PC}=1^{++}`$ it follows that the first moment of the spin-dependent non-forward gluonic coefficient must vanish. This holds true in the non-forward generalisation of the $`\overline{\mathrm{MS}}`$ factorisation scheme for handling the infra-red mass singularities. <sup>1</sup><sup>1</sup>1 See eg. the calculations of Mankiewicz et al.  and Ji and Osborne ; in the notation of Mankiewicz et al.  $`\frac{}{u}C_1^{A,g}|_{u=0}=0`$. However, there is no gauge-invariant non-forward generalisation of the popular JET and AB schemes used to describe polarised deep inelastic data. In these schemes there is an explicit gluonic contribution to the first moment of $`g_1`$ associated with the invariant contribution of $`K_+`$ (in $`A_+=0`$ gauge) to the forward matrix element of $`J_{+5}^{GI}`$ induced by a non-vanishing first moment of the spin-dependent gluonic coefficient.

The gauge dependence of the non-forward matrix elements of $`K_+`$ also means that one has to be careful with the interpretation of the spin-dependent non-forward parton distributions in terms of amplitudes to extract and then re-insert a (well defined) parton with a given momentum. Consider the non-forward spin-dependent quark distribution defined first using gauge invariant renormalisation and second via the partially conserved axial vector current. In the first case, gluonic information is intrinsically built into the definition of the “polarised quark” via the axial anomaly. In the second, the notion of “polarised quark parton” is no longer gauge invariant. Taken together, this means that one has to be careful how far one carries through the semi-classical interpretation of parton distributions. The gluonic information built into the flavour-singlet part of $`\stackrel{~}{F}_q`$ is manifest through the massive $`\eta ^{}`$ and absence of any (nearly-)massless pseudo-Goldstone axial U(1) pole in $`\stackrel{~}{E}_q`$.

The zero-modes associated with dynamical $`U_A(1)`$ symmetry breaking mean that the perturbative QCD formula (5) for deeply virtual Compton scattering may not necessarily exhaust the total cross-section. For example, in Crewther’s theory of quark-instanton interactions one would expect $`\delta (x)`$ and $`\delta (x\zeta )`$ contributions in the flavour-singlet part of $`\stackrel{~}{H}_q(x,\zeta ,t)`$ associated the transfer of incident axial-charge carried by quarks and antiquarks respectively from partons with finite momentum to zero-modes. Such zero-mode contributions are not easily reconcilable with the perturbative QCD factorisation expression (5), in particular the singularities in the leading-twist coefficient function $`\stackrel{~}{C}_q`$ in Eq.(6). This problem is not relevant to isovector pion or rho production  which is described just in terms of the isovector non-forward distributions. Any flavour-singlet zero-modes will not contribute to these processes.

## 6 Orbital angular momentum

A sum-rule  relates the form-factors appearing in $`F_q`$ and $`F_g`$ in the spin-independent part of $`𝒜_{\mathrm{DVCS}}`$ to the quark and gluonic total angular-momentum contributions to the spin of the proton. The second moments of $`F_q`$ and $`F_g`$ project out the non-forward matrix elements of the QCD energy momentum tensor:

$`P^{}|T_{q,g}^{\mu \nu }|P=\overline{U}(P^{})[A_{q,g}(t)\gamma ^{(\mu }\overline{P}^{\nu )}+B_{q,g}(t)\overline{P}^{(\mu }i\sigma ^{\nu )\alpha }\mathrm{\Delta }_\alpha /2M`$ (28)
$`+C_{q,g}(t)\mathrm{\Delta }^{(\mu }\mathrm{\Delta }^{\nu )}/M+𝒪(\mathrm{\Delta }^3)]U(P)`$
There are no massless bosons which couple to $`T^{\mu \nu }`$ in QCD, which is associated with the fact that Poincare invariance is not spontaneously broken in QCD (with corresponding Goldstone bosons). This means that the expansion in $`\mathrm{\Delta }`$ and the forward limit of the form-factors in Eq.(28) are well-defined: there are no $`\mathrm{\Delta }_\mu \mathrm{\Delta }_\nu /\mathrm{\Delta }^2`$ terms. The current associated with Lorentz transformations is

$$M_{\mu \nu \lambda }=z_\nu T_{\mu \lambda }z_\lambda T_{\mu \nu }$$
(29)
The total angular momentum operator is related to the energy-momentum tensor by

$$J_{q,g}^z=P^{},\frac{1}{2}|d^3z(\stackrel{}{z}\mathrm{x}\stackrel{}{T}_{q,g})^z|P,\frac{1}{2}$$
(30)
Substituting Eq.(28) into Eq.(30) and taking the forward limit $`\mathrm{\Delta }0`$ one obtains

$$J_{(q,g)}^z=\frac{1}{2}\left[A_{q,g}(0)+B_{q,g}(0)\right]$$
(31)
This result is Ji’s sum-rule  relating total angular momentum to the form-factors measured appearing in Eqs.(7) and (8).

If one can extract $`J_{q,g}`$ from the forward limit of $`A(t)`$ and $`B(t)`$ in hard exclusive processes, then subtracting $`S_z`$ from polarised deep inelastic scattering would give information about the orbital angular momentum $`L_z`$. Whereas the intrinsic spin $`S_z^q`$ is sensitive to the axial anomaly, the total angular momentum $`J_z^{(q+g)}`$ is not because it is measured by a conserved current. Theoretical studies  show that $`J_z^q`$ and $`J_z^g`$ are each anomaly free in perturbative QCD meaning that the axial anomaly cancels between $`S_z^q`$ and $`L_z^q`$ in perturbation theory .

This result generalises beyond perturbation theory  with the added consequence that there is no zero-mode contribution to $`J_q`$ and $`J_g`$. To see this, first consider the crossing symmetry in $`x`$ of the spin-independent quark distributions. The second moment of the charge parity plus distribution $`(q+\overline{q})(x)`$ projects out the nucleon matrix element of the energy-momentum tensor $`T_{\mu \nu }`$. Any zero mode contribution to the right-hand side of the energy-momentum sum-rule would be associated with a $`\delta ^{}(x)`$ term in $`(q+\overline{q})(x)`$. Less singular $`\delta (x)`$ terms may be induced in the spin-dependent distributions by quark instanton interactions. I know of no model which predicts a stronger singularity like $`\delta ^{}(x)`$. Phenomenologically, such a term in the spin-independent structure function $`F_2`$ would lead to a violation of the energy-momentum sum-rule for partons, which is not observed. Since there is no zero-mode contribution to $`J_z`$ it follows that any zero-mode contribution to $`g_A^{(0)}`$ is compensated by a second zero-mode with equal magnitude but opposite sign in the orbital angular momentum $`L_z^q`$. This has the practical consequence that one has to be careful how one interprets any determination of $`L_z^q`$ from DVCS through the sum-rule (31). The orbital angular momentum carried by constituent quarks and by “current quark” partons are not necessarily the same, and is distinguished by the topological term $`𝒞`$ measurable in $`\nu p`$ elastic scattering. What happens to spin and orbital angular momentum in the transition from current to constituent quarks is intimately related to the dynamics of axial U(1) symmetry breaking.

Finally, we note that in perturbative QCD one also has to be careful to quote any value of “$`L_z^q`$” with respect to the factorisation scheme used to extract “$`S_z^q`$” from polarised deep inelastic data. For example, $`\mathrm{\Delta }q_{\overline{\mathrm{MS}}}=(\mathrm{\Delta }q\frac{\alpha _s}{2\pi }\mathrm{\Delta }g)_{\mathrm{AB},\mathrm{JET}}`$ – see also Shore .

## 7 Conclusions

In summary, the axial anomaly is manifest in deeply virtual Compton scattering directly through the flavour-singlet spin-dependent part of the scattering amplitude $`𝒜_{\mathrm{DVCS}}`$ through the non-forward parton distribution $`(\stackrel{~}{F}_u+\stackrel{~}{F}_d+\stackrel{~}{F}_s)`$ and indirectly through the interpretation of information about quark orbital angular-momentum extracted from the spin-independent part of $`𝒜_{\mathrm{DVCS}}`$.

Whether QCD factorisation provides a complete description of the flavour-singlet, spin-dependent part of $`𝒜_{\mathrm{DVCS}}`$ depends on possible zero-mode contributions to $`(\stackrel{~}{F}_u+\stackrel{~}{F}_d+\stackrel{~}{F}_s)`$ which may be generated by the dynamics of axial U(1) symmetry breaking. In the absence of such zero-mode contributions perturbative QCD factorisation works for this term provided one chooses a gauge-invariant factorisation procedure such as minimal subtraction. However the resultant spin-dependent non-forward distribution will not have a simple canonical interpretation in terms of extracting and re-inserting a quark or antiquark in the target because gluonic information is intrinsic to $`(\stackrel{~}{F}_u+\stackrel{~}{F}_d+\stackrel{~}{F}_s)`$ through the anomaly. This gluonic information is manifest in the large-mass $`\eta ^{}`$ pole in $`(\stackrel{~}{E}_u+\stackrel{~}{E}_d+\stackrel{~}{E}_s)`$ associated with the flavour-singlet pseudoscalar form-factor.

QCD anomaly effects cancel between the intrinsic spin and orbital angular-momentum contributions to the total quark angular-momentum in the proton. This extends to possible zero-mode contributions associated with dynamical $`U_A(1)`$ symmetry breaking. In QCD some fraction of the intrinsic spin carried by low-energy constituent quark quasi-particles may be carried by a zero-mode in addition to the partonic contributions from finite Bjorken $`x`$. This zero-mode is compensated by a zero-mode with equal magnitude and opposite sign in the orbital angular-momentum so that the total quark angular-momentum is independent of the details of the current to constituent quark transition through $`𝒞`$ whereas the separate intrinsic and orbital contributions are not. Finally, any value of the orbital angular-momentum extracted from future DVCS experiments should also be quoted with respect to the perturbative QCD factorisation scheme used to extract the partonic intrinsic spin from polarised deep inelastic data.

Acknowledgements

It is a pleasure to thank S.J. Brodsky, A. Donnachie, P.V. Landshoff, G. Piller, A. Radyushkin and W. Weise for helpful discussions. This work was supported in part by DFG.


text:
# Statistical Approach to Gene Evolution

## Abstract

The evolution in coding DNA sequences brings new flexibility and freedom to the codon words, even as the underlying nucleotides get significantly ordered. These curious contra-rules of gene organisation are observed from the distribution of words and the second moments of the nucleotide letters. These statistical data give us the physics behind the classification of bacteria.

PACS numbers: 87.10.+e, 87.15.-v, 05.40.+j

Over the years the statistical approach to genes has become prominent. The hidden Markov models are used in the alignment routines of biological sequences . For the secondary structures of the sequences stochastic context-free and context-sensitive grammars  are applied . The recent discovery of the fractal inverse power-law correlations  in these biological chains have led to ideas that statistically these sequences have features of music and languages \[5-7\]. Languages evolve with time. The vocabulary increases; the rules that dominate get progressively optimised so the order and information content is more. The purpose of this work is to track the statistical basis of the evolution in the coding DNA sequences (CDS).

The CDS are multiple of 3-tuples, the codons. The nucleotides adenine (A), cytosine (C), guanine (G) and thymine (T) taken in groups of three work to build the amino acid chains called proteins. The word-structure of CDS is, therefore, well known. We want to study evolution in terms of these words, their distributions and the moments.

It is known that any prose does not carry all the ingredients of evolution of languages. Similarly the CDS of any gene does not have all the salient features that accompany change. The genes that are present in the whole range of organisms, from the lowest bacteria to the highest mammals, and therefore connected to fundamental life processes are normally considered to be best suited to function as evolutionary markers. With this in view we choose glyceraldehyde-3-phosphate dehydrogenase (GAPDH) CDS for its ubiquitous presence in all living beings. The enzyme it codes for catalyses one of the crucial energy-producing steps of glycolysis, the common pathway for both aerobic and anaerobic respiration.

Distribution of words is studied for languages. The frequency of words is plotted against the rank. Here the total number of occurrences of a particular word is termed its frequency. The word most frequent has rank=1, the next most has rank=2, and so on. For natural languages, the plot gives the Zipf  behaviour:

$$f_N=\frac{f_1}{N}$$
(1)
where $`N`$ stands for the rank and $`f_1`$ and $`f_N`$ are the frequencies of words of rank 1 and $`N`$ respectively. The Zipf-type approach to the study of DNA has brought methods of statistical linguistics into DNA analysis . The generalized Zipf distribution of n-tuples has provided hints that the DNA sequences may have some structural features common to languages. In this work we confine ourselves to the distribution of 3-tuples, the codons, in the CDS. The words, therefore, are non-overlapping and on the single reading frame.

The frequency-vs-rank plot of the codon words show that these distributions, given the frequency of rank 1 and the length of the sequence, are almost completely defined through the universal exponential functional form :

$$f_n=f_1.e^{\beta (N1)}$$
(2)
The parameter, called $`\beta `$, is determined by the ratio

$$\beta \frac{f_1}{L}$$
(3)
$`\beta `$ measures the frequency of rank 1 per unit length of the sequence. The exponential form (2) is to be compared to the usual Boltzmann distribution. The rank of the word is akin to energy; $`\beta `$ is analogous to inverse temperature. The relationship (3) that $`\beta `$ is frequency of rank 1 per unit length is supported well from data . The analogy between word distributions and the classical Boltzmann concepts goes deeper. A decrease in $`\beta `$, from (3), implies frequency of rank 1 per unit length goes down. In that case the vocabulary clearly increases. More words are used, thereby more states are accessed. For the GAPDH CDS we find the evolution is driving it to higher temperatures; into more freedom for words, into more randomisation. $`\beta `$ evolves monotonically.

Underneath, however, there runs a curious counterflow. Suppose we look into the nucleotides that constitute the sequence, once again in windows of size 3 and in the same reading frame. First, we ask how much order there is in the sequence. To find out we study the second moments of the letters A, C, G and T. These second moments, by themselves, do not produce any pattern. The GAPDH CDS has about 1000 bases. For each organism the proportions of A, C, G and T in the GAPDH CDS are different. This strand-bias, interestingly, masks a remarkable underlying trend.

To get there the strand-bias has to be eliminated. The order in the sequence, we assume, is its deviation from the random. We define the quantity $`X`$, a measure of this deviation, as follows:

$`X`$ = $`\frac{SecondMomentoftheBaseDistributioninGAPDHCDS}{SecondMomentoftheBaseDistributionintherandomsequencewithidenticalstrandbias}`$

Normalised as above, the effect of the strand bias is unmasked. $`X`$ values of GAPDH change monotonically with evolution. The data tells us there is an increase in persistence amongst the letters (in windows of size 3) with evolution in the CDS .

The evolution in the GAPDH CDS is then the result of these two contra trends: while words acquire greater uniformisation, the underlying letters have more order. The monotonic behaviours of $`\beta `$ and $`X`$ with evolution give us the physics behind the biological classification of bacteria.
$`\mathrm{𝐌𝐞𝐭𝐡𝐨𝐝𝐬}`$
$`\underset{¯}{WordDistributions}`$
For the codons it is known  the exponentials give somewhat better fits over the usual power laws. The exponential form, equation (2), is characterized by the parameter $`\beta `$. The quantity has some universal features in that it is almost completely determined by $`f_1`$ and the length of the CDS. The relationship

$$\beta =\frac{f_11}{L}+\frac{1}{2}.\frac{(f_11)^2}{L^2}$$
(4)
is known to fit observations on diverse genes. For the bacterial GAPDH CDS the results of $`\beta `$ are given in Table 2.
$`\underset{¯}{Moments}`$
Consider the 4-dimensional walk model  such that A, C, G and T correspond to unit steps, in the positive direction, along $`X_A`$, $`X_C`$, $`X_G`$ and $`X_T`$ axes. After $`n`$-steps if the co-ordinate of the walker is ($`n_A`$, $`n_C`$, $`n_G`$, $`n_T`$), then, clearly,

$$n=n_A+n_C+n_G+n_T$$
(5)
and $`n_i`$ ($`i`$ $``$ A,C,G,T), is the number of nucleotide of type $`i`$ in the sequence just walked.

If the sequence has $`n`$ bases, and $`n_i`$ is the number of base of type $`i`$, the strand bias of the sequence is the proportion of $`n_i`$ in $`n`$, defined as

$$p_i=\frac{n_i}{n}$$
(6)
The probability distribution for the single step in this 4-d walk is

$$P_1(x)=\underset{i}{}p_i\delta (x_i1)$$
(7)
where $`\delta `$ is the usual $`\delta `$-function of Dirac. The characteristic function of the step is the Fourier transform of equation (7),

$$P_1(k)=\underset{i}{}p_ie^{ik_i}$$
(8)
The characteristic function of $`l`$ steps

$$P_l(k)=[P_1(k)]^l$$
(9)
The second moments (i.e. the average values) of distributions may be obtained taking derivatives of $`P_l(k)`$ with respect to $`k`$. Thus for the random sequence (indicated by the subscript $`r`$) with the strand bias (6), we get the average values:

$$<n_i^2>_r=l[(l1)p_i^2+p_i]$$
(10)
$$<n_in_j>_r=l(l1)(p_i.p_j)(ij)$$
(11)
We are interested in codons, therefore, the window size $`l`$ in equations (10) and (11) is chosen to be 3. For the actual sequences we calculate $`<n_i^2>_{seq}`$ and $`<n_in_j>_{seq}`$. The quantities

$$X_D=\frac{<n_i^2>_{seq}}{<n_i^2>_r}$$
(12)
\[where D $``$ AA,CC,GG,TT\]
and

$$X_{OD}=\frac{<n_in_j>_{seq}}{<n_in_j>_r}(ij)$$
(13)
\[where OD $``$ AC,AG,AT,CG,CT,GT\]
measure the deviation of the diagonal and off-diagonal second moments of the sequence to those of the random sequence of identical strand bias respectively. Finally, we come up with an over-all averaged index, $`X`$, given by

$$X=\frac{X_D+X_{OD}}{10}$$
(14)
This $`X`$ provides a measure of the order in the sequences.

$`\mathrm{𝐎𝐛𝐬𝐞𝐫𝐯𝐚𝐭𝐢𝐨𝐧𝐬}\mathrm{𝐚𝐧𝐝}\mathrm{𝐑𝐞𝐬𝐮𝐥𝐭𝐬}`$
To set the basis for what we discuss later, we begin by recording the $`\beta `$ and the $`X`$ values of higher organisms, the eukaryotes (Table 1). We confine our discussion of the eukaryotes to three broad categories: fungi, invertebrates and vertebrates. It is known  from fossil records the oldest fungi came about 900 million years (Myr) before present (bp). The oldest fungal species, identified with certainty, are from the Ordovician period, i.e., some 500 Myr bp. The fossil records of invertebrates suggest this group came about the same time as the fungi. The vertebrates came later, about 400 Myr bp, in late Ordovician and Silurian period.

Let us look at the $`\beta `$ and the $`X`$ values of these eukaryotic groups. Fungi has the highest, followed by invertebrates, while for the vertebrates the $`\beta `$ and the $`X`$ reach minima. We conclude the $`\beta `$ and the $`X`$ decrease with evolution. The data further suggest fungi and invertebrates came about the same time and underwent parallel evolution, while the representives of the vertebrate group came later in the evolutionary line-up.

Having set the basis, let us now look at 14 bacterial species from three groups: cyanobacteria, proteobacteria (that includes vast majority of gram-negative bacteria), and the $`\mathrm{𝐵𝑎𝑐𝑖𝑙𝑙𝑢𝑠}/\mathrm{𝐶𝑙𝑜𝑠𝑡𝑟𝑖𝑑𝑖𝑢𝑚}`$ group, a type of gram-positive bacteria. Table 2 summarises the $`\beta `$ and the $`X`$ values of these samples. These bacterial groups arose during the Precambrian period of geological time-scale, but there are several schools of thought regarding their specific times of origin within this period.

We approach the bacterial GAPDH CDS with two differing statistical measures, the $`\beta `$ and the $`X`$. Interestingly, both give us almost identical trends (Figs. 1 and 2). $`\mathrm{𝐿𝑎𝑐𝑡𝑜𝑏𝑎𝑐𝑖𝑙𝑙𝑢𝑠}\mathrm{𝑑𝑒𝑙𝑏𝑟𝑢𝑒𝑐𝑘𝑖𝑖}`$, a member of the $`\mathrm{𝐵𝑎𝑐𝑖𝑙𝑙𝑢𝑠}/\mathrm{𝐶𝑙𝑜𝑠𝑡𝑟𝑖𝑑𝑖𝑢𝑚}`$ group, has the highest $`\beta `$ and $`X`$ values (Table 2). There is then a large measure of overlap between the $`\mathrm{𝐵𝑎𝑐𝑖𝑙𝑙𝑢𝑠}/\mathrm{𝐶𝑙𝑜𝑠𝑡𝑟𝑖𝑑𝑖𝑢𝑚}`$ group and the proteobacteria (Figs. 1 and 2). The extent of overlap of the $`\beta `$ values is somewhat more than that of the $`X`$. The cyanobacterial samples have the minimum values of the $`\beta `$ and the $`X`$. There is no overlap between the cyanobacterial values of the $`\beta `$ and the $`X`$ with the $`\mathrm{𝐵𝑎𝑐𝑖𝑙𝑙𝑢𝑠}/\mathrm{𝐶𝑙𝑜𝑠𝑡𝑟𝑖𝑑𝑖𝑢𝑚}`$ group. The overlap between the proteobacteria and the cyanobacteria is small. Only one proteobacterial sample, $`\mathrm{𝐵𝑟𝑢𝑐𝑒𝑙𝑙𝑎}\mathrm{𝑎𝑏𝑜𝑟𝑡𝑢𝑠}`$ has greater $`\beta `$ value than the cyanobacterial member, $`\mathrm{𝑆𝑦𝑛𝑒𝑐ℎ𝑜𝑐𝑦𝑠𝑡𝑖𝑠}`$ sp. (strain PCC 6803).

The averages of the $`\beta `$ or the $`X`$ has the maximum value in the $`\mathrm{𝐵𝑎𝑐𝑖𝑙𝑙𝑢𝑠}/\mathrm{𝐶𝑙𝑜𝑠𝑡𝑟𝑖𝑑𝑖𝑢𝑚}`$ group, followed by the proteobacteria, while the cyanobacteria samples have the lowest values. In line with our observations on the eukaryotes, we propose (Figs. 1 and 2) that the $`\mathrm{𝐵𝑎𝑐𝑖𝑙𝑙𝑢𝑠}/\mathrm{𝐶𝑙𝑜𝑠𝑡𝑟𝑖𝑑𝑖𝑢𝑚}`$ group originated some time before the proteobacterial species, but later both groups evolved in parallel. The cyanobacterial samples are of recent origin compared to these groups. The trends in the $`\beta `$ and the $`X`$ give us identical patterns that segregate the bacterial species into groups. Amusingly, the results seem to be in agreement with what is accepted so far regarding the phylogenetic relationships among these three groups . Our study of the GAPDH CDS, its word distributions, and the moments gives us the physics underlying evolution.

The decrease in $`\beta `$ with evolution for the GAPDH CDS tells us that evolution is taking the gene progressively towards higher temperatures. The $`\beta `$ value, we recall, is the frequency of rank one per unit length. Lowering of the $`\beta `$ implies less dominance of the maximum weight. In consequence, the other words enjoy greater freedom, the vocabulary increases and more states are accessed. In a sense the evolution in the GAPDH CDS mirrors Boltzmannian statistics. Even though the GAPDH CDS has evolved in a complex evolutionary regime in contact with environment, the Boltzmannian behaviour is useful. For instance, it allows us to define the word-entropy of the CDS. That gives us a measure of the information content of the words in biological chain.

At the level of the nucleotide letters A, C, G and T, the order is measured by the quantity $`X`$. As we look into the diagonal averages $`X_D`$, (12), we find it increases with evolution. For the window of size 3, this growing diagonal moment implies a rising persistent correlation. In consequence, the off-diagonal averages $`X_{OD}`$, (13), go down, decreasing antipersistence. Looked at from the letters, the sequences become less uniform and deviate more from the random sequence of identical strand bias. The order, or the information, in the arrangement of letters shows a rising trend with evolution.

Does any CDS that is an evolutionary marker evolve in ways similar to the GAPDH? We have worked with the CDS of some other glycotic enzymes, such as phosphoglycerate kinase, and found they behave similarly. Other evolutionary markers such as the ribulose-1,5-bisphosphate carboxylase/oxygenase enzyme large segment (rbcL) show similar behaviour. We use these data for biological subclassification. The CDS for ribosomal RNA is another class of sequence that is being investigated. It does not code for protein, but for RNA, and has periods other than 3. The 3 period does exist, but is not predominant.

Sequence modeling has recently become important. The fractal correlations in the sequences led to the expansion-modification system . Later came the insertion models . Here the differences in the CDS and non-coding sequences were observed and the non-coding sequences modeled. The unifying models of copying-mistake-maps  modeled both the coding and the non-coding parts. In these models the statistical features of the non-coding sequences have received emphasis. The evolutionary features of the GAPDH CDS isolates the statistical aspects that underlie evolution in coding sequences. The statistics of the word distributions and the subtle cross current of the second moments, we hope, will lead further in these efforts.
$`\mathrm{𝐀𝐜𝐤𝐧𝐨𝐰𝐥𝐞𝐝𝐠𝐦𝐞𝐧𝐭𝐬}`$
S.C. thanks Professor Anjali Mookerjee for many discussions. W.A.K. is supported by the John Fulbright foundation in the laboratory of J.C.

Electronic address: tpsc@mahendra.iacs.res.in

$`\mathrm{𝐅𝐢𝐠𝐮𝐫𝐞}\mathrm{𝐋𝐞𝐠𝐞𝐧𝐝𝐬}`$
Figure 1. The average $`\beta `$ values for the GAPDH CDS from three bacterial groups (see Table 3). The error bars indicate the standard deviation from the average values.
Figure 2. The average $`X`$ values for the GAPDH CDS from three bacterial groups (see Table 3). The error bars indicate the standard deviation from the average values.