{'fraction_non_alphanumeric': 0.03340441339650132, 'fraction_numerical': 0.010121672343422842, 'mean_word_length': 4.501645609937361, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Gauge-Invariant Formalism of Cosmological Weak Lensing

We present the gauge-invariant formalism of cosmological weak lensing, accounting for all the relativistic effects due to the scalar, vector, and tensor perturbations at the linear order. While the light propagation is fully described by the geodesic equation, the relation of the photon wavevector to the physical quantities requires the specification of the frames, where they are defined. By constructing the local tetrad bases at the observer and the source positions, we clarify the relation of the weak lensing observables such as the convergence, the shear, and the rotation to the physical size and shape defined in the source rest-frame and the observed angle and redshift measured in the observer rest-frame. Compared to the standard lensing formalism, additional relativistic effects contribute to all the lensing observables. We explicitly verify the gauge-invariance of the lensing observables and compare our results to previous work. In particular, we demonstrate that even in the presence of the vector and tensor perturbations, the physical rotation of the lensing observables vanishes at the linear order, while the tetrad basis rotates along the light propagation compared to a FRW coordinate. Though the latter is often used as a probe of primordial gravitational waves, the rotation of the tetrad basis is indeed not a physical observable. We further clarify its relation to the E-B decomposition in weak lensing. Our formalism provides a transparent and comprehensive perspective of cosmological weak lensing.

1 Introduction Significant and impressive progress has been made in observational cosmology in recent decades. However, the standard cosmology is still full of mystery, demanding further explanations both in terms of theoretical understanding and observational breakthrough. In particular, the late-time cosmic acceleration due to dark energy poses great challenges to theorists and observers alike. To tackle this pressing problem, numerous large-scale surveys have been undertaken to map the matter distribution of the Universe over the Hubble volume out to high redshifts. The current generation of these surveys include the Kilo Degree Survey (KiDS; [1]), which is a wide-field program being carried out using the VLT Survey Telescope, the Dark Energy Survey (DES; [2]), which is nominally a five year survey that will cover 5000 square degrees using the Blanco 4 m telescope, and the Subaru Hyper Suprime-Cam survey, 2 which will cover 1400 square degrees. These experiments represent a substantial improvement over previous measurements and are a significant stepping-stone to the next generation of experiments, known as stage-IV. These future programs include the Large Synoptic Survey Telescope (LSST; [3]), an 8-m dedicated ground based facility, Euclid [4] and the Wide Field Infrared Survey Telescope (WFIRST; [5]), which are two planned space based missions. Together these experiments will effectively cover the full observable extra-galactic sky and measure the shapes of roughly a billion galaxies. One of the main goals of these cosmological experiments is to measure gravitational lensing. Gravitational lensing deals with light propagation in the Universe (see, e.g., [6][7][8][9][10][11] for recent reviews). As light travels towards us from distant objects, its path is perturbed due to intervening mass. In the case when these perturbations are small, the path is typically affected by several structures along the line-of-sight. This regime is known as weak lensing or cosmic shear and leads to distortions in the observed properties of the distant sources. Since first detection in the late 1990s [12][13][14], the field of weak gravitational lensing has matured into one of the core cosmological probes. Current results include the KiDS weak lensing analysis of a 450 deg 2 survey [15], which reported a 5% precision on the parameter S 8 , a combination of σ 8 and Ω m , where Ω m is the cosmic matter density parameter and σ 8 is the amplitude of matter density fluctuations today. More recently the DES analysis of the first year data (Y1) [16], which covers an area of 1321 deg 2 , achieved 3% precision on S 8 . This steady improvement in precision will increase and should accelerate as we enter the era of stage IV experiments. However, the ambitious goals enabled by these impressive observational facilities can only be achieved, if the theoretical predictions are at the same level of precision as those set by observations.
Sachs [17,18] was the first to develop the basic foundation for the propagation of the gravitational waves, and building on this framework the early work on weak lensing was performed by Gunn [19]. With observational progress in 1990s, cosmological weak lensing received large attention, and many researchers [20][21][22][23][24][25][26] established the standard weak lensing formalism, which is further completed in the following years [27][28][29][30][31]. However, despite these theoretical and observational developments in the past decades, there has not been a complete gauge-invariant description of cosmological weak lensing. Here we derive the gauge-invariant formalism of cosmological weak lensing. The key ingredient for such a task is to identify the missing physics in the standard weak lensing formalism, which makes the predictions for the lensing observables gauge-dependent. In recent years, great attention has been paid to the relativistic description of galaxy clustering [32][33][34][35][36][37] (see [38] for review), by which the gauge ambiguities in the standard theoretical predictions are highlighted. The dominant contribution to galaxy clustering is the matter density fluctuation, and it is well known that the matter density power spectrum differs significantly near the horizon among different choices of gauge condition. It was quickly understood and argued [32] that the standard model for galaxy clustering was incomplete and our (correct) theoretical descriptions for any cosmological observables should be gauge-invariant; The key reason for such inadequacy in the standard model was that "unobservable quantities" are used to build theoretical predictions. For example, the observed redshift z is a physical observable, but our description of the radial position (or sometimes called "true redshift"z without the redshift-space distortion) is gauge-dependent. In gravitational lensing, the observed angular positionn of the source is physical, but our description of the (unlensed) "true angular position"ŝ is again gauge-dependent. The standard weak lensing formalism built in terms of such unobservable quantities is inevitably incomplete. In fact, it is well known that the lensing convergence κ in Eq. (4.19) is gauge-dependent and hence it cannot be directly associated with the physical lensing observable we measure. It is shown [39,40] that the lensing convergence we measure is the (gauge-invariant) fluctuation δD in Eq. (5.39) in the luminosity distance (or the angular diameter distance), which includes not only the standard lensing convergence κ, but also other relativistic contributions. Here we present a complete and coherent description of the gauge-invariant lensing formalism.
Establishing the proper descriptions of the observable quantities such as the observed angle and redshift is the first step, and this requires the specification of the observer, in particular, the observer frame, in which the metric is Minkowski η ab . The second step is to establish the proper descriptions of the physical quantities at the source position that will be measured by the observer, and this step also requires the specification of the source and its rest-frame, where the physical quantities such as the size and the shape of galaxies are defined. For both cases, the tetrad basis provides such a link needed to connect the light propagation in the Friedmann-Robertson-Walker (FRW) universe to the rest-frames of the observer and the source. The standard lensing formalism lacks such descriptions. In Sec. 2, we introduce the tetrad basis and present the linear-order expression for the tetrad vectors. Furthermore, there exists one more critical ingredient missing in the standard lensing formalismcheck of gauge-invariance. In the lensing literature, a gauge condition is adopted (almost exclusively the conformal Newtonian gauge), and all the calculations are performed with that gauge condition. However, by fixing the gauge condition, one loses the ability to check the gauge-invariance of theoretical predictions computed in that gauge condition. Here we perform all our calculations with the general metric representation without choosing any gauge conditions, such that we can explicitly check at each step the gauge-invariance of our expressions. This procedure greatly helps understand better which part in our theoretical descriptions can be associated with the physically meaningful quantities.
The organization of the paper is as follows: In Sec. 2, the basic observables in the observer rest-frame are expressed in relation to the photon wavevector, and a local tetrad basis is constructed to provide the connection to the photon wavevector in a FRW coordinate. In Sec. 3, we solve the photon geodesic equation, accounting for all the relativistic effects and paying particular attention to the subtleties at the observer position. The source position in a FRW coordinate is geometrically decomposed to represent the deviation from the observationally inferred position. In Sec. 4, we generalize the standard weak lensing formalism by using the full geodesic equation, we demonstrate that the standard formalism is still incomplete and gauge-dependent. In Sec. 5 we present our main results of the gauge-invariant weak lensing formalism. Using the tetrad basis at the source position, we construct the distortion matrix in terms of the physical size and shape in the source rest-frame, and we compute the lensing observables in Sec. 5.1. In Sec. 5.2, we demonstrate that the tetrad basis rotates as it is parallel transported along the photon path. However, we show that this rotation is gauge-artifact and not observable. In Sec. 5.3 we present the lensing E-B decomposition and clarify its relation to the lensing rotation. In Sec. 6, we compare our results to previous work on the lensing effects due to the gravitational waves in Sec. 6.1, on the lensing formalism by a standard ruler in Sec. 6.2, and on the lensing formalism by a Jacobi mapping approach in Sec. 6.3. We summarize our findings and discuss the implications of our new formalism for upcoming surveys in Sec. 7.
Throughout the paper, we use the Greek indices µ, ν, ρ, · · · for the space-time components in a FRW coordinate with metric g µν , in which we use the conformal time η and the Greek indices α, β, γ, · · · to represent the time and the spatial components, respectively. The Latin indices a, b, c, · · · are used to represent the internal components in a rest-frame with the Minkowski metric η ab , in which we use the proper-time t and the Latin indices i, j, k, · · · to represent the time and the spatial components. We summarize our notation convention in Table 1.

Observables in the Light Propagation
We derive the expressions for the basic quantities associated with the light propagation in terms of the metric perturbations, clarifying the difference between the FRW frame and the observer frame. In the former the light propagation is computed, and in the latter the local observables are defined and measured.

Observer rest-frame: Tetrad basis
The observers perform cosmological observations in the observer rest-frame, in which the metric is Minkowski η ab and the time direction is set by the four velocity u µ of the observer. This frame is defined only in the infinitesimal neighborhood of a given spacetime point of the observer. However, a tangent space orthogonal to the time direction u µ of the observer can be well defined, by constructing three spacelike vectors e µ i , where the index i represents three spatial directions of the observer (i = 1, 2, 3). Together with the time direction e µ 0 ≡ u µ , these four vectors are referred to as a tetrad e µ a , forming an orthonormal basis of the observer (η ab = g µν e µ a e ν b ). The Latin indices (a, b, c, · · · = 0, 1, 2, 3) are used to represent the component of a tetrad, and they are raised and lowered with the Minkowski metric.
Cosmological information is measured by observing light from a distant source, and the tetrad vectors serve as a basis for this cosmological observation. Consider a null vector k µ in a given FRW coordinate, describing the light propagation. This photon wavevector is measured in the observer rest-frame as where we expressed the components of the photon wavevector in the observer rest-frame, in terms of the observable quantities: the angular frequency ω = 2πν of the photon and the angular positionn of the source. In the observer rest-frame, a set of angles (θ, φ) is assigned to the unit directional vector n i = (sin θ cos φ, sin θ sin φ, cos θ). Trivially, these cosmological observables (e.g., ω,n, and so on) are independent of FRW coordinates, while the components of the photon wavevector k µ or the tetrad vectors e µ a are coordinate dependent. This diffeomorphism invariance of any cosmological observables was emphasized by Yoo and Durrer [41] in conjunction with their gauge-transformation properties beyond the linear order in perturbations, and it can be readily inferred from the above equation as the coordinate indices of cosmological observables are contracted with the coordinate indices of the tetrad vectors and only the internal indices of the tetrad vectors remain. It is emphasized [42] that the coordinate indices µ, ν, · · · and the internal tetrad indices a, b, · · · transform according to their own (different) symmetry groups.

Tetrad basis vectors in a FRW coordinate
In general, the four tetrad vectors in four dimensions have sixteen degrees of freedom, and ten of which are constrained by the metric tensor g µν : The remaining six degrees of freedom belong to the Lorentz symmetry at a given spacetime, i.e., three Lorentz boosts and three spatial rotations. While the tetrad formalism is completely general in choosing the remaining symmetry, the Lorentz boosts in our case are already fixed by the observer four velocity e µ 0 = u µ . This gauge choice is natural for describing the cosmological observables, and it applies to all spacetime points, defining the rest-frame of not only the observer, but also any other "observers," including the sources of our cosmological observables (see [42] for a complete description of the tetrad formalism in cosmology).
The spatial directions e µ i can be arbitrary due to the remaining symmetry in rotation in the restframe. For convenience, however, we fix the spatial symmetry by aligning the spatial tetrad directions with the FRW coordinate directions in a homogeneous universe as where δ α i is a Kronecker delta and we assumed a flat FRW universe and chose a rectangular coordinate (see Appendix A for the metric convention). With such choice, there exists no further remaining gauge symmetry in the tetrads. Indeed, this choice is implicitly made in most work in cosmology, which makes calculations simpler. However, this choice often compounds the internal and the coordinate components, which makes the calculations vulnerable to errors, when perturbations are considered. We will clarify what errors were made in previous work, with particular attention to the difference in indices -The internal components are diffeomorphism invariant and related to what we measure.
In our real universe, we have to consider the deviations from the homogeneous FRW metric. The time direction e µ 0 of the tetrad is set by the timelike four velocity field u µ as at the linear order in perturbations. The spatial directions of the tetrad can be parametrized as where δe η i and δe α i (or p j i ) are perturbations, capturing the deviation from the background. The perturbation δe α i in the spatial component of e µ i can be further split into the symmetric part S ij and the anti-symmetric part A ij as The spatial indices i, j are raised and lowered by the spatial part δ ij of the Minkowski metric. In equivalence, the perturbation δe α i can be split in terms of FRW coordinates. 3 The orthonormality condition in Eq. (2.2) constrains the time component of the spatial tetrads 8) and the symmetric part of the spatial component Despite the expression written in terms of the internal indices, the symmetric part S ij is not invariant under diffeomorphism, as the coordinate indices are contracted with the metric tensor C αβ , which changes under gauge transformations. The anti-symmetric part A ij is, however, left unconstrained by the orthonormality condition, which is why this part was unnoticed in previous work (e.g., see [36][37][38]).
To complete the derivation of the tetrad basis vectors, we consider their gauge-transformation properties. Being a four vector at every spacetime point, the tetrad vectors transform as vectors under a coordinate transformation in Eq. (A.4), and for an infinitesimal coordinate transformation by ξ µ , this relation dictates the gauge-transformation of the tetrad vectors: (2.10) 3 In complete generality, we can split the perturbation δe α i of the spatial-component of e µ i into the symmetric part S α β and the anti-symmetric part A α β as (2.7) This definition is as valid as Eq. (2.6).
where £ ξ is the Lie derivative. With the gauge-transformation of the metric perturbations in Appendix A, the anti-symmetric part has to transform as where the scalar part of the spatial transformation L α cancels and only the vector part L α remains. Given the constraint from the gauge transformation, the anti-symmetric part of the spatial tetrad vectors can be expressed in terms of the metric perturbations and the rotation of the spatial tetrad vectors as where ijk is the Levi-Civita symbol, the rotation vector Ω i captures the residual symmetry in spatial rotation and is invariant under diffeomorphism. The anti-symmetric part, often ignored in previous work, is composed of the vector perturbation C α and the spatial rotation Ω i of the tetrad vectors. Fortunately, the errors in missing the anti-symmetric part in the tetrad expressions are relatively innocuous -As a choice of internal gauge, the spatial rotation can be set zero Ω i = 0, and the perturbation calculations are performed in most cases without vector perturbations or by choosing a spatial gauge C α ≡ 0. In this work, we will keep the rotation in general, and the implication of the spatial rotation Ω i will be discussed at length in Sec. 5.2. Accounting for the symmetric and the anti-symmetric parts, we present the complete expression for the spatial tetrad vectors at the linear order [42]: Given the six degrees of freedom, we fixed the three Lorentz boosts by setting the timelike direction e µ 0 ≡ u µ as the four velocity field, but we left unspecified the remaining three spatial rotations Ω k in the spatial tetrad vectors at the perturbation level. All degrees of freedom are already fixed at the background.

Photon wavevector in a FRW coordinate
Having derived the local tetrad vectors, we are now in a position to relate the local observables to the photon wavevector k µ in the FRW frame. In the observer rest-frame, the observer measures the photon frequency ν and the angular positionn of the source. These basic observables can be used to construct a photon wavevector k a = e a µ k µ in the rest-frame as in Eq. (2.1), and they are invariant under diffeomorphism. Using the tetrad expression, we can derive the photon wavevector in a FRW coordinate 14) and the observed direction of the photon propagation in a FRW coordinate can be derived in a similar way as where the explicit expression can be readily inferred from Eq. (2.13). In the absence of perturbations, the spatial components of the photon wavevector k µ or the photon propagation direction n µ are proportional to the observed direction n i in the observer rest-frame, because the spatial tetrad vectors are by construction aligned with a FRW coordinate.
However, the presence of perturbations changes their expressions in a FRW coordinate. This is a general relativistic generalization of a Lorentz boost in special relativity, with which the observer moving with a relative velocity measures different frequency and different propagation direction. Compared to previous work [38], the difference in Eq. (2.14) arises solely due to the missing antisymmetric part A ij in the spatial tetrad vectors, and that is, the vector perturbation and the spatial rotation. We emphasize that all the quantities above are evaluated at the observer position and the expression is valid only at the observer position, because the angular frequency ω and the propagation direction n i are the quantities measured by the observer, not a field defined everywhere.
Nevertheless, it is useful to have such expressions along the photon path by generalizing the above equations. Given the observables (ω, n i ) in the observer rest-frame described by u µ , we completely fixed the photon wavevector k µ in Eq. (2.14) at the observer position by setting the orientation Ω i of the spatial tetrad e µ i . The photon wavevector is subject to the geodesic equation and the null condition, such that it is completely determined in a FRW coordinate along the null path described by the observables (ω, n i ). At any point x µ λ along the null path (parametrized by λ), we will need to specify an "observer" with four velocity u µ λ , defining the timelike direction, in which another "observation" will be performed. This observer at x µ λ will measure the frequency ω λ = − (u µ k µ ) λ in the rest-frame (different from ω at origin), and hence the "observed direction n µ λ " of this observer in a FRW coordinate is determined as in Eq. (2.15). However, as evident in Eq. (2.15), the observed direction n i λ in the observer rest-frame depends on the choice of the spatial tetrad e µ i at x µ λ , i.e., the rotation Ω i λ . This reflects the freedom to choose local coordinate directions, on which the observed angle (θ, φ) λ depends.

Photon Geodesic Path and Source Position in a FRW Coordinate
Here we solve the geodesic equation to obtain the source position in a FRW coordinate and derive the geometric distortions of the source position from the observed position. Compared to the previous work, we clarify the change in the calculations due to the missing anti-symmetric part of the spatial tetrad vectors and the spatial coordinate shift at the observer position.

Conformal transformation of the FRW metric
To facilitate the computation of the null geodesic path in a FRW coordinate, we perform a conformal transformation g µν →ĝ µν : (3.1) Since the null geodesic path (ds 2 = 0) remains unaffected by the conformal transformation, we can utilize the geodesic equation in the conformally transformed metric to derive the null path x µ . While the null path x µ (Λ) can be parametrized by any affine parameter Λ, we can physically fix the affine parameter Λ by demanding that the tangent vector along the path is the photon wavevector, 4 and Eq. (2.14) is satisfied at the observer position as the initial condition for the photon wavevector. In addition, the photon wavevector should meet the null condition 0 = k µ k µ and the geodesic equation 0 = k ν k µ ;ν at any point along the path.
With the conformal transformation, the geometry of the spacetime manifold changes, and the covariant derivatives in two different manifolds are not identical in order to satisfy their own metric compatibility: where quantities in the conformally transformed metric are represented with hat. The metric compatibility condition in two manifolds implies [43] that the covariant derivatives are related to each other with the connecting tensor aŝ Hence, the geodesic equation is not satisfied for the photon wavevector k µ with∇ µ in the conformally transformed metric. However, by re-parameterizing the photon path x µ (λ) with different affine parameter λ, we can derive the conformally transformed wavevectork µ for the same null path that satisfies the geodesic equation 0 =k ν∇ νk µ in the conformally transformed metric: where the proportionality constant C is left unconstrained in the conformal transformation, because the metric compatibility constrains only the derivative of dΛ/dλ [43]. Given the conformal transformation, the choice of the normalization C is completely free. With Eqs. (3.5) and (2.14), it appears natural to choose the normalization to fix the combination Caω that is constant everywhere in the background. 5 Therefore, we fix the normalization factor C by setting the product at the observer position [37,38,44] where the subscript o represents the observer position. The presence of perturbations makes the combination Caω vary as a function of position, while the normalization constant C is still a constant. With such condition, the conformally transformed photon wavevector can be parametrized aŝ and the perturbations to the photon wavevector at the observer position are then where we define the perturbation ∆ν in the observed frequency in the conformally transformed metric, in terms of the product because the four velocity in the conformally transformed metric isû µ = au µ andû µ =ĝ µνû ν = u µ /a, justifying the notation. Our choice of the normalization condition in Eq. (3.6) becomes ∆ν o = 0 at the observer position, but we keep the term ∆ν o in general. However, we stress that the choice of the normalization ∆ν in δn α o only affects the component in proportion to n i , while the rotation component is perpendicular to n i . With ∆ν o = 0, the perturbations to the photon wavevector at the observer position transform as where d/dλ = ∂ η − n i δ α i ∂ α is the derivative along the photon path in Eq. (3.24). Their gaugetransformation properties in general can be derived from the transformation of the photon wavevector k µ with additional constraint that the normalization condition ∆ν is imposed at one physical point p in any coordinate systems, i.e., Eq. (3.10) has the same value at the same physical point p: For our choice of the normalization condition, the fixed physical point p is the observer position. The corrections at p in Eq. (3.12) were neglected in Eq. (2.21) in [38].
In addition, we define the observed angular vector n α in a FRW coordinate as However, it should be noted that the internal index i is invariant under diffeomorphism. We also define two additional vectors θ α and φ α in terms of two unit directional vectors θ i and φ i perpendicular to n i in the observer rest-frame: Furthermore, we will use the following notation in connection to the directional vector n α : generally applicable to any vectors or tensors.

Observer position in a FRW coordinate
In a homogeneous universe, the spatial coordinate of the observer or any other observers can be set x α o = 0, and the (conformal) time coordinate is uniquely set to beη o in relation to the aget o of the Universe:η where we used bar to indicate that the coordinate position of the observer is obtained in a homogeneous universe. Indeed, this can be readily derived by considering the motion of a free-falling observer with four velocity u µ . The timelike four velocity can be parametrized in terms of the proper time τ measured by the observer in the rest-frame, and the FRW coordinates of the observer can be obtained by integrating the four velocity over the proper time: In a homogeneous universe, in which the observer four velocity is u µ = (1/a, 0), the time coordinateη o of the observer today can be obtained by setting τ i = 0 in the above equation.
In the presence of perturbations, the observer motion deviates from the static motion in a homogeneous universe. Consequently, the coordinates of the observer today also deviate fromx µ o = (η o , 0) in a homogeneous universe -the (time) coordinate lapse is [45] 18) and the (spatial) coordinate shift is [41] x where we changed the integration over the motion of the observer in terms of the proper time dτ to the integration over the time coordinate at a fixed spatial coordinate, valid at the linear order in perturbations. Using the geodesic condition of the observer motion the coordinate lapse can be further simplified as From the gauge transformation properties in Appendix A, we can readily show that these deviations in Eqs. (3.18) and (3.19) satisfy the transformation properties above. In the comoving-synchronous gauge (v = A = U α = 0), the observer position x µ o is identical to thatx µ o in a homogeneous universe. However, this is not valid in general. At the linear order in perturbations, the coordinate shift δx α o drops out in the expressions of the cosmological observables such as the luminosity distance D L in Eq. (5.39), the lensing shearγ αβ in Eq. (5.46), and the lensing rotationω in Eq. (5.47), such that there is no systematic error in missing the coordinate shift δx α o . However, the coordinate lapse δη o has significant impact [46,47], even at the linear order. The change in the observer time coordinate affects the distortion δz in the observed redshift in Eq. (3.30), because it changes the ratio of the cosmic expansion. The radial distortion δr in Eq. (3.43) is affected due to the change in the length of the null path, while the angular distortions δθ and δφ in Eq. (3.46) are not affected at the linear order. Therefore, the coordinate lapse δη o appears in the observable quantities such as the luminosity distance D L , and it was shown [46] that the absence of δη o in the luminosity distance calculations breaks the gauge invariance and the equivalence principle, causing the infrared divergences in the variance of the luminosity distance. We must emphasize that though one can set zero the coordinate lapse and shift, this choice corresponds to a gauge choice, i.e., the comoving-synchronous gauge. This means that one cannot set the lapse and shift zero and choose other gauge conditions such as the conformal Newtonian gauge than the comoving-synchronous gauge.

Photon geodesic equation and observed redshift
The photon wavevector in the conformally transformed metric trivially satisfies the geodesic equation in a homogeneous universe. In the presence of perturbations, the perturbations (δν, δn α ) to the photon wavevectork µ are constrained by the temporal and the spatial geodesic equations where we defined the derivative along the photon path with respect to the affine parameter 24) and the perturbations in the geodesic equations These perturbations were also introduced in [37,38,44], but there was a typo in Eq. (2.26) for the expression of δΓ η in [38]. Note that we already chose a rectangular coordinate for our flat FRW coordinate, in which the Christoffel symbols vanish in the background. In addition to the geodesic equation, the perturbations to the photon wavevector are subject to the null condition 0 =k µk µ : With the explicit expressions of the geodesic equations, we integrate them over the affine parameter to obtain the perturbations (δν, δn α ) λ along the photon path x µ λ : where the quantities in the square bracket are evaluated at the source and the observer positions parametrized by λ and λ o . Note that the derivative dλ along the photon path was considered only at the background level and we replaced it with the integration over the comoving distance dr, all of which are valid only when the integrands are at the linear order in perturbations. The perturbations (δν, δn α ) o at the observer position are fixed in Eq. (3.8) as the initial condition. Before we proceed to obtain the source position x µ s by integrating the geodesic equations once more over the affine parameter, we derive the expression for the observed redshift. The light emitted in the rest-frame of the source travels across the Universe, and its wavelength is stretched due to the expansion. With reference to the rest-frame wavelength or the emission frequency ω s in the source rest-frame, the observed redshift z is constructed by using the observed frequency ω o at the observer as where we used Eq. (3.10) and we defined the perturbation δz in the observed redshift. In addition to the cosmic expansion, the photon wavelength (hence the observed redshift) is affected by the peculiar velocity, the gravitational redshift, and so on, and the perturbation δz in the observed redshift captures such effects of inhomogeneities. Noting that the observer time-coordinate is η o =η o + δη o in Eq. (3.18) and the expression for the perturbation δν is Eq. (3.28), we can derive the expression for δz as (3.31) By converting to gauge-invariant variables, we isolated the gauge-dependent term in δz, and it transforms as δz = δz + HT , with which we can define a gauge-invariant variable δz χ ≡ δz + Hχ. Note that while at the observer position λ o the perturbation δz in the limit z = 0 is non-vanishing

Source position along the photon geodesic path
With the photon wavevector in the conformally transformed metric and its initial condition at the observer position, we will integrate the photon wavevector over the affine parameter to derive the source position along the photon geodesic path. The photon path is a straight line in a homogeneous universe, and the inhomogeneities in the real universe deflect the photon path from a straight line. We will begin the calculations by considering a homogeneous universe first and thereby obtaining the relation to the affine parameter set by our normalization condition in Eq. (3.6). Any position x µ λ along the photon path will be marked by the affine parameter λ, and we will use bar to indicate that the position is derived in a homogeneous universe: where we set to zero the affine parameter at the observer λ o = 0 and the position of the observer in a homogeneous universe is uniquely setx µ o ≡ (η o , 0). As a coordinate in the world-line manifold, the affine parameter is defined by the above equation as hence the spatial position becomesx The observed redshift is the only way we can assign a physically meaningful distance to cosmological objects. Since the comoving distancer is often defined in terms of a redshift parameter z, we define the affine parameter λ z and the time coordinateη z of the source in the background in terms of the observed redshift z as We have used bar for the positionx µ λ = (η,x α ) λ along the photon path to indicate that this position is evaluated in a homogeneous universe, given the observed angle n i in the observer rest-frame. We prefer to use this notation, because the position along the photon path changes in the presence of perturbations, or the mapping from a world-line manifold changes fromx µ λ to x µ λ , not because we want to introduce perturbations to a coordinate system.
With (δν, δn α ) in Eqs. (3.28) and (3.29), we can derive the expression for any position x µ λ along the path by integrating the perturbations over the affine parameter. First, we integrate the temporal part of the photon wavevector: such that the time coordinate of the source position is η s =η s + δη s . However, since the distance is more physically related to the observed redshift, we first relate the affine parameter λ s at the source position to the observed redshift as where we defined the residual deviation ∆λ s of the affine parameter λ s from λ z . Then the timecoordinate of the source position becomes η s =η z + δη s + ∆λ s ≡η z + ∆η s , (3.39) and the definition of the perturbation δz in Eq. (3.30) yields Naturally, when expressed at the observed redshift, only ∆η s will appear in the equations. 6 In the limit z = 0, ∆η s = δη o , and ∆λ s = 0. We prefer to work with (η z , ∆η s ), rather than (η s , δη s ), because only ∆η s gauge-transforms in the former combination, while (η s , δη s ) are both affected by a coordinate transformation. Note that this choice is made by convenience and we can use δη s derived in Eq. (3.37) to reproduce the relations in Eqs. (3.39) and (3.40). Compared to the time coordinateη z we infer from the observed redshift, the difference ∆η s in the source coordinate η s naturally depends on the time coordinate η o of the observer and hence the coordinate lapse δη o . However, it is independent of the rotation Ω i o of the local tetrad bases. Next, we integrate the spatial part of the photon wavevector to obtain the source position in a FRW coordinate. As mentioned, it proves convenient to express the source position x µ s around the positionx µ z inferred from the observed redshift and angle. Having computed the time distortion ∆η s in Eq. (3.40), we will compute the spatial distortion ∆x α s of the source position as Given the privileged direction n α , we decompose the spatial distortion into the radial δr distortion and the transverse distortion δx α ⊥ as where the radial distortion is [33,37,38] where we used the null condition in Eq. (3.27) and the distortion in the time coordinate in Eq. (3.40).
In [38], there was a typo in the final equation (3.5), in addition to the missing term of the coordinate shift δx α o at the observer position. The expression for δr is arranged in terms of gauge-invariant variables, isolating the gauge-dependent term n α G α s , such that δr = δr + n α L α s . Being the spatial distortion, the radial distortion δr depends on the spatial position x α o of the observer and hence the coordinate shift δx α o , but again it is independent of the rotation Ω i o of the spatial tetrad vectors (same for the distortion ∆η s in the time coordinate).
For the transverse components ∆x α ⊥ of the spatial distortions, we have to integrate δn α λ in Eq. (3.29) over the affine parameter as Further splitting the transverse distortion by using θ α and φ α in Eq. (3.13) we derive the angular distortions as [33,37,38] where we used δn α o in Eq. (3.8). For the azimuthal distortion sin θ δφ, the above equation can be used with θ α replaced with φ α . The derivative in the second integration in Eq. (3.46) cannot be pulled out as in equation (3.6) in [38] for the vector and the tensor perturbations, because it is the derivative with respect to the observed angle due to θ α : (3.47) With such expressions, we can further simplify the angular distortions as The expression for δθ is again arranged in a gauge-invariant way, thatr z δθ =r z δθ + θ α L α s . Similarly to the radial distortion, the transverse distortions depend on the observer position x α o (and the coordinate shift δx α o ). Furthermore, since they are expressed in terms of the observed angle (θ, φ) in the observer rest-frame, it is also affected by the orientation of the local tetrad basis as where for convenience we decompose the gauge-invariant vector for the rotation of the local tetrad basis as The implication is clear, such that if the local tetrad basis is rotated (Ω i o = 0) against a FRW coordinate, the angular source position would be further rotated, given the observed angle (θ, φ).
In summary, the source position x µ s , given the observed redshift z and the observed angle n i , is expressed as the sum of the positionx µ s inferred from these observables and the deviation ∆x µ s around it: because the positionx µ z inferred from the observables remains unchanged. These transformation properties are indeed satisfied in Eqs. (3.43), (3.46), and (3.40). In the limit z = 0, the source position becomes x α s = x α o = δx α o in Eq. (3.19), and hence we have δr = n α δx α o and ∆x α ⊥ = δx α ⊥o = 0. Withr z = 0 in the limit z = 0, the angular distortion (δθ, δφ) is not defined at the observer position.

Standard Weak Lensing Formalism and Its Gauge Issues
We have derived the source position in a given FRW coordinate by solving the full geodesic equation at the linear order in perturbations, including the vector and the tensor type perturbations. Here we briefly review the standard weak lensing formalism and follow the lensing formalism with the fully relativistic solution ∆x µ s in Sec. 3 to derive the expression for the gravitational lensing convergence κ, the shear γ, and the rotation ω. 7 We demonstrate that the standard weak lensing formalism is built upon the unobservable (coordinate-dependent) source position x µ s and hence it has a limitation in its relativistic generalization, which will be resolved in Sec. 5.

Short review of the standard weak lensing formalism
In classical mechanics, the gravitational interaction due to a point mass M provides a perturbation along the transverse direction to a test particle moving with the relative speed v rel : where G is the Newton's constant and b is the transverse separation (or the impact parameter). The prediction for the light deflection angleα in Einstein's general relativity is well-known to follow the same result in classical mechanics, but with additional factor two: This light deflection due to a point mass can be generalized to derive the standard weak lensing formalism by considering the gravitational potential fluctuation ψ = −GM/r due to a point mass (this indeed corresponds to the metric fluctuation α in Eq. (A.3)). The lensing potential Φ is the lineof-sight integration of the metric fluctuation, and its angular derivative gives the relation between the source angular positionŝ and the observed positionn, so called, the lens equation (see, e.g., [6,[48][49][50] for general reviews) where∇ is the angular gradient. Since the observed source position isn = (θ, φ) in a spherical coordinate, the deflection angle in the standard lensing formalism would naturally correspond to the angular distortion of the source position: Using the lens equation, the distortion matrix D (or sometimes called the amplification matrix) is defined as where I is the two-dimensional identity matrix and we defined a short hand notation for the angular derivatives of the lensing potential. The distortion matrix is conventionally decomposed into the trace, the traceless symmetric and the anti-symmetric matrices: where the trace is the gravitational lensing convergence κ and the symmetric traceless part is the lensing shear γ = γ 2 1 + γ 2 2 : Since the distortion matrix in Eq. (4.5) is symmetric, the rotation ω vanishes in the standard formalism at all orders. While the distortion matrix is defined in terms of angles, it is often assumed in literature that the line-of-sight direction is along z-axis (n/ /ẑ, i.e., θ = 0), and two angles are aligned with x-y plane. 8 In such a setting, consider two small angular vectors at the source position subtended respectively by dθ and dφ at the observer position The solid angle at the source subtended by these two angular vectors is then related to the solid angle at the observer as and hence the gravitational lensing magnification µ is then For this reason, the distortion matrix is often called the inverse magnification matrix. Using the Poisson equation in cosmology, the gravitational lensing convergence can be computed in terms of the matter density fluctuation δ m in the comoving gauge as where we used and ignored the boundary terms. The standard lensing formalism is based on the lens equation and the lensing potential in Eq. (4.3). However, the source angular positionŝ = (θ + δθ, φ + δφ) is gauge-dependent, and the lensing potential that is responsible for the angular distortion (δθ, δφ) is ill-defined. Indeed, we already know that 2ψ in Eq. (4.3) should be (α χ − ϕ χ ) to match the leading terms for δθ in Eq. (3.46) and the Poisson equation in Eq. (4.12) is indeed an Einstein equation with ψ there replaced by −ϕ χ . 9 Furthermore, there exist no contributions from the vector and the tensor perturbations in the standard lensing formalism. Finally, while the derivations in this subsection assume no linearity, all formulas of the standard lensing formalism turn out to be valid only at the linear order in perturbations.

Relativistic generalization of the standard weak lensing formalism
Despite the issues listed in Sec. 4.1, the standard lensing formalism can be readily generalized and be put in the general relativistic framework. Gravitational lensing exclusively deals with the distortion in the source angular positionŝ = (θ + δθ, φ + δφ), compared to the observed angular position n = (θ, φ). So we can generalize the standard lensing formalism by replacing the deflection angle −∇Φ due to the lensing potential by the angular distortion (δθ, δφ) with all the relativistic effects taken into account. However, already apparent in the calculations in Sec. 3, the source position x µ s is unobservable and gauge-dependent. Derived based on the unobservable quantities, the distortion matrix and the lensing observables in the generalized lensing formalism will be gauge-dependent. Here we prove this by explicitly computing the distortion matrix and the lensing observables with the angular distortion (δθ, δφ) in Eq. (3.46).
First, let us consider the "generalized" lens equation Despite this incompleteness, we proceed with this simple relativistic generalization for now. Consider an extended source in the sky that appears subtended by an infinitesimal size (dθ, dφ) in angle centered at the observed angle (θ, φ). Since the angular position of the source isŝ = (θ + δθ, φ + δφ) in a spherical coordinate, the infinitesimal size of the source is then where we computed the variation of the source position and ignored the variation along the radial direction n i . Note that these angular sizes (dθ, dφ) should not be confused with the angular distortion (δθ, δφ). Expressing the angular extent of the source along the observed angular directions θ i and φ i , the distortion matrix can be derived as In this way, we generalized the standard weak lensing formalism by using the angular distortion (δθ, δφ) computed from the geodesic equation. The distortion matrix computed in Eq. (4.18) is exact for an infinitesimal size of the source, provided that the angular distortion (δθ, δφ) is also exact (computed at all orders in perturbations). However, we defined the angular extent (ds θ , ds φ ) of the source along the observed angular directions, rather than the angular directions set by the source position s i = (θ+δθ, φ+δφ), i.e., the left-hand side of Eq. (4.18) should be (ds θ , sin s θ ds φ ). This is indeed the ambiguity in the standard weak lensing formalism, and we show that the distortion matrix defined this way reproduces the standard results for the convergence and the shear. However, it is clear that the distortion matrix is not physically well-defined in the standard weak lensing formalism. Given the distortion matrix in Eq. (4.18) and the angular distortion (δθ, δφ) in Eq. (3.46), it is straightforward to compute the lensing observables. First, the gravitational lensing convergence is where all the terms are arranged in terms of gauge-invariant variables except two terms multiplied with G α , such that the gravitational lensing convergence gauge-transforms as There exist typos for the expression κ in [38], in addition to the spatial shift δx α o of the observer position. The absence of the rotation Ω i o in Eq. (4.19) represents that the gravitational lensing convergence is independent of rotation of the local tetrad basis. This is indeed true for the rotation along the line-of-sight direction, but not for the other directions beyond the linear order in perturbations.
Most works on weak lensing in literature consider only the scalar contribution and adopts the conformal Newtonian gauge, in which G α ≡ 0. While there exist the vector and the tensor contributions to κ, one can safely assume that the initial conditions are devoid of any vector or tensor contributions at the linear order in perturbations. Under such assumptions, the gravitational lensing convergence is The first term is the velocity contribution at the observer, and the second term is the contribution due to the coordinate shift of the observer position. These two contributions are often missing in literature. While the first term is sometimes considered in the calculations (e.g., [40]), the second contribution was never considered in literature. Since the conformal Newtonian gauge leaves no gauge freedom, the expression of κ cN is gauge-invariant. 10 However, as expected in Eq. (4.19), the lensing convergence κ (or κ cN ) is not the correct lensing observable we measure in surveys. The next lensing observable is the gravitational lensing shear, or the traceless symmetric part of the distortion matrix. According to the decomposition of the distortion matrix in Eq. (4.6), two independent shear components can be computed as The expressions for the gravitational lensing shear are more complicated due to their rotational properties. Since the shear components γ 1 and γ 2 are confined to a two-dimensional plane perpendicular to the line-of-sight direction, they transform under rotation around the line-of-sight direction, but the shear components are invariant under 180 degree rotation. This implies that they can be expressed in terms of the spin ±2 objects ±2 γ as where a function s f (θ, φ) with spin s transforms under a clock-wise rotation of axes by Ξ as (see, e.g., [51][52][53][54]) Using the basis vectors m i ± with spin ∓1 under rotation that describe objects with spin ±1 we can construct the basis matrices with spin ∓2 and their variants where we used the internal indices for their expressions, as they are defined in the observer restframe. However, as in Eq. (3.13), we will project the internal indices into a FRW coordinate with the background relation δ α i . The spin ±2 shear components ±2 γ can then be expressed as where the shear matrix is [40] It is apparent that the shear matrix γ αβ (or the shear components γ 1 , γ 2 ) is gauge-dependent, and they transform as The shear matrix is independent of the rotation Ω i o , because it only contributes to the second-order shear matrix. The expression for the shear matrix in Eq. (4.31) is again incomplete and gaugedependent. However, when the conformal Newtonian gauge is adopted and only the scalar perturbations are considered, this expression corresponds to the correct expressionγ αβ in Eq. (5.46). However, when the tensor modes are considered, the shear matrix γ αβ contains the tensor contribution at the observer, but not at the source position. This missing contribution is called the FNC term [55] (in relation to the metric shear [56]), and the absence of the FNC term at the source position also indicates that Eq. (4.31) is incomplete.
One remaining component of the lensing observables is the rotation ω. The anti-symmetric part of the distortion matrix in Eq. (4.18) can be obtained in a similar way as The rotation naturally changes in proportion to Ω n o , or the rotation of the local tetrad basis along the line-of-sight direction. In fact, the physical rotation should be against the local tetrad basis, absorbing the rotation Ω n o of the local basis. The expression for the rotation is gauge-dependent due to the presence of G α , and it transforms as 2ω = 2ω + n · ∇ × L . (4.34) Furthermore, the rotation is non-vanishing in the presence of the vector and the tensor perturbations, while it vanishes in their absence. This is also known as the Skrotsky effect [57] or the gravitational Faraday effect. However, this contribution to the rotation along the line-of-sight is artificial, arising from the rotation of the tetrad basis against the global coordinate. Indeed the tetrad basis vectors are parallel transported along the geodesic, and hence its contribution to the rotation vanishes, as the lensing images are always compared against the local tetrad basis. We come back to this issue in Sec. 5.2.

Gauge-Invariant Formalism of Weak Lensing
Here we present the full gauge-invariant formalism of cosmological weak lensing and resolve the issues associated with the lensing rotation.

Geometric approach to weak lensing
In Sec. 3 we have derived the source position x µ s in a FRW coordinate, given the observed redshift z and the angular position n i in the rest-frame of the observer. The source position is different from the inferred positionx µ z , and this difference is geometrically decomposed as the radial distortion δr, the angular distortion (δθ, δφ), and the time distortion ∆η s of the source position. Here we extend the standard weak lensing formalism by applying the geometric approach and checking the gaugeinvariance of the lensing observables.
In order to develop a gauge-invariant formalism of weak lensing, we need to separate observable quantities and unobservable quantities in weak lensing. The first and the foremost is the source angular position s i that is not observable. While we measure the source position n i in the rest-frame, the source position x µ s is constructed in a FRW coordinate by tracing the photon path backward in time, and its geometric deviations from the inferred position are gauge-dependent as derived in Sec. 3. So the source angular position in a FRW coordinate can be legitimately expressed as Consider again an extended source in the sky of the observer that appears subtended by an infinitesimal angular size (dθ, dφ) as in Sec. 4. The angular size of the source in a FRW coordinate is then +φ α sin θ dφ 1 + cot θ δθ + ∂ ∂φ δφ + dθ ∂ ∂θ (sin θ δφ) − n α δθ dθ + sin 2 θ δφ dφ .
Compared to Eq. (4.17), we notice the presence of the additional term along the line-of-sight. Furthermore, the source position in a FRW coordinate is indeed specified in terms of its four-dimensional position: η s =η z + ∆η s , x α s =x α z + δr n α + ∆x α ⊥ =r z s α + δr n α , (5.4) Therefore, the source size in a FRW coordinate that would appear subtended by angular size (dθ, dφ) at the fixed observed redshift z is then dx α s =r z ds α + d (δr) n α + δr dn α ≡r z (dn α + ∆s α ) , where the first term represents the source size in the absence of any perturbations and the second term with ∆s α is defined to capture any deviations due to the perturbations We will not need to consider the variation in the time coordinate ∆η s at the linear order. Compared to Eq. (4.17), the radial distortion δr also contributes to the source size ∆s α in a FRW coordinate. As we further proceed, it will become clear that the radial component in proportion to n α will be projected out, so that ignoring this part in Eq. (4.17) causes no systematic errors at the linear order in perturbations. However, note the presence of the radial distortion δr in proportion to θ α and φ α , compared to Eq. (4.17).
In Section 2, we set up the tetrad basis at the observer position to relate the local observables to the photon wavevector k µ in a FRW coordinate. By solving the geodesic equation in Section 3, the photon wavevector is traced back to the source position. To properly relate the source size dx µ s in a FRW coordinate to the physical size dx i s of the source in its rest-frame, we need to construct the tetrad basis in the source rest-frame. The expressions for the tetrad basis vectors at the source position are almost identical to those at the observer position in Eqs. (2.4) and (2.13), except the fact that all the quantities are evaluated at the source position and the velocity vector U α represents the velocity of the source. Finally, we need to define a two-dimensional plane in the rest-frame of the source, perpendicular to the light propagation. A small area in this plane such as a circle or an ellipse will appear distorted at the observer, and its observed shape compared to that in the source rest-frame will determine the physical lensing observables such as the magnification (convergence), the shear, and the rotation that we can measure First, the photon propagation direction measured by the source would be expressed in a FRW coordinate as where we defined the perturbation to the photon propagation at the source position ∆n α s ≡ δn α s + U α s − ∆ν s n α .
Again, the photon wavevector k µ at the source position is fully determined (up to Ω i o ) given the observation at the observer position. Once the source velocity is specified as an "observer" (in fact, an "emitter"), the photon propagation direction n µ s (or ∆n α s ) measured by this observer is then fully determined, as indicated above.
Next, this photon propagation direction n µ s in a FRW coordinate needs to be projected into the rest-frame of the source to define a two-dimensional area perpendicular to the light propagation. In the source rest-frame, the "observed angle n i s " measured by an "observer" at the source position is then n i s = g µν e µ i n ν s g αβ e α i n β + O(2) = n i + 2C αβ n β δ αi + ∆n α δ i α + n α δe α j δ ji (5.10) where the integration represents the perturbation contribution of the photon propagation along the line-of-sight, the terms at the source and the observer positions are due to the misalignment of their rest-frames and FRW coordinates, and the expression is fully gauge-invariant. This equation is indeed the correct lens equation, which generalizes over Eqs. (4.3) and (4.15).
The observed angle n i s at the source (or the emission angle) in the rest-frame is identical to the observed angle n i at the observer position in the background due to the construction of our local tetrad basis. Even in the absence of any perturbations, these two angles can be different, simply by setting up the local coordinates differently. In perturbation theory, all these differences are, however, shifted to the perturbations, while the background remains unaffected, and indeed the difference in the observed angles is proportional to the difference in rotation of the local tetrad bases at the source and the observer positions: At the moment, the rotations Ω i of the local tetrad bases at both positions are undetermined. By parameterizing the source angle with perturbations (∆θ, ∆φ) in the source rest-frame as we can express the light propagation direction in the source rest-frame and set up two orthonormal bases θ i s and φ i s that are perpendicular to the propagation direction and form an orthonormal basis where we defined the perturbation vectors ∆θ i s and ∆φ i s for the notational simplicity and note that (θ, φ) in Eq. (5.12) is the same observed angle at the observer position. Given the explicit expression for n i s in Eq. (5.10), we can compute the perturbations to the source angle (θ s , φ s ) as and the expression for sin θ ∆φ is equivalent with θ i replaced by φ i . These perturbations in angle are gauge-invariant, and they contain the rotation of the local tetrad bases Having obtained two orthonormal basis vectors, we are in a position to compute the physical size in the source rest-frame that is expressed as dx µ s in a FRW coordinate and would be measured at the observer position with angular size (dθ, dφ). In reality, the source position would be measured as the observed angle (θ, φ) at the observed redshift z, and the physical size of a standard ruler in the source rest-frame would be subtended by its angular size (dθ, dφ), as opposed to our reverse construction of the standard ruler. So by measuring the angular size with prior knowledge of the standard ruler, we can construct the lensing observables, and compared to the standard lensing formalism, the lensing observables obtained this way are physically well-defined. The source size dx µ s in a FRW coordinate can be projected into the rest-frame along the two orthogonal directions in the plane perpendicular to the light propagation direction: The changes in the observed angles at the source position due to the distortion of the local frame contribute to the physical sizes dn α ∆θ α s + δe α i θ i = dθ δe θ θ + sin θ dφ δe φ θ + cos θ ∆φ , (5.21) where we used the short-hand notation δe θ φ ≡ θ α δe α j φ j for instance. Putting it altogether, the physical sizes in Eq. (5.18) can be simplified as where we used the same short-hand notation for C αβ and noted that δe θ θ = −C θθ and δe φ φ = −C φφ . Noting that the angular diameter distance in the background isD A = a zrz and the scale factor a s at the source is related to the perturbation δz in the observed redshift we can define the distortion matrixD (with hat) as dθ sin θ dφ , (5.25) and the elements of the matrix can be read off aŝ δθ − cos θ δφ + 2 C θφ + δe φ θ + cos θ ∆φ , (5.28) where the elements of the distortion matrix D (without hat) in Eq. (4.18) are put in the parenthesis to facilitate the comparison. There exist notable differences in Eq. (5.25) that are physically unambiguous, when compared to the distortion matrix D: The source position in the sky appears to the observer as the observed angle (θ, φ) at the observed distanceD A (z) (or at the observed redshift z). Given the knowledge of the physical sizes and the orientation (dL θs , dL φs ) of a standard ruler, their "observed angular size" (ds θ , ds φ ) in the absence of perturbations would be dL θs dL φs ≡D A ds θ sin θ ds φ , (5.30) and hence the distortion matrix defined above is indeed the lensing distortion matrix, providing the relation between the angular sizes at the source and at the observer: Note that their observed angular position is (θ, φ), not (θ + δθ, φ + δφ), naturally resolving the ambiguity in Eq. (4.18).
To ensure these expressions are correct and physically well-defined, we explicitly verify the gauge-invariance of the distortion matrix. Any physically well-defined quantities should be gaugeinvariant at the linear order [41]. For a general coordinate transformation in Eq. (A.4), we already derived in Sec. 3 how the geometric distortions (δr, δθ, δφ) of the source position transform. Given that ∆φ in Eq. (5.16) is gauge-invariant, it is relatively straightforward to show that each component of the distortion matrixD is indeed gauge-invariant, by using how the remaining components transform Having verified the gauge-invariance of the distortion matrix, we will decompose the distortion matrix and derive the lensing observables. First, we consider the convergence of the distortion matrix. Consider two physical separation vectors dL i dθ and dL i dφ in the source rest-frame that would appear to the observer subtended by dθ and dφ, respectively. Using the distortion matrix in Eq. (5.25), we can derive these two separation vectors in the source rest-frame: Note that we have two separation vectors and two angular vectors in the source rest-frame. The physical area dA spanned by these two separation vectors in the source rest-frame is where the relation is non-perturbative, given the distortion matrix. This relation is indeed the definition for the angular diameter distance D A : where δD is the fluctuation in the angular diameter distance. 11 By expressing the physical area in terms of the observed angles in the absence of perturbations (or the source angles), dA s = ijk n i s dL j θs dL k φs =D 2 A dΩ s , dΩ s = sin θ ds dθ θ ds dφ φ − ds dφ θ ds dθ φ , (5.37) 11 The fluctuation δD in the angular diameter distance is identical to the fluctuation in the luminosity distance, because the luminosity distance is DL = (1 + z) 2 DA. the determinant of the distortion matrix is indeed related to the observed magnification µ as At the linear order in perturbations, the determinant of the distortion matrix is where the convergence κ in the standard lensing formalism is given in Eq. (4.19) and we used The convergenceκ (with hat) of the distortion matrixD is then Therefore, the diagonal componentκ (with hat) of the lensing observables is related to the physical magnification (or the luminosity distance), and it is indeed gauge-invariant, as opposed to the gauge-dependent convergence κ (without hat) in Eq. (4.19). The expression is independent of the rotations Ω i at the source and the observer positions, because it only involves the ratio of the physical and the inferred areas. The detailed properties of the gauge-invariant luminosity distance fluctuation are studied in [46,47]. Next, we compute the gravitational lensing shear of the distortion matrixD. Two components of the lensing shear can be readily read off from the distortion matrix as where the shear components γ 1 and γ 2 in Eqs. (4.22) and (4.23) are shown in the parenthesis. The additional contributions in the shear componentsγ 1 andγ 2 arise due to the change of a FRW frame to the source rest-frame: which exactly cancel the gauge-dependent terms in γ 1 and γ 2 and add the tensor contribution C αβ at the source position (or the FNC term at the source [55]). In the same way in Sec. 4, we can construct the spin ±2 shear components ±2γ =γ 1 ± iγ 2 and derive the shear matrix Compared to Eq. (4.31), the shear matrixγ αβ (with hat) is gauge-invariant, and there exist additional tensor contribution at the source position, or the FNC term [55]. This contribution arises because the physical length is defined in the source rest-frame, not in the FRW frame. Indeed, its presence is necessary to prevent the infrared-divergences of the shear matrix (see Sec. 6.1). Finally, we derive the last remaining lensing observable of the distortion matrixD, or the rotation (with hat) where the rotation ω (without hat) from D in Eq. (4.33) is shown in the parenthesis. The additional contribution to the rotationω from the frame change is which cancels the gauge-dependent part in ω (without hat). Compared to Eq. (4.33), the rotationω has extra terms that involves the rotation at the source position, i.e., Ω n s and ∆φ. Any physical rotation should be measured against the local tetrad bases, and hence it naturally involves the difference (Ω n o − Ω n s ) in orientation. In short, the standard lensing formalism lacks the specification of the source rest-frame, in which the physical size and shape are defined. By setting up a local tetrad basis at the source position, we fixed the standard weak lensing formalism and explicitly verified the gauge-invariance of the lensing observables.

Physical rotation of the images vs rotation of the tetrad basis
We have computed the rotationω as one of the lensing observables in Sec. 5.1 by comparing the orientation of the physical lengths L θs and L φs along two orthogonal directions in the source rest-frame to the orientation they appear in the observer rest-frame. However, as apparent in Eq. (5.47), the rotationω depends on the orientations Ω s and Ω o of the local rest-frames of the source and the observer around the light propagation direction. Technically, if we change Ω o , the observed angle (θ, φ) changes or rotates. But in reality, the fact that we assign some numbers to the observed light propagation direction (θ, φ) means that we already fixed the orientation Ω o of our local coordinate. However, the orientation Ω s of the source remains unspecified, and the lensing observableω depends on this unspecified orientation. In other words, by rotating the orientation of the local coordinate in the source rest-frame, we can set the rotation to zero. The question then naturally arises, "is the rotation physical and measurable?" Let's revisit how the lensing observables are defined physically. The lensing convergenceκ (or the determinant at the linear order) can be defined as the ratio of the physical area in the source rest-frame to the area inferred from the observed angular size. So it is independent of the orientations at the source and the observer positions. The lensing shear can be defined as the observed ellipticity of a circular object in the source rest-frame. Hence the shear is again independent of the orientation of the source, while it depends on the orientation of the observer coordinate. Though Eq. (5.46) is independent of Ω o at the linear order, the shear matrix rotates in general as we rotate the local coordinate, and two shear componentsγ 1 andγ 2 represent the shear amplitude and orientation (the latter of which depends on the orientation of the observer coordinate). To define the rotation in a physical way, we need to synchronize the orientations at the source and the observer positions, from which any deviation can be assigned to rotation. In a curved spacetime, this can be achieved by parallel transporting the local tetrad basis along the photon propagation direction. Indeed, this is the only physically meaningful way, as the parallel transport is path-dependent.
One subtlety associated with this procedure is that when the local tetrad vectors e µ a are parallel transported to the source position, the transported timelike vectorȇ µ 0 (with check) at the source position is different from the source velocity u µ s . To construct the parallel transported basisê µ a in the source rest-frame, we need to Lorentz boost the transported basisȇ µ a with the source velocity vector u µ s . Indeed, one can show [42] that any vectors A i ⊥ and B i ⊥ defined in a plane perpendicular to a photon propagation direction n i of a null vector in Eq. (2.1) satisfy where those with hat represent the components after the Lorentz boost (see, also, [58]). In our case, the photon propagation directionn i constructed from n µ s with the transported timelike vectorȇ µ 0 changes ton i due to the Lorentz boost, and indeed we already derived the photon propagation directionn i = n i s in Eq. (5.10). For our purposes, the plane at the observer is defined in terms of two orthonormal vectors θ µφ µ (or θ i -φ i in the rest-frame), and the transported plane spanned byθ µ s -φ µ s defines an "oriented" plane perpendicular to the photon propagation direction n µ s but with the transported velocityȇ µ 0 . Equation (5.49) implies that the physical size and the shapeÂ i ⊥ defined in a plane perpendicular to the photon propagation direction n µ s with the correct source velocity u µ s is equivalent to those defined in a plane perpendicular to n µ s withȇ µ 0 . Therefore, the source size dx µ that appears subtended by the angular size (dθ, dφ) can be projected into any of the two planes for computing the lensing observables. This argument is indeed crucial to the Jacobi mapping method in Sec. 6.3 that is often implicitly assumed (see, e.g., [39] for the discussion).
Given the photon propagation direction n µ s with the correct source velocity u µ s , we computed the photon propagation direction n i s in Eq. (5.10) and the plane θ i s -φ i s perpendicular to n i s . To compute the physical rotation, first we need to parallel transport the local tetrad basis vectors θ µ and φ µ to the source position, then the transported vectorsθ µ s andφ µ s need to be Lorentz boosted to the source rest-frame with u µ s . Finally, we need to subtract the angle Θ from the rotationω in Eq. (5.47) that the Lorentz boosted vectorsθ µ s andφ µ s make against θ µ s and φ µ s , respectively (note that the plane spanned byθ µ s -φ µ s is identical to that by θ µ s -φ µ s , but their individual directions are not aligned). This procedure is illustrated in Figure 1.
Fortunately, the second step is trivial at the linear order in perturbations, as the Lorentz boost by U α s of the source is already at the linear order. The spatial tetrad vectors after the Lorentz boost areê identical in their components at the linear order. Therefore, the directional vectors ofθ µ s andφ µ s in the source rest-frame given the tetrad basis e µ i arê where a particular attention needs to be paid to the difference in notation of ∆θ i s and ∆φ i s (without hat) in Eq. (5.14). So the angle Θ between two directional vectors is The rotation of the local coordinate in the source rest-frame by the angle Θ will affect the distortion matrixD in Eq. (5.25) as R(Θ)D, where R is the rotation matrix. This needs to be contrasted by In particular, the two-dimensional plane spanned byθ µ s -φ µ s needs to be Lorentz boosted to the source rest-frame to form the correct and oriented plane denoted byθ µ s -φ µ s . This plane is identical to the plane spanned by θ µ s -φ µ s , but two sets of the basis vectors make an angle Θ.
the usual rotation of coordinates, where the distortion matrix changes as R(Θ)DR t (Θ). The latter arises when we rotate the local coordinates both at the source and the observer positions by the same angle Θ, while the former in our case rotates only the local coordinate at the source position. Given this change, only the rotation among the lensing observables is affected aŝ ω →ω + Θ , (5.54) at the linear order in perturbations, while the convergenceκ and the shearγ αβ remain unaffected. Therefore, our task boils down to computing the transported vectorsθ µ s andφ µ s to derive the angle Θ. The condition for the parallel transport of the tetrad basis vectors is where we used the conformal transformation relation in Eq. (3.5). At the background, we can readily check that two spatial tetrad vectors perpendicular to the light propagation are trivially transported as while the propagation direction vector and the four velocity are parallel transported in a non-trivial way, even at the background, i.e., the simple parametrizations e µ n = 0, do not satisfy the parallel transport condition. For our purposes, we will focus on the two spatial tetrad vectorsθ µ andφ µ . Integrating the parallel transport condition, we can derive the perturbation δȇ µ θ of the spatial tetrad vectors parametrized in Eq. (2.5) as where the initial condition (δe α θ ) o is set by Eq. (2.13). Therefore, the directional vectorsθ i s andφ i s that are parallel transported and Lorentz boosted arê where the last terms δe θ,φ i in both expressions are the perturbations of the spatial tetrad vectors at the source position we used for computing the lensing observables, not the parallel transported tetrad vectorsȇ i θ,φ . So, the angle between two basis vectors is then where we simplified it with an integration by part and the boundary terms are The rotation angle Θ is the rotation of the local tetrad basis vectors (θ µ -φ µ ), as we parallel transport them to the source position, and it is indeed the rotation angleω in Eq. (5.47). However, we should emphasize that the rotation of the tetrad basis is a "rotation" against an arbitrary global FRW coordinate and there is no way to compare the orientations at two different spacetime points without parallel transporting the basis. In this regard, the so-called Skrotsky effect is not physical. Therefore, by setting the local coordinate at the source position to that of the observer parallel transported and Lorentz boosted, the misalignment angle Θ is zero, and the orientation of the local coordinate at the source position is completely determined as Putting this in Eq. (5.47), we find that the physical lensing rotation iŝ with no further arbitrariness left. We want to emphasize that this is the only physically meaningful way to define the rotation and the physical rotationω vanishes only at the linear order. In other words, the lensing rotation is not a coordinate artifact (but the Skrotsky effect is), and it is non-zero beyond the linear order in perturbations.

E-B decomposition and its relation to the rotation
The traceless symmetric shear tensor γ αβ (or γ 1 and γ 2 ) was decomposed in terms of the helicity eigenstates ±2 γ of spin s = ±2. However, as introduced in [59], the shear tensor can be decomposed in terms of two potentials Φ E and Φ B (e.g., see [28,31,60,61]). Following the approach developed in [60], the lens equation can be generically written as 64) where ij is the Levi-Civita in two-dimensional angular basis, the E-mode potential Φ E represents the usual lensing potential, and the B-mode potential Φ B represents the pseudo vector of the deflection angle. 12 Given the expression for the angular position of the source, the distortion matrix in Eq. (4.18) can be derived as in Sec. 4.2: and the lensing observables can be directly read off from the distortion matrix as: These relations of the lensing observables to the E-and B-mode potentials can be compactly represented as [28,31,60,61] 75) 12 Here our notation convention differs in the following quantities with minus sign from those quantities defined in [60]: The shear matrix in [28,31] is indeed the ellipticity matrix, and hence their shear matrix corresponds to 2γ αβ in our notation, giving rise to a factor 2 difference in the E-B decomposition.
where the notation for the two potentials γ E and γ B in [28] is γ E = ∇ 2 Φ E /2 and γ B = ∇ 2 Φ B /2. The spin ±2 shear components ±2 γ are also related to the two potentials as The two potentials are uniquely defined up to the transformation that supplies the same lensing observables under the conditions: The lens equation written in terms of unit angular vector has two degrees of freedom, and these two degrees of freedom are captured by two independent potentials Φ E and Φ B . As derived above, all the lensing observables are therefore expressed in terms of Φ E and Φ B : The convergence κ is exclusively described by Φ E , and the rotation ω is also exclusively by Φ B , while the shear matrix γ ij receives the contributions of both Φ E and Φ B . To compute the two potentials, we need to use the lens equation and the lensing observables. From the lens equation (4.3) in the standard lensing formalism, the two potentials can be read off as where the E-mode potential Φ E is the projected lensing potential, sourced by the density fluctuations, while the B-mode potential is absent Φ B = 0 to all orders. In Sec. 4.2 we generalized the standard weak lensing formalism by solving the geodesic equation and accounting for the contributions from the scalar, vector, and tensor perturbations, and we derived the lens equation (4.15) and the lensing observables in Eqs. (4.19), (4.22), (4.23), and (4.33). The two potentials in the generalized lensing formalism can then be derived as where∇ −2 is the inverse (angular) Laplacian operator and we used the useful relations (and their inverse relations)∇ the source positions. The B-mode potential Φ B is also excited in the presence of the vector and the tensor perturbations, giving rise to the lensing rotation 2ω = −∇ 2 Φ B in Eq. (4.33). However, as emphasized in Sec. 4.2, all the lensing observables in this simple generalization of the standard formalism are gauge-dependent, and consequently the two potentials derived above are also gaugedependent: In summary, while the E-B decomposition is a useful tool, the physical lensing observables (κ, γ ij ,ω, orD) are not fully described by Φ E and Φ B alone, and hence the relation in Eqs. (5.70)−(5.73) breaks down for the physical lensing observables.

Comparison to Previous Work
Here we compare our gauge-invariant lensing formalism to previous work in literature. While there exists an extensive work in lensing, there are relatively few papers that treat the weak lensing observables in a fully relativistic framework. In most cases, previous work adopts the conformal Newtonian gauge and computes only the scalar contributions. Few papers considered the vector and the tensor contributions, and even fewer papers checked the gauge-invariance of the lensing observables. In this section, we first compare our results in the presence of tensor modes, highlighting the need for the tensor contributions at the source and the observer positions. Then we discuss the Cosmic Ruler papers [40,55] that computed the lensing observables with a different approach and provided comprehensive accounts of the lensing observables. While in literature there exist several work [62][63][64][65][66][67] that computed the lensing observables beyond the linear order, we will focus on our comparison to those only at the linear order in perturbations.

Contributions of the primordial gravity waves to the lensing observables
The primordial gravity waves in the early Universe provide the most natural way to excite the tensor modes on large scales at the linear order in perturbations. The presence of the tensor perturbations C αβ affects the lensing observables in many subtle ways. The photon path is deflected (in Sec. 3), as it propagates along the line-of-sight direction; The rest-frames of the source and the observer are affected (in Sec. 5), changing the relation between the photon wavevector in the FRW frame and the observables in the observer rest-frame. These tensor contributions to gravitational lensing have been considered in the past [40,55,56,68]. Here we briefly compare our calculations to the previous work.
Dodelson, Rozo, & Stebbins [56] considered the tensor contributions to the lensing observables. Focusing only on the tensor modes, they solve the geodesic equation, and derive the light deflection due to the tensor modes. Upon identifying their notation the source position r in their Equation (3) can be arranged to match the angular distortion, when only the transverse directions are considered: Compared to our angular distortion in Eq. (3.48), is exactly the rotation ω (without hat) in Eq. (4.33), except the rotation Ω i o = 0 and the minus sign due to the different sign convention. However, as we showed, this rotation ω is a coordinate artifact due to a rotation of the tetrad basis along the light propagation, when compared against a global coordinate. Interestingly, by imposing the consistency condition among the lensing observables, they introduced the tensor contribution at the source position (so called the metric shear), which is what we need to fix the gauge-dependent shear matrix γ αβ when we transform from the FRW frame to the source rest-frame.
Recently, Adamek, Durrer & Tansella [68] considered the tensor contributions h ij from the primordial gravitational waves to the lensing observables, and they used the relativistic cosmological simulations [69,70] to compare the gravitational wave signals in the lensing observables to the tensor contributions generated by the nonlinear evolution of large scale structure in the late time. They solved the geodesic equation and derived the deflection angle α a , which corresponds to our angular distortion α a → (δθ, sin θ δφ) , h ij → 2C αβ . (6.7) Their Equation (2.5) is equivalent to our Eq. (3.46), but without the tensor contribution at the observer position. The absence of the tensor contribution causes the infrared divergences in the power spectrum, which they regulate by introducing a counter term. These tensor contributions at the source and the observer positions are not constrained by the geodesic equation, as they arise in defining the rest-frames. However, their significance is apparent in both cases [56,68].
In application to CMB polarization, Dai [71] computes the rotation of CMB polarization and argues that it is maximally correlated with the lensing rotation, affecting the CMB polarization power spectrum. However, the lensing rotation by the vector and the tensor vanishes at the linear order, and the rotation in [71] is a coordinate artifact. Furthermore, since polarization is parallel transported, the rotation of polarization is zero to all orders, and its correlation to the lensing rotation beyond the linear order is therefore zero.

Cosmic Ruler approach: lensing by a standard ruler
It is well known from the early days in lensing literature that the geodesic equation and its infinitesimal deviation vectors provide the essential ingredients for describing the lensing phenomena (see, e.g., [18,19]). However, as emphasized throughout this work, the photon wavevector obtained by solving the geodesic equation needs to be related to the observables in the observer rest-frame and the physical size and shapes in the source rest-frame. This requires setting up the local tetrad bases at the observer and the source positions. In this regard, a comprehensive and pioneering work has been done [40,55] under the name of the "Cosmic Ruler," with which we briefly compare our calculations.
The cosmic ruler approach assumes a standard ruler in the rest-frame of the source and computes its relation to the observables. Imagine a stick of a known length L (or the standard ruler) and measure the light from this stick. In particular, we measure the angular positions n i 1 and n i 2 and the observed redshifts z 1 and z 2 of the two end points of the stick: where two separation vectors are defined as In the limit the size of the ruler goes to zero (small ruler), we consider a case, in which two end points of the ruler are at the same observed redshift (z 1 = z 2 ), but with slightly different angular positions (n i 2 = n i 1 + dn i ; transverse ruler). We will refer to this limit as the "small transverse ruler limit for lensing," which is relevant for our purposes. 13 In this limit, these separation vectors correspond to Given these positions in a FRW coordinate, the (known) length of the standard ruler can be computed by projecting two end points to the rest-frame of the source as where the projection tensor for the source position x µ 1 is This projection tensor ensures that when contracted, any four vectors are projected in the rest-frame of the source. Since the projection and any other operations in this case involve two end points with two different redshifts, there always exists ambiguity regarding the choice of the operation point, i.e., either at x µ 1 or x µ 2 . The observed size L z of the standard ruler will be inferred based on the observed angular size and the redshift (again we use z 1 ) as where the orientation of the ruler is further decomposed [40] along the line-of-sight direction δx and the transverse direction δx ⊥ (by using n 1 ). Using Eq. (6.12), we can compute the size of the standard ruler at the linear order in perturbations and relate it to the observed size L z as where δη = −δx is assumed [40] and the deviation vector ∆x β s is expanded as Having derived the relation between the physical size L and the observed size L z and measured six observables (two redshifts and two angular positions), we can determine six physical quantities associated with the distortion of the standard ruler as [40] 1 where we defined the inferred comoving size L 2 c ≡ δx 2 + δx 2 ⊥ . The scalar C and the two-component transverse vector B α are related to the radial distortion. For our purposes of deriving the lensing observables, we take the small transverse ruler limit (δx = 0), and the symmetric matrix A αβ plays the equivalent role of the distortion matrixD. The symmetric matrix can be further decomposed in terms of the trace M and the spin-2 shear γ 1 & γ 2 as where the components of the matrix is written in a Cartesian coordinate aligned with (θ, φ, n) and the magnification scalar M and the shear matrix γ ij are Identifying their notation ∆x → δr , ∆x i ⊥ →r z δθ θ i + sin θ δφ φ i , (6.20) we can derive the relation of the ruler observables to our lensing observables as In their calculations [40], the photon wavevectork µ is related to the observed angle and redshift with an implicit assumption ∆ν o = 0 in our calculation. Furthermore, they assume that the observer position is x µ o =x µ o , i.e., the coordinate lapse δη o and the spatial shift δx α o are ignored. At the linear order in perturbations, the magnification M and the shear matrix γ ij are independent of the spatial shift δx α o , but the magnification depends on the coordinate lapse δη o . Its absence is the source for the infrared divergence in the variance of the luminosity distance (see [46]). 14 The tetrad basis at the source position is essential to our approach in relating the photon wavevector to the physical size and the shape in the source rest-frame. This part was replaced in the cosmic ruler by the projection tensor H µν in Eq. (6.12). In fact, the projection operation has to be performed not only to the restframe of the source, but also to a plane perpendicular to the photon propagation. So, using the photon direction n µ in Eq. (5.7) at the source position, the correct projection tensor can be derived as However, at the linear order, the additional term from n µ n ν is always multiplied by the transverse deviation δx µ ⊥ , and hence these contributions vanish asn α s ∝ n α . Furthermore, in the cosmic ruler, the orientation of the standard ruler is not specified, and hence the rotation is absent in A αβ by construction. While the rotationω at the linear order indeed vanishes, this is not true beyond the linear order. A straightforward extension of the cosmic ruler to include the lensing rotation would be to consider two (distinguishable) standard rulers in the source rest-frame, e.g., a red stick and a blue stick with known sizes and orientations in the source rest-frame. Despite these points for improvement, the cosmic ruler approach [40,55] provides a very comprehensive account of the lensing observables with an explicit check of gauge-invariance and the consideration of the source rest-frame, in which the standard ruler is defined. Moreover, their calculations of the magnification and the shear matrix are fully consistent with our calculations, and their work has been further extended [55,72] to include the intrinsic shear in a gauge-invariant way, which is beyond the scope of our current work.

Jacobi mapping approach
The Jacobi mapping approach provides a different way to formulate the lensing observables (see, e.g., [26,39,58,66,67]). Consider an infinitesimal field ξ µ (or the Jacobi field) that connects two photon paths parametrized by the same affine parameter. 15 The Jacobi field can be expressed in terms of the local tetrad basis parallel transported along the photon path as whereȇ µ A (A = θ, φ) forms the two-dimensional plane perpendicular to the light propagation and e µ 0 is the parallel transported velocity u µ λ . This basis is often referred to as the Sachs basis [17] and is more convenient than any other tetrad basisȇ µ a that are parallel transported, because its spatial 14 However, see their arXiv version (1204.3625v3) of the Cosmic Ruler paper [40], where they included the coordinate lapse δηo in the expressions, exactly to prevent the infrared divergence of the theoretical predictions. The spatial shift δx α o is not included, but again its contribution drops out in the lensing observables. 15 Only in this subsection, we use ξ µ to denote the Jacobi field. This should not be confused with an infinitesimal coordinate transformation in Eq. (A.4). basis vectors are aligned with the light propagation direction n µ λ and two orthonormal directionsȇ µ A . Since the photon path is described by the null vector k µ , the infinitesimal Jacobi field can be further simplified as where we used Eq. (2.15) to obtain the relation ξ 0 = −ξ n . The geodesic deviation equation for the Jacobi field is where R µ νρσ is the Riemann tensor. By projecting the Jacobi field into the Sachs basis, the governing equation for the Jacobi field can be derived [44] as The solution of the Jacobi field can be obtained by doubly integrating the source term on the righthand side, and it can be expressed in terms of the linear map, or the Jacobi map D AB as Since the Jacobi field ξ A is the physical separation in a plane perpendicular to the photon propagation direction in the rest-frame of an "observer" with u µ λ , the Jacobi field ξ A at the source position can be related to the physical length and the shape in the source rest-frame, once we correct the difference between the parallel transported velocity u µ λ and the source velocity u µ s . However, this correction is trivial as described in Eq. (5.49). Therefore, the distortion matrix that relates the physical size and the shape in the source rest-frame to the observed angles in the observer rest-frame can be readily constructed by using the Jacobi map [73]. In particular, the anti-symmetric part of the Jacobi map will be related to the lensing rotation. Given the governing differential equation for the Jacobi map, the solution is in general not symmetric (hence non-vanishing rotation), though the source term R AB is symmetric. However, the background solution is trivial D AB ∝ δ AB , and hence the linear-order solution is symmetric with vanishing rotation, in full agreement with our result. The Jacobi mapping approach provides a non-perturbative description of the lensing observables. However, since its formalism is based on the parallel transported tetrad basis, its specification is often absent and hence the relation to the observables is left unspecified. In contrast, our geometric approach provides this missing link, and in particular our explicit calculation of the parallel transport of the tetrad basis illuminates the physical origin of the vanishing rotation, which appears mysterious in the Jacobi mapping approach.

Discussion and Summary
We have presented a gauge-invariant formalism of cosmological weak lensing, accounting for all the relativistic effects associated with the lensing observables at the linear order in perturbations. Without choosing a gauge condition, we have solved the geodesic equation for the light propagation and constructed the tetrad bases for the rest-frames of the observer and the source. The last step is important in establishing the relation of the photon wavevectors to the observed angle and redshift in the observer rest-frame and the physical size and shape in the source rest-frame. We have demonstrated that the standard weak lensing formalism can be naturally generalized to account for all the relativistic effects associated with the light propagation by using the solution of the geodesic equation. However, without specifying the observer and the source rest-frames, the standard lensing formalism is shown to be deficient and gauge-dependent. Using the tetrad bases at the observer and the source positions, we have improved the standard lensing formalism and derived the lensing observables such as the gravitational lensing convergence, the lensing shear, and the rotation in a gauge-invariant way. With full generality, the derivations are lengthier than when a gauge condition is adopted, but by choosing a gauge condition, one loses a way to check the gauge-invariance of the lensing observables. Indeed, very little attention has been paid to this aspect in lensing literature, and some of the gauge issues in the standard weak lensing in Sec. 4 could have been readily spotted by comparing numerical results in another gauge conditions, for instance, the synchronous gauge. We emphasize that it is important to derive equations with the general metric condition and explicitly verify their gauge-invariance. This procedure provides a great sanity check of nonlinear perturbation calculations [41].
The key point in developing the gauge-invariant formalism of weak lensing and improving upon the standard lensing formalism is to identify the rest-frames, in which observable and physical quantities are defined. The tetrad basis vectors e µ a form an orthonormal basis with the Minkowski metric and connect the local rest-frame to the FRW frame, providing the exact ingredient for our purpose. In cosmology, there exists a privileged timelike direction set by the observer four velocity e µ 0 = u µ , which defines the rest-frame of the observer. Furthermore, the spatial directions e µ i in this rest-frame are used to measure any directional quantities in cosmology; for example, the observed propagation direction of the photon wavevector k µ in the FRW frame is measured in the rest frame by contracting against the spatial tetrads as in Eq. (2.1). Any measurements in the observer rest-frame are indeed expressed in terms of diffeomorphism scalars constructed by contracting against e µ a , and this condition automatically guarantees that any observables measured in the rest-frame are gaugeinvariant at the linear order in perturbations [41]. For this reason, it is difficult to overemphasize the significant role of the tetrad bases in cosmology that is unfortunately often neglected in literature. It is indeed more advantageous to go beyond the tetrad bases at the observer and the source positions and consider a tetrad field that consists of the observer families in describing physical observables, as it will lead to a fully nonlinear formalism in cosmology without coordinates [42].
With four tetrad vectors, there exist 16 degrees of freedom, 10 out of which are constrained by the metric g µν = η ab e a µ e b ν at each point. By setting the privileged direction e µ 0 = u µ , we fix three degrees of freedom associated with the boosts, but three remain unconstrained. This remaining freedom is unconstrained by the metric tensor, and it represents the rotation of the spatial tetrad vectors. While this anti-symmetric part of the spatial tetrad in Eqs. (2.6) and (2.11) is often neglected in literature, it is necessary to include the anti-symmetric part to ensure that the spatial tetrad vectors transform as four vectors in a FRW coordinate. Fortunately though, the missing part at the linear order results in the systematic errors, only in the vector perturbations and the rotation of the spatial tetrad vectors that are often subdominant. In the same spirit, the rest-frame of the source is as important as the rest-frame of the observer. For our purpose, the tetrad basis at the source position needs to be used to establish the rest-frame of the source, in which the physical size and shape of the source galaxy are defined. In conjunction with the tetrad basis at the observer position, this part is crucial for deriving the gauge-invariant expressions for the lensing observables. It is well-known [56,68] that the tensor perturbations yield the infrared divergences in the lensing shear and this pathology was fixed by introducing a counter term. However, we have shown that this correction term naturally arises when we transform the FRW frame to the source rest-frame (this correction term is also known as the FNC term [55], 16 or the metric shear [56,74]). Such correction terms due to the frame change exist not only in tensor perturbations, but also in scalar and vector perturbations.
In addition to the frame change, there exists another ingredient often missing in the perturbation calculations -the observer position is different from the positionx µ o = (η o , 0) in a homogeneous universe, as the observer drifts away from the background path in the presence of perturbations. This deviation is characterized by the coordinate lapse δη o in Eq. (3.18) and the coordinate shift δx α o in Eq. (3.19). These perturbations are uniquely determined, and their presence is necessary for deriving the gauge-invariant expressions. These deviations vanish in the comoving gauge, but they cannot be set zero with other gauge conditions. At the linear order, however, the coordinate shift δx α o drops out in all the lensing observables, but the coordinate lapse δη o contributes to the lensing convergence (or the luminosity distance). Their absence in the calculation of the variance of the luminosity distance is shown [46,47] to be the cause of the infrared divergences.
In the presence of the vector and the tensor perturbations, the tetrad basis vectors "rotate" even at the linear order in perturbations as they are parallel transported along the photon path, when they are seen from a global FRW coordinate aligned with the local tetrad basis at the observer position. This rotation of the tetrad basis, often known as the Skrotsky effect [57], translates into the lensing rotation, and its potential measurements may be considered as a probe of the primordial gravitational waves, though its constraining power is expected to be low [56]. However, we have shown that the Skrotsky effect is an artifact of using a global FRW coordinate for observables, and the only physically meaningful way to compare two points in curved space is to parallel transport the basis vectors. We stress that even in the presence of the vector and tensor perturbations the lensing rotation vanishes (but only) at the linear order in perturbations, when compared to the basis parallel transported to the source position. However, this point should not be confused with the statement that the vanishing rotation implies the vanishing lensing B-mode; Instead, the vector and the tensor perturbations contribute to the lensing shear (or non-vanishing B-mode), while the physical rotation is zero, as described in Sec. 5.3.
Indeed, one way to measure the lensing rotation even with a single source was already discussed [75,76] in the past by using the polarization measurements. Though the polarization vector is perpendicular to the photon propagation as the lensing images, it belongs to the tangent space of the central geodesic, while lensing images are extended in space, albeit infinitesimally small. Therefore, polarization is parallel transported (to all orders in perturbations), and hence its measurement can be used to infer the base direction to synchronize the local coordinates at the observer and the source positions, in the absence of any significant magnetic field along the path. 17 Under the assumption that the morphology of the source is aligned with its polarization, one can infer the lensing rotation from a single system, when combined with the polarization measurements.
In summary, the physical lensing observables can be found in Eq. (5.41) for the lensing convergence, in Eq. (5.46) for the lensing shear, and in Eq. (5.47) for the rotation. Compared to the standard lensing formalism, there exist additional relativistic effects in all the lensing observables. The physical rotation is zero (ω = 0) in Eq. (5.47), because the spatial orientation Ω n s of the source frame indeed cancels the remaining terms (see Sec. 5.2). The lensing convergenceκ we measure is indeed the luminosity distance, not the usual (coordinate) convergence κ in Eq. (4.19). Compared to the standard lensing formalism, the additional velocity contributions in the luminosity distance, sometimes referred to as the Doppler lensing, are significant and already measured in current surveys (see, e.g., [47,77]), while the contributions of the gravitational potential or the primordial gravity waves are small, demanding special techniques to be measured in the upcoming surveys (see, e.g., [78]). Looking to the future, we believe that the relativistic effects in large scale structure will provide nearby region is fully contained in the tetrad field. 17 In principle, multi-frequency observations can decode the Faraday rotation due to the magnetic fields.
a great opportunity to probe the nature of gravity and understand the physical mechanism of the perturbation generation in the early Universe. Our gauge-invariant lensing formalism will be essential in providing correct theoretical predictions for the upcoming surveys (see, e.g., [79]). and the vector perturbations transform as and the gauge-invariant variables 15) corresponding to the Bardeen variables: 16) With our interest in the observer four velocity, we introduce a gauge-invariant variable V α by combining the scalar and the vector gauge-invariant variables Given the observed angle and frequency in the observer rest-frame, the photon wavevector k µ in Eq. (2.14) can be constructed in a FRW coordinate. Under the coordinate transformation in Eq. (A.4), this photon wavevector should transform as a four vector, and this transformation property constrains how the conformally transformed wavevectork µ in Eq.