Augusto Rupérez Micola

Albert Banal-Estañol
email: a.banal-estanol@city.ac.uk

Derek W Bunn

Laurent David

Tobias Kretschmer

Olivier Massol

Phanish Puranam

Stéphane Tchung-Ming

Xavier Vives

Herbert Dawid

Giorgio Fagiolo

Incentives and Coordination in Vertically Related Energy Markets *

Keywords: Agent-based modelling, energy markets, reward interdependence. JEL codes: C63, L22, L97

We present an agent-based model of a multi-tier energy market. We show how reward interdependence between strategic business units within a vertically integrated firm can increase its profits in oligopolistic energy markets. The effects are shown to be distinct from those of the raising rivals' costs model. In our case, higher prices relate to the nature of energy markets, which facilitate the emergence of financial netback effects.

A c c e p t e d M a n u s c r i p t 1 Introduction

This paper studies the vertical relationships between gas and electricity markets. Vertical relationships are those that involve an exchange between sequential stages of the value chain. In the energy industry, gas is an important input for electricity generation, and therefore wholesale natural gas and electricity markets are vertically interrelated. The same is true for wholesale and retail electricity markets since retailers buy electricity from wholesalers (Stern 1998).

Vertical interactions may involve separate firms or different strategic business units (SBUs) within the same firm (Gulati et al. 2005). Vertical integration is widespread among European energy firms. Gas producers often own gasfired power plants, and many electricity firms consist of generation and retail SBUs (Finon and Midttun 2004). Mergers between gas and electricity firms are relatively new in the United States, but they are occurring at a rapid pace (Hunger 2003). Moreover, the merger activity appears to be accelerating as competition opportunities expand, incentive regulation diffuses more widely, and regulators have become less hostile to mergers (Joskow 2000).

Several streams of literature have studied the advantages as well as the disadvantages of vertical integration. Industrial economists have extensively analysed whether vertically related firms could benefit from foreclosing nonintegrated rivals (see Rey and Tirole 2004 for a recent survey). [START_REF] Bushnell | Vertical arrangements, market structure, and competition: An analysis of restructured U.S. electricity markets[END_REF] Ordover et al. (1990), for example, show that a vertically integrated firm in a bilateral duopoly may benefit from disadvantaging its downstream rival. [START_REF] Chatterjee | Vertical strategies and market structure: A systematic risk analysis[END_REF] The argument runs as follows. If the upstream unit ceases to sell in the input market, its upstream rival will face less competition and raise prices. Higher input prices increase the costs of the non-integrated downstream rival. This firm is forced to reduce production and increase prices. As a result, the vertically integrated firm can increase profits by raising both its end-user market share and price. [START_REF] Choi | Vertical foreclosure with the choice of input specifications[END_REF] Studies of vertical relationships in energy markets (e.g. Granitz and Klein 1996, Bushnell et al. 2005) often explain their findings using this foreclosure argument. However, its logic depends crucially on the firm's ability either to internalise transactions or to set an internal transfer price that is different from the external (input) market price. In practice, wholesale energy markets are often compulsory, so trading internalisation is not feasible. Moreover, the standard energy market mechanism is the uniform price auction, which seems to [START_REF] Bushnell | Vertical arrangements, market structure, and competition: An analysis of restructured U.S. electricity markets[END_REF] Management scholars have identified several other motives for firms to integrate vertically (Harrigan 1984 and 1986), including the reduction of transaction costs (Williamson 1975, Mahoney 1992), the reduction of corporate risk (Chatterjee et al. 1992) and the elimination of the double marginalisation inefficiency (see e.g. Gaudet and van Long 1996). 2 Salinger (1988), on the other hand, analyses whether vertical integration leads to higher or lower prices in a successive oligopoly setting. [START_REF] Choi | Vertical foreclosure with the choice of input specifications[END_REF] Ordover et al. assumed that the vertically integrated firm can commit to limit its supplies to the downstream rivals and that the upstream competitor can charge only linear prices. Several authors have relaxed some of their assumptions. In particular, Choi and Yi (2000) and Ma (1997) dispense in different settings with the commitment assumption, although not with the linear pricing one. Allowing for a broader set of tariffs, Hart and Tirole (1990) show that vertical integration may also benefit the integrated firms if the upstream unit is more efficient than its competitors.

Page 3 of 19

A c c e p t e d M a n u s c r i p t make differences between internal and external prices impossible at the outset. Thus, two of the main resting points of the foreclosure or raising rivals' costs logic are often not present in energy markets.

This paper introduces an agent-based simulation model of natural gas and electricity markets. We show how corporate incentives linking the strategic business units within a vertically integrated firm can increase prices in energy markets. The effects are shown to be different from those arising in the standard models of vertical foreclosure. In our case, higher prices are related to the existence of financial netback effects in energy markets.

Wholesale natural gas is often priced against wholesale electricity prices, which, in turn, are usually set with reference to retail tariffs. Industry players refer to this financial process as netback or spark spread pricing. The use of netback pricing is well documented in the energy markets literature (for an early discussion see Moxnes 1987) and suggests a sequential relationship between energy markets in which gas and electricity markets are cleared in a down-toupstream sequence. This is in contrast to the Ordover et al.'s "physical flow" formulation, widely adopted by the foreclosure literature, whereby upstream prices are determined before those downstream.

Trading is not a simple sequential procedure in most industry supply chains, as it is often iterative and tentative before deals are struck, so that modeling it either way is an abstraction. However, in the energy case, there are good reasons to support netback pricing rather than the physical flow formulation. First, retail prices are generally fixed for longer periods of time than upstream prices. For example, in the UK, retail consumers cannot change supplier during the first twenty eight days of signing a contract, but upstream prices vary every half hour. Therefore, retailers and generators take downstream prices as given when trading in the wholesale market. [START_REF] Finon | Reshaping European Gas and Electricity Industries[END_REF] Second, electricity is non-storable, consumption and generation have to be balanced at all times, but end-users are free to choose volume. As a consequence, information on retail consumption volumes has to flow upwards in the supply chain via prices. [START_REF] Gaudet | Vertical integration, foreclosure and profits in the presence of double marginalisation[END_REF] In this context, it seems plausible that widely used vertical incentives should play a role on the determination of energy prices. The general reward system of an organisation influences the behavioural choices of its members. Bonuses tied to overall profits create incentives for cooperative behaviour both between individuals (Zander and Wolfe, 1964; Wageman and Baker, 1997) and across de-A c c e p t e d M a n u s c r i p t partments within a firm (Petersen 1992, Kretschmer and Puranam forthcoming). For individuals, the more interdependent the task, the more interdependent the reward system should be (Wageman and Baker) because it results in a positive relationship between effectiveness of the integrative devices and organisational performance (Lawrence and Lorsch 1967). For firms, the importance of cooperation between SBUs grows with their interdependence (Gulati and Singh 1998) and the higher the inter-unit synergies, the more useful the collaborative incentives are (Kretschmer and Puranam). Collaborative incentives, however, not only encourage cooperation but may also enhance free riding. Indeed, rewards based on aggregate profits hinder the identification of individual performances. As a consequence, individuals have more incentives to shirk hoping that the others will compensate (Holmstrom 1982, Petersen 1992). [START_REF] Granitz | Monopolization by "Raising Rivals' Costs": The Standard Oil case[END_REF] Despite the importance of collaborative incentives, the existing literature provides no guidance as to how they should be given to sequential SBUs in vertically integrated energy firms. In order to fill this gap, we consider a setting consisting of two sequential, multiple-unit, compulsory, uniform-price auctions representing a wholesale and a retail energy market. Although quite realistic, this complex trading environment presents a manifold of non-Pareto ranked Nash equilibria (von der Fehr and Harbord 1993). To achieve predictions, we adopt an inductive selection method based on the adaptive theory of reinforcement learning put forward by Roth and Erev (1995).

The agent-based simulations show that coordination overcomes the potential disadvantages of broad collaborative incentives due to the large interdependences between energy markets. More importantly, our results uncover a simple but powerful mechanism to exert vertical market power. Using collaborative incentives that link the reward to the performance of the different SBUs, vertically integrated firms induce higher prices and achieve higher profit. These observable outcomes are similar to those of the foreclosure argument. Closer inspection, however, reveals that our downstream SBU behaves less competitively, increasing downstream prices at the expense of market share rather than taking advantage of the rivals' higher costs, as in the foreclosure argument. Moreover, the upstream SBU behaves more competitively and benefits from the higher prices downstream.

The remainder of the paper is organised as follows. Section 2 outlines the agent-based simulation model. Section 3 presents the results, which are extended in Section 4. A short discussion follows in Section 5.

The Computational Model 2.1 General Setting

The model incorporates key features of energy markets in the short-run. Consider two sequential, oligopolistic markets, a "wholesale" and a "retail" market.

Page 5 of 19 A c c e p t e d M a n u s c r i p t

They can represent a wholesale and a retail electricity market. Alternatively, one can think of wholesale gas and electricity markets, whereby gas sold in the gas market is an input for the generation of electricity. In the wholesale market, there are A wholesalers that sell energy to B retailers. These, in turn, re-sell this energy in the end-user market. [START_REF] Gulati | Adaptation in vertical relationships: Beyond incentive conflict[END_REF] Marginal costs are assumed to be constant throughout and normalised to 0 for simplicity. There are no transmission constraints or storage.

Firms are assumed to be capacity constrained. Generators, for example, are limited by the aggregated maximum capacity of their power plants. Wholesale gas sales are also typically constrained by long-term contracts signed between gas producers and shippers, together with high penalties to discourage overutilization of the network. As a stylization to isolate the effects emerging from vertical incentives, we focus on the situation in which capacity is equal across tiers (i.e. there are no particular bottlenecks in the supply chain). Moreover, firms in each tier are assumed to be identical. Hence, if one denotes market capacity as K, the individual wholesaler capacity is K w = K A and that of retailer

K r = K B .

Market Rules

Goods are traded repeatedly along the value chain (Figure 1). In a given round t, two uniform price auctions take place, first at the retail and then at the wholesale level, following the netback pricing sequence. In each market i, i ∈ {r, w}, trading occurs as follows. Suppliers simultaneously submit single price bids at which they are willing to sell (up to) their capacity, starting from 0 and up to P i (t), a maximum level for each market described below. An independent auctioneer determines the uniform market price P i (t) by intersecting the ad hoc supply function with the corresponding demand curve, Q i (t). The independent auctioneer assigns full capacity, q i j (t) = K i , to the m sellers with bids below the market price; the remaining capacity, q i j (t) = Q i (t) -mK i , to the sellers with bids equal to the market price;8 and zero sales, q i j (t) = 0, to those submitting above the market price. The price and the individual quantities are then communicated independently to each supplier.

We follow the usual simplification of modeling retail demand as fully inelastic. The literature has established the extremely low price elasticity of short term electricity demand, originating among others from the lack of real-time metering systems (e.g. Stoft 2002). In particular, the inelastic market demand Q r (t) is drawn, independently in each round, from a uniform distribution in [ Qrε, Qr + ε], where Qr is the expected end-user demand and ε accounts for the small uncertainty typical in day-ahead forecasting. We assume that there always is overall overcapacity, Qr + ε < K, consistent with the normal operations A c c e p t e d M a n u s c r i p t of most de-regulated energy markets. [START_REF] Harrigan | Formulating vertical integration strategies[END_REF] Hence, our model excludes blackouts due to extreme weather or technical failure.

Possible retail prices are bounded between (0, Ψ], with Ψ being the maximum reasonable end-user price, and therefore P r (t) = Ψ. [START_REF] Harrigan | Matching vertical integration strategies to competitive conditions[END_REF] Retail commitments are honoured with purchases in the wholesale market. The retailers' aggregated demand curve in the wholesale market is equal to the market demand if the price is below the retail price, Q w (t) = Q r (t) if P w (t) ≤ P r (t), and zero otherwise, Q w (t) = 0 if P w (t) > P r (t). Accordingly, wholesalers submit bids bounded between (0, P r (t)] and P w (t) = P r (t).

Profits for each firm type are

π w a (t) = P w (t) q w a (t) for a = 1...A (1) π r b (t) = [P r (t) -P w (t)] q r b (t) for b = 1...B. (2)
There are apparently no other multi-tier simulations driven by netback principles in the energy modelling literature. This new method is the paper's main methodological contribution.

Vertical Integration and Reward Interdependence

In the basic model, it is assumed that a wholesaler (without loss of generality, a = 1) and a retailer (b = 1) are vertically integrated in that they belong to Page 7 of 19 A c c e p t e d M a n u s c r i p t the same organisational structure (i.e. the same firm). Trading is compulsory and firms are not allowed to price discriminate or set an internal transfer price. However, the vertically integrated firm can influence its traders' decisions with incentives depending not only on their own SBU profits, but also on those of the other SBUs,

Ω w 1 (t) = (1 -α)π w 1 (t) + απ r 1 (t) and (3) Ω r 1 (t) = (1 -α)π r 1 (t) + απ w 1 (t), (4)
where α = {0, .01, .02, ..., .5} parameterises the "reward interdependence" (Wageman and Baker) between the two vertically related SBUs. A small α represents narrow incentives, which become broader for growing α. Note that for α = 0, Ω w 1 (t) = π w 1 (t) and Ω r 1 (t) = π r 1 (t), and SBUs trade as if they were independent. For tractability and realism, the model is restricted to the case in which SBUs are rewarded predominantly on the basis of their own performance (i.e. α ≤ 0.5).

Managers in the non-integrated firms do not have reward interdependences, so their incentives are correlated to their own unit performance,

Ω i j (t) = π i j (t) for i = {w, r} and j 6 = 1. (5)

Bidding and Behavioural Learning

The feasible price offer domain for each firm is approximated by a discrete grid consisting of a fixed number of possible actions (independent of t). In each trading period, suppliers choose among S i possible prices, equally spaced between the minimum and the maximum price offer, (0, P i (t)]. Hence, the set of possible actions, A i at a tier i in a period t is given by

A i (t) = s * Ã P i (t) S i ! for s = 1, ..., S i . (6)
Notice that in the wholesale market, the set of possible prices changes over time.

In both markets, actions with a lower s are more competitive or closer to the marginal costs. Each trader plays each possible action with a given likelihood or "propensity", r i j,s . The probability that an agent j plays an action s is given by its propensity divided by the sum of the propensities of all possible actions,

p i j,s (t) = r i j,s (t) P S i s=1 r i j,s (t)
.

Propensities for all actions are initialised to the firms' maximum profit, r i j,s (1) = ΨK i for all s and i, so that all actions have the same initial probability, p i j,s (1) = 1 S i for all s and i.

6

A c c e p t e d M a n u s c r i p t

At the end of each round, traders reinforce the selected action, k, through an increase in its propensity equivalent to the performance, Ω i j (t). Moreover, actions that are similar, k -1 and k + 1, are also reinforced, by Ω i j (t) * (1δ) where 0 < δ < 1 ("persistent local experimentation" in the terminology of Roth and Erev). All propensities are discounted by γ ("gradual forgetting"), and actions whose probability falls below a certain threshold are removed from the space of choice ("extinction in finite time"). Summarising, the pre-extinction propensities r i0 j,s are

r i0 j,s (t) = ⎧ ⎨ ⎩ (1 -γ) r i j,s (t -1) + Ω i j (t) if s = k (1 -γ) r i j,s (t -1) + (1 -δ) Ω i j (t) if s = k -1 or s = k + 1 (1 -γ) r i j,s (t -1) if s 6 = k -1, s 6 = k and s 6 = k + 1,
and the final propensities, corrected by the extinction feature, are

r i j,s (t) = r i0 j,s (t)I { r i0 j,s (t) S i s=1 r i0 j,s (t) >μ} ( 8
)
where I is an indicator function that takes value 1 if the condition between brackets is satisfied and zero otherwise.

Simulation Parameters

In the first instance, industry structure is simplified to consist of two wholesalers and three retailers (i.e. A = 2; B = 3). Remember that the two vertically integrated entities are denoted by a = 1 and b = 1. Total capacity is set to K = 300 so that the individual capacity of a wholesaler is K w = 150 and of a retailer K r = 100. Expected market demand is Qr = 240 and ε = 5, and hence there is an expected excess capacity of 20% with about 5% uncertainty in the day-ahead forecasted demand, approximately of the magnitude observed in energy markets. The end-user reasonable price ceiling is set at Ψ = 200.

We study fifty-one reward interdependence cases α = {0, .01, .02, .03, ..., .50} with 50 simulation runs consisting of 500 periods each. The data consists of averages for the last 200 periods for each run and case. These represent long term stationary values to which the three markets converge, based on the Roth and Erev reinforcement parameters.11 ' [START_REF] Holmstrom | Moral hazard in teams[END_REF]

Results

This section describes the simulation results. We report first the results on market prices and firm profits. We then turn to the analysis of firm behaviour, Page 9 of 19

Market Prices and Profits

Figures 2 and 3 report the 95% mean confidence intervals of simulated prices in the two sequential markets. They present the relationship between reward interdependence (α, on the horizontal axis) and the wholesale and retail prices (P w and P r , on the vertical axis), respectively. As shown by Figure 2, wholesalers coordinate on higher prices as α grows. The lack of reward interdependence between the integrated wholesaler and retailer (α = 0) results in an average wholesale price of approximately 95 monetary units, which increases to about 102 units for overall performance-based incentives (α = .50). That is, if the SBUs of the vertically integrated firm traded as if they were independent, wholesale market prices would be lower.

Figure 3 shows that retail prices are also increasing in the level of reward interdependence. When α = 0, the simulation produces an average price of about 132; for α = .50 it is about 148. As a result, the expected absolute size of the resource rent shared by the two tiers (= P r Q r ) also increases. The proportion accruing to the retailers, however, does not increase. Retail prices increase only to compensate the higher wholesale prices. Increasing wholesale prices with constant production costs results in a higher proportion of the rent staying with the wholesale duopoly. As a consequence, α influences both the profits of the vertically integrated firm and its two SBUs (Figure 4). Although broader vertical incentives do not increase the profits of its retail unit, they increase those of its wholesale unit. The resulting overall corporate profit (π w 1 + π r 1 ) is clearly increasing from about 14,300 monetary units to over 16,000. Thus, there is a positive relationship between vertical reward interdependence across the two strategic business units, the firm's overall profit and the concentration of profit upstream.

Summarising, as in Ordover et al., wholesale and retail prices increase as a result of vertical integration. Indeed, the vertically integrated firm achieves higher profit by using rewards based on overall firm performance. With this type of reward interdependence, prices are higher than if the two business units had been independent.

In order to provide an explanation for these results one needs to explore the individual behaviour of firms, to which we turn next.

Price-setting Behaviour

In a uniform auction, prices are determined by the last bid on the merit order, and therefore being a marginal seller is a necessary condition to influence prices. Hence, the frequency with which each firm takes the marginal price-setting position on the supply merit order represents one dimension of its predisposition to exert market power (at the expense of market share).

The vertical axes in Figures 5 and6 provide averages over the frequencies with which each market player sets the price in the 200 end-of-simulation periods averaged across the 50 simulation runs. An increase in the reward interdependence in the vertically integrated firm creates two simultaneous effects. On the one hand, the proportion of trading periods in which the wholesale unit sets prices goes down from about 50% when α = 0 to less than 40% when α = 0.5 (Figure 5). As in Ordover et al., vertical integration provides incentives for the wholesale unit to concede market power upstream and, as a result, for the competitor to push up the wholesale price.

On the other hand, the proportion of trading periods in which the integrated retail unit is price-setting increases from 33% when α = 0 to about 50% for α = 0.5 (Figure 6). As a result, there is a negative relationship between α and the retail unit's market share, which is at odds with Ordover et al.'s findings.

Latent Intensity of Competition

A similar behaviour can be identified in the S i end-of-simulation individual latent probability distributions from which agents choose bids, depicted in Figures 7 and8. The concentration of probabilities is largely invariant across a large number of periods once the market reaches convergence, so the distributions on the last trading period are an indication of the firms' long-term mixed strategies.

On the horizontal axes, strategies are identified with numbers ranging from 1 for the more competitive to 100 for the highest possible bid. Cumulative probabilities for the two tiers are calculated on the vertical axes for each element of the strategy space. The curves summarise the cumulative bidding A c c e p t e d M a n u s c r i p t probabilities for α = 0 and α = 0.5, averaged across the 50 simulation runs. In all cases, probabilities concentrated on lower and higher strategies result in the agents behaving more and less competitively, respectively. That is, curve movements to the upper-left and lower-right corners are indications of more and less competition, respectively. [START_REF] Hunger | Analyzing gas and electric donvergence mergers: A supply curve is worth a thousand words[END_REF] The figures offer a number of insights linking individual probability distributions to market outcomes. First, the wholesalers' distributions are very similar under α = 0 but become different under α = 0.5 (Figure 7). Reward interdependence incentives in the vertically integrated firm have the effect of making the integrated wholesaler's bids more competitive (upper-left movement). Moreover, the rival wholesaler's prior is slightly less competitive (movement to the lower-right). Market power is concentrated on the non-integrated wholesaler. Second, probability distributions on the retailer side are similar for all firms under α = 0 (Figure 8). However, when α = 0.5, the integrated retailer tends to exert more market power than its competitors (movement to the lower-right).

Firm Learning, Behaviour and Market Outcomes

Through their dynamic trading interaction, firms learn to prioritise those bidding strategies that achieve higher payoffs and choose them more often. Each firm's price setting frequency is related to the strategy reinforcement, and once marginal supply and demand patterns are established, price regularities follow. The results suggest a link between reward interdependence, firm learning, trading behaviour and market outcomes.

As reward interdependences increase, the retail SBU submits higher bids, sets prices more often and increases them markedly, albeit at the expense of market share. However, wholesale prices are bounded below retail prices due to the netback procedure. A higher retail price shifts up the wholesale demand curve and, thus, provides more scope for higher wholesale prices. Then, the wholesale SBU trades more often on the baseload, making it easier for its opponent to increase prices. As a consequence, the integrated firm improves its upstream market share and profitability, and this has a positive effect on the total profits of the vertically integrated firm.

Our results yield prices that are superficially equivalent to those of foreclosure, which, by construction, cannot occur in our setting. We uncover a different mechanism to exert vertical market power in which the downstream SBU behaves less competitively, gives up market share and pushes up retail prices, rather than benefit from its' rivals' higher costs, as in Ordover et

Extensions

Generators may buy their inputs in wholesale markets one step up in the value chain. For example, gas sold in the wholesale market may be an input for the generation of the electricity traded in the wholesale electricity market that, in turn, is sold to end users in the retail market. This section extends the previous analysis from two to three tiers. We also analyse whether the results depend on the position of the vertically integrated SBUs in the value chain and on market concentration levels.

Three-tier Markets

Consider three sequential energy markets, with the same characteristics as before. [START_REF] Joskow | Deregulation and regulatory reform in the U.S. electric power sector[END_REF] Suppose that in the gas market there are two upstream natural gas shippers that sell gas to three electricity generators. These generators buy gas to produce electricity and sell it in the wholesale electricity market to four electricity retailers, who re-sell in the end-user market. New simulations analyse reward interdependences between a gas shipper and a generator plus those between a gas shipper and a retailer and a generator and a retailer, keeping all other parameters constant. Table 1 summarises results for α = 0 and α = .5 in these new simulations.

Similar to the previous results, reward interdependences between a natural gas shipper and an electricity generator result in an increase in both wholesale gas and electricity markets. Wholesale gas prices, P g , increase from an average of 59 for α = 0 to 63 when α = .5 (+6.77%) whereas wholesale electricity prices, P e , increase from 83 to 93 (+12.04%). Interestingly, there is a negligible change in retail electricity prices, P c (+0.88%), and therefore end-users would not be hurt here. Market prices seem to indicate that retailers are particularly damaged from reward interdependences between a gas shipper and a generator.

Two complementary findings emerge from the simulated consequences of reward interdependences linking the fortunes of a natural gas shipper and a retailer. First, retail prices increase from an average of 124 for α = 0 to 144 when α = .5 (+16.12%), compared to a negligible change when α links a shipper and a generator. Reward interdependence induces the retailer to implement a less competitive-mixed strategy through which the firm becomes price setting more often with the retail price increasing steadily.

Secondly, wholesale prices increase, too. The effect on retail prices moves up in the value chain. It translates into higher wholesale electricity and gas prices from about 83 to 95 (+14.45%) and 59 to 67 (+15%), respectively. Reward interdependences hurt end users in this case. Although the generators' probability priors and trading behaviour change little, the integrated shipper becomes more competitive and facilitates the exertion of market power by its upstream rival, which leads to a wholesale electricity-gas difference similar to the initial case. Moreover, by construction, higher retail prices widen the range within which wholesale electricity prices are determined. With a constant number of bidding strategies, S i , similar generators' probability distributions result in higher wholesale prices as the retail price increases.

The simulations assuming integration between a generator and a retailer show increases in wholesale electricity and retail prices (+6.02% and +9.67%). The downstream unit manages to leverage its overall revenue and implement a more collusive mixed strategy that expands the base for higher wholesale electricity prices. The upstream unit, on the other hand, bids more often on the base-load part of the supply curve, with the other generators setting (higher) prices. It is interesting to note how the non-integrated gas shippers capture a large proportion of the rent through a higher wholesale natural gas increase of 6.77%.

Symmetric Market Structure across Tiers

In the previous analysis, market structure remained progressively less concentrated. The vertical foreclosure literature suggests that higher upstream concentration is an important element of the firms' ability to exert vertical market power. We check this result in our setting, through a symmetric-across-tiers simulation.

Table 2 summarises market prices for the case in which there are two firms in each tier. Vertical integration between a shipper and a generator leads, as in the asymmetric case, to higher wholesale electricity and natural gas prices but not higher end user prices. Higher prices with increasing α, therefore, hinge on Page 17 of 19 A c c e p t e d M a n u s c r i p t the market clearing sequence rather than on market structure asymmetries. [START_REF] Kretschmer | Exploiting Synergies in Complex Organizations: When do Collaborative Incentives Create Value? Forthcoming in Organization Science[END_REF]

Discussion

The existence of interdependences between vertically related SBUs has become a bedrock in the business strategy literature, yet we know relatively little about the forces shaping these interdependences and their effects. A set of agentbased simulations identifies one such effect in the de-regulated energy industry. Reward interdependences between SBUs lead to trading coordination and higher prices. Under tight reward interdependence structures, vertically integrated firms increase the scope for upstream profits. This leads to strategic behaviour with the superficial appearance of foreclosure, but based on a quite different principle. This paper adds to the preceding literature in at least three ways.

Methodologically, we are not aware of any other multi-tier energy simulations driven by netback principles in the literature. The financial dependence between electricity and natural gas markets does not seem to be captured by the classical models of vertical foreclosure, where causal pricing relationships are sequential from the upstream to the downstream segment.

The second contribution relates to the literature on the sources of vertical market power in the energy industry (e.g. Bushnell et al. 2005, Granitz and Klein 1996, Kühn and Machado 2004). The simulations suggest a new mechanism to solve the puzzle of how vertical market power appears in some energy markets where it should not. The downstream unit submits higher bids and increases retail prices. As a result, the range of possible upstream prices also increases. Then, the upstream unit stays on the baseload part of the supply curve, increases its market share and allows its opponent to increase wholesale prices. Hence, vertical market power emerges in a compulsory, uniform price auction without trading internalisation or price discrimination.

Thirdly, the research identifies a link between internal incentive structures, SBU behaviour and firm performance. Reward interdependence has been shown to be an instrument leading to market power via higher vertical SBU coordination. Ways in which reward interdependences can be articulated include direct bonuses and stock options, and casual evidence indicates that these are widespread in the energy industry. It is interesting to note that such reward interdependence contracts are internal to the firm and, hence, normally fall outside the scope for regulatory intervention. Whether firms use them explicitly as a way of aligning their interests to those of their SBU employees is an interesting question for future empirical work.

Page 18 of 19

Figure 1 :

 1 Figure 1: Sequential clearing of the market simulation.

  A c c e p t e d M a n u s c r i p t

Figure 2 :

 2 Figure 2: Wholesale prices

A c c e p t e d M a n u s c r i p t

Figure 3 :

 3 Figure 3: Retail prices.

Figure 4 :

 4 Figure 4: Profit, vertically integrated firm.

  A c c e p t e d M a n u s c r i p t

Figure 5 :

 5 Figure 5: Wholesalers' price-setting frequencies.

Figure 6 :

 6 Figure 6: Retailers' price-setting frequencies.

  A c c e p t e d M a n u s c r i p t

Figure 7 :

 7 Figure 7: End of simulation distribution of strategies of wholesalers.

Figure 8 :

 8 Figure 8: End of simulation distribution of strategies of retailers.

  s c r i p t

Table 1 :

 1 al. The identification of this mechanism is this paper's main economic policy contribution. Average Prices for a market with two gas shippers, three electricity generators, and four electricity retailers.

	Page 15 of 19

Table 2 :

 2 Average Prices for a market with two gas shippers, two electricity generators, and two electricity retailers.

In contrast, in other markets, one can argue that upstream producers fix their prices before downstream firms compete. In the typical example of a producer and a retailer, the retailer buys from the producer and then competes in the downstream market. As such, the downstream firm would take the upstream price as given (the item is already bought).

In practice there are several cycles of this. Ahead of real-time, retailers will make forecasts of demand and seek to contract power to cover it. Generators will make forecasts of output and purchase fuels to cover it. In real time, end-user demand is discretionary but at a preset fixed price. Retailers will have to cover imbalances between their contract and actual offtake at the spot electricity imbalance price, which is actually set slightly ex post depending upon what the system operator has to purchase in real time to ensure system stability. Similarly, there will be real-time energy balancing in the gas network adapting to what the power stations have had to do.

Broad incentives could also obstruct learning since it is more difficult to identify the most successful business strategies.

Although relevant in the medium term, we do not deal with entry and exit of firms, variation in expected end-user demand or capacity expansion.

In case of a tie, the selling firm is selected randomly.

For example, the UK energy system includes a reserve margin of about 20% of expected peak demand.

This upper price ceiling can be understood as a limit triggering regulatory intervention or the cost of alternative, expensive, peaking load fuels to which the system administrator could switch at short notice if prices exceed Ψ.

For each level of reward interdependence (α), the one-lag with trend Augmented Dickey-Fuller test-statistics for both wholesale and retail price series are lower than -10. Given that the 95% critical value is -3.43, the null hypothesis of non-stationarity is always rejected.

[START_REF] Holmstrom | Moral hazard in teams[END_REF] The Roth and Erev parameters used throughout are γ = 0.01, δ = 0.5 and μ = 0.0005.

Notice that how competitive the strategies are is determined with respect to the maximum bid, and therefore more competitive strategies could translate into higher absolute bids.

See the working paper version of the paper for more details on the model and results, as well as for additional figures, for this three-tier case.

Results under two symmetric firms in two tiers and three symmetric firms in three tiers are qualitatively equivalent to the case of two symmetric firms in three tiers.


Kazuya Koyama

Non-Gaussianity of quantum fields during inflation

In this review, we discuss how non-Gaussianity of cosmological perturbations arises from inflation. After introducing the in-in formalism to calculate the n-point correlation function of quantum fields, we present the computation of the bispectrum of the curvature perturbation generated in general single field inflation models. The shapes of the bispectrum are compared with the local-type non-Gaussianity that arises from non-linear dynamics on super-horizon scales.

I. INTRODUCTION

There is currently a great deal of interest in the statistical properties of primordial perturbations from inflation, because measurements of any non-Gaussianity will improve by about an order of magnitude over the next few years, for example with the Planck [1]. This will provide a key way to discriminate between the many models of inflation. Although single field models of slow-roll inflation typically generate a small level of non-Gaussianity [2,3], there may be an observable level generated in many alternative models of early universe .

We are interested in the primordial curvature perturbation on uniform density hypersurfaces, ζ, on large scales, which is directly related to temperature anisotropies in Cosmic Microwave background (CMB). The power spectrum of ζ is defined as

ζ(k 1 )ζ(k 2 ) = (2π) 3 δ (3) (k 1 + k 2 )P ζ (k 1 ), (1)
where ζ(k) is a Fourier component of ζ(t, x i ). If ζ obeys Gaussian statistics, the power spectrum determines all statistical quantities. However, if the distribution function of ζ deviates from Gaussian statistics, we need to specify higher order statistics. The first non-trivial statistics is the bispectrum defined by

ζ(k 1 )ζ(k 2 )ζ(k 2 ) = (2π) 3 δ (3) (k 1 + k 2 + k 3 )B ζ (k 1 , k 2 , k 3 ). ( 2
)
Currently observations of the CMB have concentrated on constraining the 3-point function (bispectrum) [59][60][61].

The most instructive way to understand how non-vanishing bispectrum appears from non-linearities is to use the so-called delta-N formalism [62][63][64][65]. See [66] for a concise review of the delta-N formalism. This is based on the separate universe approach [67][68][69]. This considers each super-Hubble scale patch to be evolving like a separate Friedman-Robertson-Walker universe which is locally homogeneous. By patching these regions together we can track the evolution of the curvature perturbation on large scales just by using background quantities.

The number of e-foldings, N , given by

N = t fin tini H(t)dt , (3)
is evaluated from an initial flat hypersurface to a final uniform-density hypersurface. The perturbation in the number of e-foldings, δN , is the difference between the curvature perturbations on the initial and final hypersurfaces. We wish to calculate primordial perturbations, hence we pick a final uniform density hypersurface to be at a fixed time during the standard radiation dominated era, for example during primordial nucleosynthesis. The initial time is arbitrary provided it is after the Hubble exit time of all relevant scales. It is often convenient to pick this time to be shortly after Hubble-exit.

Let us consider a model where a scalar field φ determines the expansion of the Universe. The scalar field acquires quantum fluctuations δφ under horizon scales which cause fluctuations in e-foldings in each super-Hubble scale patch. Then ζ can be written as ζ(t, x i ) = N ,φ (t, t * )δφ(t * , x i ) + 1 2 N ,φφ (t, t * ) δφ(t * , x i ) 2 -δφ(t * , x i ) 2 + ...,

where t * denotes the horizon crossing time and N ,φ = dN/dφ and N ,φφ = d 2 N/dφ 2 . Using Eq. ( 4), the bispectrum of the curvature perturbation is calculated as

ζ(k 1 )ζ(k 2 )ζ(k 3 ) = N 3 ,φ δφ(k 1 )δφ(k 2 )δφ(k 3 ) + N ,φφ N 2 ,φ 2
δφ(k 1 )δφ(k 2 )[δφ δφ](k 3 ) + perms, (5) where [δφ δφ] denotes a convolution

[δφ δφ](k) = d 3 q (2π) 3 δφ(q)δφ(k -q). (6)
A diagrammatic approach to compute the higher order correlation function is developed in Ref. [70].

There are in general two contributions to the bispectrum of the curvature perturbations. One arises from the bispectrum of quantum fluctuations of the scalar field generated under horizon scales, δφ(k 1 )δφ(k 2 )δφ(k 3 ) . The other is coming from the non-linear evolution of the scalar field on super-horizon scales determined by the second derivative of N with respect to the field, N ,φφ . The latter contribution exists even if the field perturbations at the horizon crossing are Gaussian. This contribution is often called the local type non-Gaussianity as this arises from a non-linear local relation between the curvature perturbation ζ and the field perturbations δφ in real space. On the other hand, the contribution from the non-linearity of quantum fields depends on non-linear interactions under horizon scales and its k-dependence is in general very different from the local type non-Gaussianity. In this review, we derive this contribution in several inflation models and compare its shape with the local type non-Gaussianity.

The structure of this review is as follows. In section II, we review the in-in formalism to calculate the npoint function of quantum fields. In section III, the effective action for higher order perturbations is derived in two gauges that are necessary to compute the interaction Hamiltonian at third order. We calculate the bispectrum of the curvature perturbation in slow-roll inflation and k-inflation. In section V, we compare the shape of the bispectrum in k-inflation models with that in the local-type non-Gaussianity. Section VI is devoted to conclusions.

II. QUANTUM CORRELATIONS IN THE IN-IN FORMALISM

In this section, we review the in-in formalism to calculate the n-point function of quantum fields by following Ref. [71]. Consider a general Hamiltonian system, with canonical variables φ(x, t) and conjugates π(x, t) satisfying the commutation relations

φ(x, t), π(y, t) = iδ 3 (x -y) , φ(x, t), φ(y, t) = π(x, t), π(y, t) = 0 , (7)
and the equations of motion

φ(x, t) = i H[φ(t), π(t)], φ(x, t) , π(x, t) = i H[φ(t), π(t)], π(x, t) . (8)
The Hamiltonian H is a functional of the φ(x, t) and π(x, t) at fixed time t, which according to Eq. ( 8) is independent of the time at which these variables are evaluated. We assume the existence of a time-dependent classical number solution φ 0 (x, t), π 0 (x, t), satisfying the classical equations of motion and we expand around this solution, writing

φ(x, t) = φ 0 (t) + δφ(x, t) , π (x, t) = π 0 (t) + δπ(x, t) . (9)
Since classical numbers commute with everything, the fluctuations satisfy the same commutation rules (7) as the total variables. Now, although H generates the time-dependence of φ(x, t) and π(x, t), it is H rather than H that generates the time dependence of δφ(x, t) and δπ(x, t) where H[δφ(t), δπ(t); t] is the sum of all terms in H[φ 0 (t) + δφ(t), π 0 (t) + δπ(t)] of second and higher order in the δφ(x, t) and/or δπ(x, t):

δ φ(x, t) = i H[δφ(x, t), δπ(x, t); t], δφ(x, t) , δπ(x, t) = i H[δφ(x, t), δπ(x, t); t], δπ(x, t) . ( 10
)
This then is our prescription for constructing the time-dependent Hamiltonian H that governs the timedependence of the fluctuations: expand the original Hamiltonian H in powers of fluctuations δφ and δπ, and throw away the terms of zeroth and first order in these fluctuations. It is this construction that gives H an explicit dependence on time. It follows from Eq. (10) that the fluctuations at time t can be expressed in terms of the same operators at some very early time t 0 through a unitary transformation

δφ(t) = U -1 (t, t 0 )δφ(t 0 ) U (t, t 0 ) , δπ(t) = U -1 (t, t 0 )δπ(t 0 ) U (t, t 0 ) , (11)
where U (t, t 0 ) is defined by the differential equation

d dt U (t, t 0 ) = -i H[δφ(t), δπ(t); t] U (t, t 0 ), (12)
and the initial condition U (t 0 , t 0 ) = 1. In the application that concerns us in cosmology, we can take t 0 = -∞, by which we mean any time early enough so that the wavelengths of interest are deep inside the horizon. From now on, we will omit the x dependence of δφ and δπ to simplify the notation.

To calculate U (t, t 0 ), we now further decompose H into a kinematic term H 0 that is quadratic in the fluctuations, and an interaction term H I :

H[δφ(t), δπ(t); t] = H 0 [δφ(t), δπ(t); t] + H I [δφ(t), δπ(t); t] , (13)
and we seek to calculate U as a power series in H I . To this end, we introduce an "interaction picture": we define fluctuation operators δφ I (t) and δπ I (t) whose time dependence is generated by the quadratic part of the Hamiltonian:

δ φI (t) = i H 0 [δφ I (t)
and the initial conditions δφ I (t 0 ) = δφ(t 0 ), δπ I (t 0 ) = δπ(t 0 ). Because H 0 is quadratic, the interaction picture operators are free fields, satisfying linear wave equations. It follows from Eq. ( 14) that in evaluating H 0 [δφ I , δπ I ; t] we can take the time argument of δφ I and δπ I to have any value, and in particular we can take it as t 0 , so that H 0 [δφ I (t), δπ I (t); t] = H 0 [δφ(t 0 ), δπ(t 0 ); t] , but the intrinsic time-dependence of H 0 still remains. The solution of Eq. ( 14) can again be written as a unitary transformation:

δφ I a (t) = U -1 0 (t, t 0 )δφ a (t 0 )U 0 (t, t 0 ) , δπ I a (t) = U -1 0 (t, t 0 )δπ a (t 0 )U 0 (t, t 0 ) , (15)
with U 0 defined by the differential equation

d dt U 0 (t, t 0 ) = -i H 0 [δφ(t 0 ), δπ(t 0 ); t] U 0 (t, t 0 ) ( 1 6 )
and the initial condition U 0 (t 0 , t 0 ) = 1 . Then from Eqs. ( 12) and ( 16) we have

d dt U -1 0 (t, t 0 )U (t, t 0 ) = -iU -1 0 (t, t 0 )H I [δφ(t 0 ), δπ(t 0 ); t]U (t, t 0 ) .
This gives

U (t, t 0 ) = U 0 (t, t 0 )F (t, t 0 ) , d dt F (t, t 0 ) = -iH I (t)F (t, t 0 ) , F(t 0 , t 0 ) = 1 . ( 17
)
where H I (t) is the interaction Hamiltonian in the interaction picture:

H I (t) ≡ U 0 (t, t 0 )H I [δφ(t 0 ), δπ(t 0 ); t]U -1 0 (t, t 0 ) = H I [δφ I (t), δπ I (t); t] ( 1 8 )
The solution of equations like ( 17) is well known (see for example [START_REF] Peskin | ????[END_REF])

F (t, t 0 ) = T exp -i t t0 H I (t) dt , (19)
where T indicates that the products of H I s in the power series expansion of the exponential are to be timeordered; that is, they are to be written from left to right in the decreasing order of time arguments. The solution for the fluctuations in terms of the free fields of the interaction picture is given by Eqs. (11) and (19). Then expectation values of some product A(t) of field operators are obtained as

A(t) = T exp i t t0 H I (t) dt A I (t) T exp -i t t0 H I (t) dt , (20)
where A(t) is any δφ(x, t) or δπ(x, t) or any product of the δφs and/or δπs, all at the same time t but in general with different space coordinates, and A I (t) is the same product of δφ I (x, t) and/or δπ I (x, t). Also, T denotes anti-time-ordering: products of H I s in the power series expansion of the exponential are to be written from left to right in the increasing order of time arguments. It is more convenient to use a formula equivalent to Eq. ( 20):

A(t) = ∞ N =0 i N t -∞ dt N tN -∞ dt N -1 • • • t2 -∞ dt 1 H I (t 1 ), H I (t 2 ), • • • H I (t N ), A I (t) • • • . ( 21
)
In the following section, we apply this formula to calculate the bispectrum of quantum fields generated during inflation.

III. NON-LINEAR COSMOLOGICAL PERTURBATIONS

In this section, we calculate the action for higher order cosmological perturbations that is necessary to compute the interaction Hamiltonian at the third order in the in-in formalism. This calculation is pioneered by Ref. [3] and extended to general inflation models by Ref. [18,19,24]. Here we review the derivation of the higher order action by following Ref. [24,29].

A. Inflation models

To set up our notation, let us first review the formalism in [START_REF] Garriga | [END_REF] where a general Lagrangian for the inflaton field is considered. The Lagrangian is of the general form

S = 1 2 d 4 x √ -g M 2 pl R + 2P (X, φ) , ( 22
)
where φ is the inflaton field and

X = -(1/2)g µν ∂ µ φ∂ ν φ. The reduced Planck mass is M pl = (8πG) -1/2
and the signature of the metric is (-1, 1, 1, 1). The energy of the inflaton field is

E = 2XP ,X -P , (23)
where P ,X denote the derivative with respect to X. Suppose the universe is homogeneous with a Friedmann-Robertson-Walker metric

ds 2 = -dt 2 + a 2 (t)dx 2 3 . (24)
Here a(t) is the scale factor and H = ȧ/a is the Hubble parameter of the universe. The equations of motion of the gravitational dynamics are the Friedmann equation and the continuity equation

3M 2 pl H 2 = E , Ė = -3H(E + P ) . (25)
It is useful to define the "speed of sound" c s as (see [74,75] for a definition of the sound speed)

c 2 s = P ,X E ,X = P ,X P ,X + 2XP ,XX , (26)
and some "slow variation parameters" as in standard slow roll inflation

= - Ḣ H 2 = XP ,X M 2 pl H 2 , η = ˙ H , s = ċs c s H . ( 27
)
These parameters are more general than the usual slow roll parameters (which are defined through properties of a flat potential, assuming canonical kinetic terms), and in general depend on derivative terms as well as the potential. For example, in Dirac-Born-Infeld (DBI) inflation the potential can be steep, and kinetically driven inflation can occur even in absence of a potential. We also note that the smallness of the parameters , η, s does not imply that the rolling of inflaton is slow. When we refer to the slow-roll expansion, we assume that all the three slow variation parameters are small. The primordial power spectrum is derived for this general Lagrangian in [START_REF] Garriga | [END_REF]

P ζ (k) = 1 36π 2 M 4 pl E 2 c s (P + E) = 1 8π 2 M 2 pl H 2 c s , (28)
where the expression is evaluated at the time of horizon exit at c s k = aH and

P ζ (k) = k 3 P ζ (k)/2π 2 . The spectral index is n s -1 = d ln P ζ (k) d ln k = -2 -η -s . (29)
In order to have an almost scale invariant power spectrum, we need to require the 3 parameters , η, s to be very small, which we will denote simply as O( ). We note that in inflationary models with standard kinetic terms the speed of sound is c s = 1. In the case of DBI inflation, the speed of sound can be very small. In the case of arbitrary c s , Eqs. ( 28) and ( 29) for the power spectrum and its index at leading order is still valid as long as the variation of the sound speed is slow, namely s 1. In the following we set M pl = 1.

B. Effective action for higher order perturbations

Now in the general setup described by the action (22), we need to expand the action up to the cubic order in perturbations to obtain the third order interacting Hamiltonian. For this purpose,it is convenient to use the ADM metric formalism [76]. The ADM line element reads

ds 2 = -N 2 dt 2 + h ij dx i + N i dt dx j + N j dt , ( 30
)
where N is the lapse function, N i is the shift vector and h ij is the 3D metric. The action (22) becomes

S = 1 2 dtd 3 x √ hN (3) R + 2P + 1 2 dtd 3 x √ hN -1 E ij E ij -E 2 . ( 31
)
The tensor E ij is defined as

E ij = 1 2 ḣij -∇ i N j -∇ j N i , (32)
and it is related to the extrinsic curvature by K ij = N -1 E ij . ∇ i is the covariant derivative with respect to h ij and all contra-variant indices in this section are raised with h ij unless stated otherwise. The Hamiltonian and momentum constraints are respectively

(3) R + 2P -2π 2 N -2 P ,X -N -2 E ij E ij -E 2 = 0, ∇ j N -1 E j i -∇ i N -1 E = πN -1 ∇ i φP ,X , (33)
where π is defined as

π ≡ φ -N j ∇ j φ. (34)
We decompose the shift vector N i into scalar and intrinsic vector parts as

N i = Ñi + ∂ i ψ, (35)
where ∂ i Ñ i = 0, here indices are raised with δ ij . Before we consider perturbations around our background let us count the number of degrees of freedom (dof) that we have. There are five scalar functions, the field φ, N , ψ, deth and h ij ∼ ∂ i ∂ j H, where H is a scalar function and deth denotes the determinant of the 3D metric. Also, there are two vector modes Ñ i and h ij ∼ ∂ i χ j , where χ j is an arbitrary vector. Both Ñ i and χ j satisfy a divergenceless condition and so carry four dof. Furthermore, we also have a transverse and traceless tensor mode γ ij that contains two additional dof. Because our theory is invariant under change of coordinates we can eliminate some of these dof. For instance, a spatial reparametrization like

x i = xi + ∂ i ˜ (x, t) + i (t) (x, t)
, where ˜ and i (t) are arbitrary and ∂ i i (t) = 0, can be chosen so that it removes one scalar dof and one vector mode. A time reparametrization would eliminate another scalar dof. Constraints in the action will eliminate further two scalar dof and a vector mode. In the end we are left with one scalar, zero vector and one tensor modes that correspond to three physical propagating dof. In this review, we are primarily interested in a scalar degree of freedom.

In order to identify this scalar degree of freedom, we need to fix a gauge. There are two commonly used gauges. In the next subsection we derive the higher order action in these gauges.

C. Non-linear perturbations in the comoving gauge

In the comoving gauge, the scalar degree of freedom is the so-called curvature perturbation ζ and the inflaton fluctuations vanish. The 3D metric is perturbed as

δφ = 0, h ij = a 2 e 2ζ ĥij , ĥij = δ ij + γ ij + 1 2 γ ik γ k j + • • • (36)
where det ĥ = 1, γ ij is a tensor perturbation that we assume to be a second order quantity, i. e.

γ ij = O(ζ 2 ).
It obeys the traceless and transverse conditions γ i i = ∂ i γ ij = 0 (indices are raised with δ ij ). ζ is the gauge invariant scalar perturbation. In (36), we have ignored the first order tensor perturbations (1) γ ij GW . This is because any correlation function involving this tensor mode will be smaller than a correlation function involving only scalars, see results of [3].

We expand N and N i in power of the perturbation ζ

N = 1 + α 1 + α 2 + • • • , (37)

Ñi = Ñi

(1)

+ Ñi

(2)

+ • • • , (38)
ψ = ψ 1 + ψ 2 + • • • , (39)
where α n , Ñi

(n)
and ψ n are of order ζ n . In order to compute the effective action to order O(ζ 3 ), as pointed out in [3], in the ADM formalism one only needs to consider the perturbations of N and N i to the first order O(ζ). This is because their perturbations at order O(ζ 3 ) such as α 3 will multiply the constraint equation at the zeroth order O(ζ 0 ) which vanishes, and the second order perturbations such as α 2 will multiply a factor which vanishes by the first order solution. So the first order solution for N and N i is enough for our purpose. Therefore our task is simplified. In order to expand the action (22) to quadratic and cubic order in the primordial scalar perturbation ζ, we only need to plug in the solution for the first order perturbation in N and N i and do the expansion. Now, the strategy is to solve the constraint equations for the lapse function and shift vector in terms of ζ and then plug in the solutions in the expanded action up to third order. At first order in ζ, a particular solution for equations ( 33) is [3,18]:

α 1 = ζ H , Ñi (1)
= 0,

ψ 1 = - ζ H + χ, ∂ 2 χ = a 2 c 2 s ζ. ( 40
)
The second order action is

S 2 = dtd 3 x a 3 c 2 s ζ2 -a (∂ζ) 2 . ( 41
)
The third order action is [3,18,24]

S 3 = dtd 3 x -aζ (∂ζ) 2 -a 3 (Σ + 2λ) ζ3
H 3 + 3a 3 c 2 s ζ ζ2 + 1 2a 3ζ - ζ H ∂ i ∂ j ψ 1 ∂ i ∂ j ψ 1 -∂ 2 ψ 1 ∂ 2 ψ 1 - 2 a ∂ i ψ 1 ∂ i ζ∂ 2 ψ 1 . (42)
Here we defined two parameters following [18]

Σ = XP ,X + 2X 2 P ,XX = H 2 c 2 s , (43)
λ = X 2 P ,XX + 2 3 X 3 P ,XXX . (44)

D. Non-linear perturbations in the uniform curvature gauge

In this gauge, the inflaton perturbation does not vanish and the 3D metric takes the form

φ(x, t) = φ 0 + δφ(x, t), h ij = a 2 ĥij , ĥij = δ ij + γij + 1 2 γik γk j + • • • (45)
where det ĥ = 1 and γij is a tensor perturbation that we assume to be a second order quantity, i.e., γij = O(δφ 2 ). It obeys the traceless and transverse conditions γi i = ∂ i γij = 0 (indices are raised with δ ij ). We expand N and N i in powers of the perturbation δφ(x, t)

N = 1 + α 1 + α 2 + • • • , (46)

Ñi = Ñi

(1)

+ Ñi

(2)

+ • • • , ( 47
) ψ = ψ 1 + ψ 2 + • • • , (48)
where α n , Ñi

and ψ n are of order δφ n and φ 0 (t) is the background value of the field. At first order in δφ, a particular solution for equations (33) is [3,77]:

α 1 = 1 2H φ0 δφP ,X , Ñi (1)
= 0,

∂ 2 ψ 1 = a 2 c 2 s d dt - H φ δφ . (49)
The second-order action is given by

S 2 = dtd 3 xa 3 P ,XX X 0 δ φ 2 + 2X 0 α 2 1 -2 φ0 α 1 δ φ + P ,Xφ φ0 δφ δ φ -2X 0 α 1 δφ + 1 2 P ,φφ δφ 2 +P ,X 1 2 δ φ 2 -φ0 α 1 δ φ + X 0 α 2 1 -a -2 1 2 (∂δφ) 2 + φ0 ∂ i δφ∂ i ψ 1 -3H 2 α 2 1 + P ,φ α 1 δφ -2a -2 Hα 1 ∂ 2 ψ 1 , (50)
where X 0 = φ2 0 /2. The third-order action is obtained as

S 3 = dtd 3 xa 3 P ,XX 1 2 φ0 δ φ 3 + X 0 α 1 -4 δ φ 2 + 5 φ0 α 1 δ φ -4X 0 α 2 1 +a -2 - 1 2 φ0 δ φ(∂δφ) 2 + X 0 α 1 (∂δφ) 2 -2X 0 δ φ -φ0 α 1 ∂ i δφ∂ i ψ 1 +P ,Xφ 1 2 δφ δ φ 2 -φ0 α 1 δφ δ φ + X 0 α 2 1 δφ -a -2 1 2 δφ(∂δφ) 2 + φ0 δφ∂ i δφ∂ i ψ 1 +P ,XXX X 0 1 3 φ0 δ φ 3 + X 0 α 1 -2 δ φ 2 + 2 φ0 α 1 δ φ - 4 3 X 0 α 2 1 +P ,XXφ X 0 δφ δ φ 2 -2 φ0 α 1 δφ δ φ + 2X 0 α 2 1 δφ + P ,Xφφ 1 2 φ0 δ φ -X 0 α 1 δφ 2 +P ,X α 1 - 1 2 δ φ 2 + φ0 α 1 δ φ -X 0 α 2 1 -a -2 1 2 α 1 (∂δφ) 2 + δ φ -φ0 α 1 ∂ i δφ∂ i ψ 1 + 1 2 P ,φφ α 1 δφ 2 + 1 6 P ,φφφ δφ 3 + 3H 2 α 3 1 + 2a -2 Hα 2 1 ∂ 2 ψ 1 + 1 2 a -4 α 1 (∂ 2 ψ 1 ) 2 -∂ i ∂ j ψ 1 ∂ i ∂ j ψ 1 . (51)

E. Relation between gauges

The gauges used in the previous two sections are of course related by a gauge transformation. Introducing a new variable ζ n defined by ζ n = -Hδφ/ φ0 , ζ in the comoving gauge is related to δφ in the flat gauge as [3]

ζ = ζ n + f (ζ n ), (52)
where

f (ζ) = η 4c 2 s ζ 2 + 1 c 2 s H ζ ζ + 1 4a 2 H 2 -(∂ζ)(∂ζ) + ∂ -2 (∂ i ∂ j (∂ i ζ∂ j ζ)) + 1 2a 2 H (∂ζ)(∂χ) -∂ -2 (∂ i ∂ j (∂ i ζ∂ j χ)) .
(53) On large scales where ζ n becomes constant and we get

ζ = ζ n + η 4c 2 s ζ 2 n . (54)
We can show that this is nothing but the expression for ζ obtained in the delta-N formalism using the relations

N ,φ = - H φ0 , N ,φφ = Ḣ φ2 0 - φ0 H φ3 0 , (55)
and the definition of η in Eq. ( 27). There are two ways to calculate the bispectrum of ζ on large scales. One is to calculate the bispectrum of ζ n = -Hδφ/ φ0 in the flat gauge and apply the delta-N formalism. It is also possible to calculate the bispectrum of ζ in the comoving gauge. We will use both approaches in the next section.

IV. BISPECTRUM OF CURVATURE PERTURBATION

In this section, we first calculate the three point function for the field perturbations in the in-in formalism using the cubic order action obtained in the previous section. We only consider leading order terms in slow-roll expansions. In the following, we consider two inflation models, k-inflation and standard slow-roll inflation. Then using the delta-N formalism, we derive the bispectrum of the curvature perturbation. We follow the calculations in Ref. [19,30,31]. We also discuss a method to calculate it directly in the comoving gauge presented in Ref. [24].

A. Bispectrum of quantum fields in k-inflation

First let us consider models with non-standard kinetic terms. This is known as k-inflation. In this case, the leading order terms in the slow-roll expansion in the action in the flat gauge ( 50) and ( 51) are given by

S 2 = dtd 3 x a 3 P ,X 2
1 c 2 s δ φ 2 - 1 a 2 (∂δφ) 2 , (56)
S 3 = dtd 3 x P ,XX φ 2 + P ,XXX φ3 6 a 3 δ φ 3 -P ,XX φ 2 a δ φ(∂δφ) 2 . ( 57
)
The perturbations in the interacting picture are promoted to quantum operators like

δφ(τ, x) = 1 (2π) 3 d 3 kδφ(τ, k)e ik•x , δφ(τ, k) = u(τ, k)a(k) + u * (τ, -k)a † (-k). (58)
a(k) and a † (-k) are the annihilation and creation operator respectively, that satisfy the usual commutation relations

a(k 1 ), a † (k 2 ) = (2π) 3 δ (3) (k 1 -k 2 ), [a(k 1 ), a(k 2 )] = a † (k 1 ), a † (k 2 ) = 0. ( 59
)
At leading order the solution for the mode functions is given by

u(τ, k) = H 2c s P ,X 1 k 3/2 (1 + ikc s τ ) e -ikcsτ . (60)
Using the in-in formalism Eq. ( 21), the vacuum expectation value of the three point operator in the interaction picture is written as [3,71]

δφ(t, k 1 )δφ(t, k 2 )δφ(t, k 3 ) = -i t t0 d t δφ(t, k 1 )δφ(t, k 2 )δφ(t, k 3 ), H I ( t) , (61)
where t 0 is some early time during inflation when the field's vacuum fluctuation are deep inside the horizons, t is some time after horizon exit. If one uses conformal time, it's a good approximation to perform the integration from -∞ to 0 because τ ≈ -(aH) -1 . H I denotes the interaction Hamiltonian and it is given by H I = -L 3 , where L 3 is the Lagrangian obtained from the action (57). using the solution for the mode function and commutation relations for the creation and annihilation operators, we get

δφ(k 1 )δφ(k 2 )δφ(k 3 ) = -(2π) 3 δ (3) (k 1 + k 2 + k 3 ) H 4 √ 2 c 2 s (P ,X ) 3/2 1 Π 3 i=1 k 3 i A k-inf φ (k 1 , k 2 , k 3 ), (62)
where

A k-inf φ = - 3λ Σ k 2 1 k 2 2 k 2 3 K 3 + 1 c 2 s -1 k 2 1 k 2 • k 3 K 1 + k 2 + k 3 K + 2 k 2 k 3 K 2 + 2 cyclic terms (63) = 1 c 2 s -1 - 2λ Σ 3k 2 1 k 2 2 k 2 3 K + 1 -c 2 s c 2 s   - 1 K i>j k 2 i k 2 j + 1 2K 2 i =j k 2 i k 3 j + 1 8 i k 3 i   . (64)

B. Bispectrum of quantum fields in slow-roll inflation

In slow-roll inflation with a standard kinetic term P (X) = X -V (φ), leading-order terms in the third-order action (51) are given by

S 3 = dt dx 3 a 3 - 1 a 2 δ φ∂ψ∂δφ - 1 4H φδφ( δ φ) 2 - 1 a 4 1 4H φδφ(∂δφ) 2 ,
where ψ was defined in Eq. ( 49) and is given to leading order by

∂ 2 ψ = - a 2 2H φ δ φ. ( 65
)
The three point function can be calculated in the same way. We obtain

δφ(k 1 )δφ(k 2 )δφ(k 3 ) = -(2π) 3 δ (3) (k 1 + k 2 + k 3 ) H 4 √ 2 1 Π 3 i=1 k 3 i A stand φ (k 1 , k 2 , k 3 ), ( 66
)
where

A stand φ (k 1 , k 2 , k 3 ) =   - 1 8 i k 3 i + 1 8 i =j k i k 2 j + 1 K i>j k 2 i k 2 j   . (67)

C. Bispectrum of curvature perturbation

Now we can apply the delta-N formalism to calculate the bispectrum of the curvature perturbation. We define the bispectrum of the curvature perturbation as

ζ(k 1 )ζ(k 2 )ζ(k 3 ) = (2π) 7 δ 3 (k 1 + k 2 + k 3 )(P ζ ) 2 1 i k 3 i A ζ (k 1 , k 2 , k 3 ) , ( 68
)
where P ζ is given by Eq. ( 28). In k-inflation models, in the small sound speed limit and at the leading order in slow-roll expansion, the relation between the curvature perturbation and the field perturbation is simply given by ζ = -Hδφ/ φ = P ,X δφ/2 . Then the three point function for ζ is given by Eq. ( 68) where

A ζ = A k-inf φ
given by Eq. ( 64). In standard slow-roll inflation case, the relation between ζ and ζ n can be written as ζ = ζ n + ηζ 2 n /4. Then the the bispectrum of the curvature perturbation is given by Eq. ( 68) where A ζ is given by

A ζ = A stand φ + η 1 8 i k 3 i . (69)

D. Computation in comoving gauge

It is also possible to calculate the bispectrum of the curvature perturbation in comoving gauge and this gives a very useful consistency check. Here we follow Ref. [24] and see how this works.

In fact the cubic effective action in (57) looks like order O( 0 ) in the slow variation parameters while in the previous section, we find that the bispectrum is suppressed by slow-roll parameters in slow-roll inflation.

In slow-roll inflation, as emphasized and demonstrated in Ref. [3], one can perform a lot of integrations by parts and cancel terms of order O( 0 ) and O( ). The resulting cubic action is actually of leading order O( 2) in slow roll parameters. A similar analysis can be performed for the general Lagrangian in Ref. [18]. Except for terms that are proportional to 1 -c 2 s or λ, the rest of the terms can be cancelled to the second order O( 2) and the cubic order action Eq. ( 42) can be rewritten as

S 3 = dtd 3 x -a 3 Σ(1 - 1 c 2 s ) + 2λ ζ3 H 3 + a 3 c 4 s ( -3 + 3c 2 s )ζ ζ2 + a c 2 s ( -2s + 1 -c 2 s )ζ(∂ζ) 2 -2a c 2 s ζ(∂ζ)(∂χ) + a 3 2c 2 s d dt η c 2 s ζ 2 ζ + 2a (∂ζ)(∂χ)∂ 2 χ + 4a (∂ 2 ζ)(∂χ) 2 + 2f (ζ) δL δζ 1 , ( 70
)
where χ is defined in Eq. ( 40) and in the last term

δL δζ 1 = a d∂ 2 χ dt + H∂ 2 χ -∂ 2 ζ . ( 71
)
Here ∂ -2 is the inverse Laplacian, δL/δζ| 1 is the variation of the quadratic action with respect to the perturbation ζ, therefore the last term which is proportional to δL/δζ| 1 can be absorbed by a field redefinition of ζ. It can be easily shown that the field redefinition that absorbs this term is

ζ → ζ n + f (ζ n ), ( 72
)
where f (ζ) is given by Eq. ( 53). This is nothing but the relation between ζ and ζ n obtained from the guage transformation between the flat gauge and comoving guage. One then computes the vacuum expectation value of the three point function in the interaction picture in the same way. We get Eq. ( 68) with

A = 1 c 2 s -1 - 2λ Σ 3k 2 1 k 2 2 k 2 3 2K 3 + 1 c 2 s -1   - 1 K i>j k 2 i k 2 j + 1 2K 2 i =j k 2 i k 3 j + 1 8 i k 3 i   + c 2 s   - 1 8 i k 3 i + 1 8 i =j k i k 2 j + 1 K i>j k 2 i k 2 j   + η c 2 s 1 8 i k 3 i + s c 2 s   - 1 4 i k 3 i - 1 K i>j k 2 i k 2 j + 1 2K 2 i =j k 2 i k 3 j   . (73)
In k-inflation, the first two terms are the leading order contributions in slow-roll expansions. The remaining terms are O( ). Note that one should take into account O( ) corrections from the leading order contribution in order to obtain a full expression up to O( ) in k-inflation. In standard slow-roll inflation, the first two terms and the last term vanishes and c s = 1.

V. SHAPES OF BISPECTRUM

In this section, we compare the prediction of the bispectrum in k-inflation with the local-type non-Gaussianity. The discussions in this section are based on Refs. [78,79].

The bispectrum in the local-type non-Gaussianity is often characterized by

ζ = ζ n + 3 5 f local NL (ζ 2 n -ζ 2 n ), (74)
where ζ n obeys Gaussian statistics. Originally the f local NL parameter was introduce to parametrise a nonlinearity in the curvature perturbation Φ in the Longitudinal gauge which is related to ζ as Φ = (3/5)ζ. Note that Ref. [3] uses a different sign convention for f local NL from WMAP papers (see [80] for the latest result). Here we follow the definition used in WMAP papers. Eq. ( 74) is nothing more than the expression for ζ in the delta-N formalism. In slow-roll inflation f local NL is O( ) but models like curvaton predicts f local NL larger than one. The bispectrum of curvature perturbation is given by

ζ(k 1 )ζ(k 2 )ζ(k 3 ) = (2π) 3 δ 3 (k 1 + k 2 + k 3 )(P ζ ) 2 F (k 1 , k 2 , k 3 ), (75)
where

F local (k 1 , k 2 , k 3 ) = (2π) 4 3 10 f local N L 1 k 3 1 k 3 2 + 1 k 3 2 k 3 3 + 1 k 3 3 k 3 1 . (76)
Eq. ( 74) describes (at leading order) the most generic form of non-Gaussianity which is local in real space. This form is therefore expected for models where non-linearities develop outside the horizon. This happens for all the models in which the fluctuations of an additional light field, different from the inflaton, contribute to the curvature perturbations we observe. In this case non-linearities come from the evolution of this field outside the horizon and from the conversion mechanism which transforms the fluctuations of this field into the curvature perturbations. Both these sources of non-linearity give a non-Gaussianity of the form (74) because they occur outside the horizon. Examples of this general scenario are the curvaton models [81], models with fluctuations in the reheating efficiency [8,82] and multi-field inflationary models [6].

Being local in position space, Eq. ( 74) describes correlation among Fourier modes of very different k. It is instructive to take the limit in which one of the modes becomes of very long wavelength [3], k 3 → 0, which implies, due to momentum conservation, that the other two k's become equal and opposite. The long wavelength mode ζ k3 freezes out much before the others and behaves as a background for their evolution. In this limit F local is proportional to the power spectrum of the short and long wavelength modes

F local ∝ 1 k 3 3 1 k 3 1 . ( 77
)
This means that the short wavelength 2-point function ζ k1 ζ -k1 depends linearly on the long wavelength mode

ζ k3 ζ k3 ζ k1 ζ -k1 ∝ ζ k3 ζ -k3 ∂ ∂ζ k3 ζ k1 ζ -k1 . (78)
From this point of view we expect that any bispectra will reduce to the local shape (76) in the degenerate limit we considered if the derivative with respect to the long wavelength mode does not vanish.

In standard single field slow-roll inflation, as pointed out in [3], different points along the background wave are equivalent to shift in time along the inflaton trajectory, so that the derivative with respect to the background wave is proportional to the tilt of the scalar spectrum. This can be explicitly checked in the full expression of the 3-point function (Eqs. ( 67) and ( 69)):

F stand (k 1 , k 2 , k 3 ) = (2π) 4 i k 3 i     - 1 8 i k 3 i + 1 8 i =j k i k 2 j + 1 K i>j k 2 i k 2 j   + η 1 8 i k 3 i   , (79)
In the limit k 3 → 0 Eq. ( 79) goes as

F stand (k 3 → 0) ∝ (2 + η) 1 k 3 3 k 3 1 = -(n s -1) 1 k 3 3 k 3 1 . ( 80
)
As expected the tilt in the spectrum n s fixes the degenerate limit of the 3-point function. Note however that expression (79) is not of the local form (76) but contains contributions which are important for nondegenerate triangles. If we compare expression ( 76) and ( 79) and neglect the different shape dependence, we see that standard single-field inflation predicts f local NL of order of the slow-roll parameters. We have seen that the degenerate limit k 3 → 0 describes the effect of a slowly-varying long-wavelength perturbation on the 2-point function of short wavelength modes. In many models, the correlation is much weaker in this limit than in the local model (76). Physically this means that the correlation is among modes with comparable wavelength which go out of the horizon nearly at the same time. In this case the 3-point function in the degenerate limit is suppressed by powers of k 3 with respect to the behaviour of Eq. ( 77). We have correlation among modes of comparable wavelength in all models in which the non-Gaussianity is generated by derivative interactions: these interactions become exponentially irrelevant when the modes go out of the horizon because both time and spatial derivatives become small, so that all the correlation is among modes freezing almost at the same time.

FIG. 1: Left: Plot of the function F (1, k2/k1, k3/k1)(k2/k1) 2 (k3/k1) 2 predicted by the DBI models [10]. Right: Difference between the above plot and the analogous one (top of fig. 2) for the factorizable equilateral shape used in the analysis. From [79].

K-inflation is a typical example for these type of models. The three point function is obtained in Eq. ( 73) and given by

F (k 1 , k 2 , k 3 ) = (2π) 4 i k 3 i   1 c 2 s -1 - 2λ Σ 3k 2 1 k 2 2 k 2 3 2K 3 + 1 c 2 s -1   - 1 K i>j k 2 i k 2 j + 1 2K 2 i =j k 2 i k 3 j + 1 8 i k 3 i     . (81)
In a model of inflation based on the DBI action,

P (φ, X) = -f (φ) -1 1 -2Xf (φ) + f (φ) -1 + V (φ), (82)
λ is given by

λ = Σ 2 1 c 2 s -1 . ( 83
)
Thus the first term vanishes in Eq. (81). Unfortunately, a function F is not factorizable, so it is not easy to perform an optimal analysis using CMB observations. However, it is a very good approximation to take a factorizable shape function F which is close to Eq. ( 81) and perform the analysis for this shape. In the limit k 1 → 0 with k 2 and k 3 fixed, all the equilateral functions diverge as k -1 1 [78] (while the local form eq. ( 76) goes as k -3 1 ). The factorizable function that satisfies this condition is given by

F (k 1 , k 2 , k 3 ) = (2π) 4 9 10 f equil NL - 1 k 3 1 k 3 2 - 1 k 3 1 k 3 3 - 1 k 3 2 k 3 3 - 2 k 2 1 k 2 2 k 2 3 + 1 k 1 k 2 2 k 3 3 + (5 perm.) , (84)
where the permutations act only on the last term in parentheses. In figure 1 we study the equilateral function predicted in DBI inflation (81). In the second part of the figure we show the difference between this function and the factorizable one used in our analysis. We see that the relative difference is quite small. The same remains true for other equilateral shapes (see [78] for the analogous plots for other models).

In figure 2, we compare this function with the local shape. The dependence of both functions under a common rescaling of all k's is fixed to be ∝ k -6 by scale invariance, so that we can factor out k -6 1 for example. Everything will now depend only on the ratios k 2 /k 1 and k 3 /k 1 , which fix the shape of the triangle in momentum space. For each shape we plot F (1, k 2 /k 1 , k 3 /k 1 )(k 2 /k 1 ) 2 (k 3 /k 1 ) 2 ; this is the relevant quantity if we are interested in the relative importance of different triangular shapes. The square of this function gives the signal to noise contribution of a particular shape in momentum space [78]. We see that for the function (84), the signal to noise is concentrated on equilateral configurations, while squeezed triangles with one side much smaller than the others are the most relevant for the local shape.

FIG. 2: Plot of the function F (1, k2/k1, k3/k1)(k2/k1) 2 (k3/k1) 2 for the equilateral shape used in the analysis (left) and for the local shape (right). The functions are both normalized to unity for equilateral configurations

k 2 k 1 = k 3 k 1 = 1. Since F (k1, k2, k3) is symmetric in its three arguments, it is sufficient to specify it for k1 ≥ k2 ≥ k3, so k 3 k 1 ≤ k 2 k 1 ≤ 1 above.
Moreover, the triangle inequality says that no side can be longer than the sum of the other two, so we only plot F in the triangular region 1

-k 2 k 1 ≤ k 3 k 1 ≤ k 2 k 1 ≤ 1
above, setting it to zero elsewhere. From [79].

VI. CONCLUSION

In this review, the bispectrum of curvature perturbation is calculated using the in-in formalsim and the delta-N formalism. There are two distinct contributions to the bispectrum. One is coming from a non-linear relation between the curvature perturbation ζ and quantum fluctuations of a field at the horizon crossing. In this case the non-linearities come from the evolution of this field outside the horizon. Being local in a position space, the shape of the bispectrum is highly non-local in a Fourie space having a maximum signal for the squeezed configuration k 3 k 1 , k 2 . The other contribution is coming from the bispectrum of quantum fields generated under horizon scales. In models in which the non-Gaussianity is generated by derivative interactions such as DBI inflation and k-inflation models, we have correlation among modes of comparable wavelength in all models and these interactions become exponentially irrelevant when the modes go out of the horizon because both time and spatial derivatives become small, so that all the correlation is among modes freezing almost at the same time. Then the bispectrum has a peak at the equilateral configuration

k 1 ∼ k 2 ∼ k 3 .
In order to put constraints on the bispectrum from CMB observations, it is necessary to construct an estimator that uses a model prediction for the bispectrum as an template. For the local-type non-Gaussianity and the equilateral non-Gaussianity, the constraints obtained in WMAP 7-year results are [80] -10 < f local NL < 74, -214 < f equil NL < 266.

at 95% confidence level. Recently, it has been found that Large Scale Structure (LSS) can give a similar level of constraints on the local type non-Gaussianity from the scale-dependent bias effects on the halo power spectrum [83] while this effect is absent in the equilateral non-Gaussianity [84,85]. The constraint on f local NL is obtained from SDSS as -29 < f local NL < 70 [86] and combining it to WMAP 7-year results, we get -5 < f local NL < 59 [80]. Thus currents observations are consistent with Gaussian primordial curvature perturbations but future experiments such as Planck will give much tighter constraints f local NL ∼ O(1) and we may be able to detect a deviation from Gaussianity which has a huge impact on early universe models.

There are a lot issues that are not covered by this review. We will mention some of the issues here:

• Trispectrum In this review, we concentrated on the leading order non-Gaussianity, i.e. the bispectrum but it has been recognized that the trispectrum could give a useful information to distinguish between many possible models that predict large non-Gaussianity. For the local non-Gaussianity, we can easily extend the model by expanding ζ up to the third order ζ = ζ n + f NL ζ 2 n /2 + g NL ζ 3 n /6. The trispetrum is characterized by two parameters τ NL = f 2 NL and g NL [17,28,77]. The constraints on these parameters are rather weak, -631 < f 2 NL < 717 and -3.80 < g NL /10 6 < 3.88 from WMAP 5-year results at 95% confidence level [87]. But again the Planck will improve these significantly [88]. The full trispectrum in DBI inflation at the leading order in small sound speed limit has been obtained [25,[89][90][91]. Unlike the bispectrum, there are still two degrees of freedom even for the equilateral configurations k 1 ∼ k 2 ∼ k 3 ∼ k 4 and also the form of the trispectrum is too complicated to be used for the estimator.

It is necessary to develop approximations for the shape of the trispectrum in DBI inflation as is done for the bispectrum.

• Multi-field inflation

In single field inflation models with a standard kinetic term, the resulting non-Gaussianity is small suppressed by slow-roll parameters. However, in multi-field models, it is possible to have large local non-Gaussianity due to non-linear dynamics of fields outside horizon. The delta-N formalism is easy to be extended to multi-field models and there have been extensive study of non-Gaussianity in multi-field models using the delta-N formalism (see for example [92] and references therein). Multi field effects are also important in DBI inflation. In DBI inflation, fluctuations along the entropy directions of the fields that are orthogonal to the field trajectory have the same sound speed as the adiabatic fluctuations along the field trajectory [30][31][32]93]. If the trajectory makes a turn in a field space, this converts the entropy perturbations to the curvature perturbations. Although the bispectrum is enhanced by this conversion, the enhancement of the power spectrum is stronger and f equil N L becomes smaller in multi-field models, which help ease stringent constraints on DBI inflation models in string theory [31]. It has been shown that the trispectrum is enhanced for a given f equil N L [94,95]. It is also possible that the multi-field effects modify the bipsectrum for quantum field at the horizon crossing. In the so-called quasi-single inflation models [96], the entropy perturbations develops large non-Gaussianity. The conversion of entropy perturbations to the curvature perturbation can happen near the horizon crossing and during this transition the shape of the bispectrum can be modified in a non-trivial way.

There are many other possibilities to get large non-Gaussianity of quantum fields such as a feature in inflaton potentials [26,97]. All these models predict distinct shapes of the bispectrum and trispectrum. In the future, we may be able to exploit CMB and LSS data to distinguish between many possible early universe models via non-Gaussianity.

Acknowledgments

KK is supported by the UK's Science & Technology Facilities Council, the European Research Council and Research Councils UK.