# SafeCritic: Collision-Aware Trajectory Prediction

Tessa van der Heiden
<EMAIL_ADDRESS>&Naveen Shankar Nagaraja*
<EMAIL_ADDRESS>
&Christian Weiß*
<EMAIL_ADDRESS>&Efstratios Gavves
<EMAIL_ADDRESS>BMW Group, Landshuter Str. 26, 85716 Unterschleißheim, Germany
University of Amsterdam, Science Park 904, 1012 WX Amsterdam, Netherlands

###### Abstract

Navigating complex urban environments safely is a key to realize fully
autonomous systems. Predicting future locations of vulnerable road users, such
as pedestrians and cyclists, thus, has received a lot of attention in the
recent years. While previous works have addressed modeling interactions with
the static (obstacles) and dynamic (humans) environment agents, we address an
important gap in trajectory prediction. We propose SafeCritic, a model that
synergizes generative adversarial networks for generating multiple “real”
trajectories with reinforcement learning to generate “safe” trajectories. The
Discriminator evaluates the generated candidates on whether they are
consistent with the observed inputs. The Critic network is environmentally
aware to prune trajectories that are in collision or are in violation with the
environment. The auto-encoding loss stabilizes training and prevents mode-
collapse. We demonstrate results on two large scale data sets with a
considerable improvement over state-of-the-art. We also show that the Critic
is able to classify the safety of trajectories.

_K_ eywords Pedestrian trajectory prediction $\cdot$ generative adversarial
networks $\cdot$ reinforcement learning

##  1 Introduction

Forecasting future locations of humans is an essential task to develop safe
autonomous agents that interact with people. Autonomous systems such as
household and factory robots, security systems and self driving cars need to
know where people will be given their previous movement to plan safe paths.

Predicting people’s future movement is challenging as the target location of a
person cannot be known. The static and dynamic layout of the environment, be
it a road network or an indoors map complicates the forecasting. Static
obstacles place restrictions beyond their narrow definition of spatial extent
(pedestrians will tend not to walk between closely situated parked cars).
Considering dynamic obstacles like other pedestrians, bicycles or moving cars,
one must not only minimize direct future collisions but also take into account
hidden social rules (keeping to the right, respecting personal space). This
becomes more complex since there are multiple unknown target locations each of
the dynamic obstacles might be heading to. Consequently, predicting all the
possible and plausible future trajectories is a multi-modal problem.
Respecting the hidden rules and addressing the challenges is key for systems
that are not just accurate but also safe.

Computationally tractable generative models has led to a great number of
methods [1, 2, 3] that model the multi-modality of future trajectory
prediction. While these models focus on producing accurate future
trajectories, they do not take into account the safety of the trajectories
explicitly in terms of collisions. As a result, the generated paths might be
unsafe, wherein they collide with other road users, or even infeasible, i.e.,
they go through static scene elements like buildings, cars or trees.
Consequently, these trajectories cannot easily be trusted in practice.

In this work we argue that one should take into account safety, when designing
a model to predict future trajectories. The problem is that future
trajectories are generated online and cannot be predetermined. Thus, applying
standard and differentiable supervised learning is not straightforward. Even
if one could minimize collisions with static elements and the ground truth
trajectories, the model must learn to minimize collisions with the generated,
and therefore unseen, future trajectories. Thus, in the absence of a clear
learning objective, we rely on Generative Adversarial Networks (GAN) [4] for
the generative modelling and cast the learning within the Inverse
Reinforcement Learning framework [5]. Inverse Reinforcement Learning guides
the model towards optimizing an unknown reward function and shares many
similarities with GANs. Here, our unknown reward function accounts for
learning safe and “collision-free” future trajectories that have not even been
generated yet.

This work has three contributions: (a) to the best of our knowledge, we are
the first to propose generating “safe” future trajectories, minimizing
collisions with static and dynamic elements in the scene, (b) the generated
trajectories are not only safer but also more accurate. The reason is that the
model learns to prune out unlikely trajectories and generalize better. In
fact, we show that the remaining few collisions in the trajectories strongly
correlate with ground truth, thus showing that the model learns to also
anticipate collisions. and (c), we rely on attention mechanism for spatial
relations in the scene to propose a compact representation for modelling
interaction among all agents. The representation, which we refer to as
attentive scene representation, is invariant to the number of people in the
scene and thus is scalable.

##  2 Related work

Future trajectory prediction. Predicting future trajectories of road users has
recently gained an increased attention in the literature. In their pioneering
work Alahi _et al_[6] define future trajectory prediction as a time series
modelling problem and proposed SocialLSTM, which aggregates information of
multiple agents in a single recurrent model. While SocialLSTM is a
deterministic model, Gupta _et al_[1] proposed SocialGAN to learn to generate
multi-modal future trajectories. To make sure that novel future trajectories
are sampled, the encoder and decoder were trained using a Conditional
Generative Adversarial Network [7, 8] learning framework, leading to
impressive results.

Extending on SocialGAN, SoPhie [2] incorporate semantic visual features
extracted by deep networks, combined with an attention mechanism. Inspired by
[2], DESIRE [3] also relies on semantic scene context and introduces a ranking
and refinement module to order the randomly sampled trajectories by their
likelihood. Further, [9] develop R2P2, which is a non-deterministic model
trained by modified cross-entropy objective to balance betwee the distribution
and demonstrations. SocialAttention [10] connects multiple RNNs in a graph
structure to efficiently deal with the varying amount of people in a scene.

Future trajectory prediction models beyond the deep learning have also been
explored. The Social Force [11] is one of the first models of pedestrian
trajectory prediction that was successfully applied in the real-world [12].
Interacting Gaussian Process (GP) based methods [13] model the trajectory of
each person as a GP with interaction potential terms that couples each
person’s GP. An extensive array of comparisons is available in [14], based on
the TrajNet [15] benchmark.

Here, we adopt a modified Conditional Generative Adversarial Network framework
[7, 8] to generate novel future trajectories from observed past ones. Compared
to [1, 2, 2, 3] our focus is to generate trajectories that are not just
accurate but also lead to minimum collisions and thus are safe. Safe
trajectories are different from trajectories that try to imitate the ground
truth [9], as the latter may lead to implausible paths, _e.g_ , pedestrians
going through walls.

Generative Adversarial Networks + Reinforcement Learning. To model safe
trajectories we propose a model that combines Generative Adversarial Networks
Generative Adversarial Networks [4] and Inverse Reinforcement Learning [5].
[4] are generative models that rely on the min-max optimization of two
competing sub-models, the generator and the discriminator. As standard GANs do
not have an encoder, Conditional GANs [7, 8] include an encoder, which can be
used to generate novel samples given observations. Finn _et al_[16] show that
GANs resemble Inverse Reinforcement Learning [5] in that the discriminator
learns the cost function and the generator represents the policy. Inspired by
[17, 18] who showed the potentials of using external simulators as critic
objectives [19], we propose to use a critic network to model collision
minimization by Inverse Reinforcement Learning. To prevent mode collapse, an
additional auto-encoding loss is considered.

##  3  SafeCritic

Figure 1: The SafeCritic network consisting of a Generator trained by
Discriminator and Critic  networks. The Generator observes each trajectory
together with the dynamic and static scene and computes a future trajectory.
The Critic and discriminator score the predictions. The orange boxes are the
three different losses for the Generator. The Dynamic Scene Representation are
the hidden states of each neighbor $j$ and embedded into a 3D tensor for each
target agent $i$. The Static Scene Representation are segmented points around
each target. The attention of the person in red is shown on the static and
dynamic scene.

###  3.1 Overview

Our goal is a model that predicts accurate, diverse and safe, _i.e_ collision-
free, future trajectories for all relevant road users in a driving scene. We
denote the past trajectories of all $n$ road users participating in the scene
with $\mathbf{X}=\\{X_{1},X_{2},...,X_{n}\\}$. Each trajectory is a sequence
of vectors $X_{i}=[\overline{x}_{i}^{t}]_{t=-T{obs}}^{0}$ for $T_{obs}$ past
consecutive time steps. Each element is a $\overline{x}_{i}^{t}=(\delta
x_{i}^{t},\delta y_{i}^{t})$, a 2D coordinate displacement vector for agent
$i$, with $\delta x_{i}^{t}=x_{i}^{t+1}-x_{i}^{t}$. We also denote the ground
truth future trajectories also containing 2D displacement vectors as elements
with $\mathbf{Y}=\\{Y_{1},Y_{2},...,Y_{n}\\}$ and the predicted trajectories
with $\mathbf{\hat{Y}}$. The time span of the prediction is
$[T_{0},\ldots,T_{\text{pred}}]$. To reduce notation clutter, we drop the
subscript and superscripts $i$ and $t$ when they can be inferred from the
context.

The proposed SafeCritic uses a conditional Generative Adversarial Network [1],
as well as a Critic Network [20], [21],[22], to generate novel and safe future
trajectories for the road users in the scene, see Figure 1 for an overview.
Next, we describe the particular design of SafeCritic.

###  3.2 Generating Trajectories with SafeCritic

The conditional GAN of SafeCritic comprises of two modules, a _Generator_ ,
and a _Discriminator_ , as seen in the left of Figure 1. The generator of the
conditional GAN is composed of two sub-networks.

The first sub-network is a $\theta$-parameterized encoder
$\mathbf{u}=g^{enc}(\mathbf{X};\theta)$ that receives as input the past
trajectories of all road users and maps them to a common embedding space
$\mathbf{u}$. The purpose of the encoder is to aggregate the current state of
all road users jointly. The aggregated state of all road users is then passed
to the second sub-network, which is a $\omega$-parameterized decoder
$\hat{\mathbf{Y}}=g^{dec}(\mathbf{u},F_{d},F_{s},z;\omega)$, where
$F_{d},F_{s}$ are the dynamic and static aggregated features that will be
detailed in 3.4, and $z\sim N(0,1)$ is the random noise vector for sampling
novel trajectories. The purpose of the decoder is to generate novel future
coordinates $\hat{\mathbf{Y}}$ for all $n$ road users simultaneously.

The generator is followed by the discriminator, which is implemented as a
third $\phi$-parameterized sub-network
$D=d(\mathbf{X},\mathbf{Y},\mathbf{\hat{Y}}|\phi)$. The discriminator receives
as inputs the past trajectories $\mathbf{X}$, as well as the ground truth
future trajectories $\mathbf{Y}$ and the predicted trajectories
$\mathbf{\hat{Y}}$. Being adversaries, the generator and the discriminator are
trained in an alternative manner [4]: the generator generates gradually more
realistic trajectories, while the discriminator classifies which trajectories
are fake thus motivating the generator to become better.

###  3.3 Critic for Obtaining Safe Future Trajectories

While the generator and the discriminator modules are trained to generate
plausible future trajectories, the proposed trajectories are no provisions to
avoid possible collisions in the predicted future trajectories. To this end,
we incorporate a Critic module responsible for minimizing the static and
dynamic future collisions, as shown in the right of Figure 3.

A natural objective for the Critic module is a supervised loss that minimizes
the number of future dynamic and static collisions. Unfortunately, integrating
such an objective to the model is not straightforward. The number of non-
colliding trajectories existing in data sets is overwhelmingly higher than for
colliding trajectories. Moreover, having a generative model means that the
generated trajectories will not necessarily imitate the ground truth ones.

To this end, we rely on the off-policy actor-critic framework [23] to define
the objective of Critic . The actor is the trajectory generator decoder
$g^{dec}(\cdot)$, taking a displacement action $\hat{\mathbf{Y}}$. The Critic
then is a $\psi$-parameterized sub-network
$v=o(\mathbf{\hat{Y}},F_{s},F_{d};\psi)$ that predicts how likely an eminent
collision is or not, given $\hat{\mathbf{Y}}$. The action value $v$, is the
score for the future action causing a collision. Another choice would have
been to simply return the reward signal to the generator. However, in that
case we would need to rely on high-variance gradients returned by REINFORCE.

Generator For training SafeCritic, we must take into account the opposing
goals of the generator and the discriminator, as well as the regularizing
effect of the critic. Following standard GAN procedure we arrive at a min-max
training objective for the generator and the discriminator [4]. As the
discriminator quickly learns to tell the real and generated samples apart,
training the generator is harder than training the discriminator, leading to
slow convergence or even failure [24].

Inspired from [1], we introduce an additional auto-encoding loss to the ground
truth for the best generated sample, $\mathbf{\hat{Y}}$,
$\operatorname*{arg\,min}_{\theta,\phi,\omega}\sum_{i}\left\lVert{\mathbf{Y}_{i}-\mathbf{\hat{Y}}_{i_{K}}}\right\rVert_{2}^{2}$.
An additional benefit of the auto-encoding loss is that it encourages the
model to avoid trivial solutions and mode collapse, as well as to increase the
diversity of future generated trajectories.

Critic For the training of the Critic we employ a version of the Deep
Deterministic Policy Gradient algorithm [25], more specifically, the DDPG
variant [26] that uses a sequence of observations. Regarding the training
objective of the Critic sub-network, we simply opt for minimizing the squared
error when a collision was predicted and when it occurred:
$\operatorname*{arg\,min}_{\psi}\mathbb{E}_{z\sim
p_{data}(z)}\left\lVert{o(\hat{\mathbf{Y}},F_{s},F_{d};\psi)-R_{\mathbf{\hat{Y}}}}\right\rVert_{2},$
where $R_{\mathbf{\hat{Y}}}$ indicates the presence or not of a collision and
can be computed by any collision checking algorithm, here [27]. It is
important to note that following [6, 2, 1], all future predictions by the
generator and the Critic are performed sequantially, one step at a time. The
Critic takes as input the generated trajectory from $g^{dec}$ up to time $t$,
while the collision checker algorithm computes $R_{\mathbf{\hat{Y}}}$ given
the generated trajectory at the next time step, $t+1$.

The maximum reward is attained when the future trajectories cause minimum
number of collisions. In the end, the Critic sub-network acts as a regularizer
to the rest of the network, guiding the generator to generate novel
trajectories that cause mimimal collisions.

###  3.4 Attentive Scene Representation

To generate novel, accurate, diverse and safe future trajectories, a good
representation of the scene as well as the rest of the agents in the scene is
required. What is more, this representation must attend to the relevant
elements in the vicinity of every agent.

To this end, we propose the following tensor representations $F_{d},F_{s}$ for
capturing possible dynamic interactions between the agents themselves, as well
as the static interactions between the agents and the environment,
respectively, see Fig. 1. Both representations for dynamic and static
interactions follow the same design. We center a grid around the $i$-th agent
$G_{i}=\\{g_{rc}\\}$, where $r$ and $c$ run over the height and width of the
grid. We then check if a given grid location is occupied by either another
agent or a static element. If the $rc$ grid location is occupied by an agent
$j$, we then set the respective entry in $F_{d,i}$ to the hidden state of that
agent, namely $F_{d,i}^{rc}=h_{j}$. Similarly for static features, if the $rc$
grid location is occupied by a static element we set the respective entry in
$F_{s,i}$ to the one-hot vector representation of that static element.

Last, we employ an attention mechanism to prioritize certain elements in the
latent state representations $F_{d},F_{s}$ of the $i$-th agent. The attention
is applied spatially, such that specific coordinates on the grid receive
higher importance than others.

###  3.5 Architecture

For the encoder we opt for a Long short-term memory (LSTM) network [1]. The
observed trajectories are first embedded into a multilayer perceptron (MLP) in
order to obtain a fixed length vector,
$\mathbf{e}=\texttt{MLP}(\mathbf{X};W_{emb})$, where $\mathbf{X}$ are the
displacement coordinates of the pedestrians in the scene at all time steps
$t=[-T_{obs},...,0]$. The LSTM then updates its hidden states sequentially.
The last hidden state of the LSTM is concatenated with the static and dynamic
attentive representations $F_{d},F_{s}$ and the noise input $z$ and then
passed to the decoder. For the decoder we rely also on a separate LSTM network
that _does not_ share weights with the encoder. The LSTM decoder generates the
hidden states of each person for $t=[0,...,T_{\text{pred}}]$ predicted time
steps. The final hidden state is converted into a displacement vector and
added onto the previous prediction $\hat{\mathbf{Y}}^{t-1}$,
$\hat{\mathbf{Y}}^{t}=\hat{\mathbf{Y}}^{t-1}+\texttt{MLP}(\mathbf{h}_{dec}^{t};W_{T})$
where $W_{T}$ is a learnable transformation matrix. The critic takes $K$
predictions $\hat{\mathbf{Y}}_{i}^{t}$ of agent $i$, together with the state
in which that agent was $F_{s}^{t},F_{d}^{t}$ and assigns a reward for each
time-step. In this way, it learns to predict the likelihood of collision at
$t+1$, which is dependant on the motion of agent $i$, but also of another
agent $j$. If there are multiple approaching neighbors, chances of a collision
are higher and since we want to reward avoiding collisions especially in these
difficult situations, we do not normalize over the number of agents. We
implement the Critic as an LSTM and an MLP and train it with Mean Squared
Error between the number of collisions.

##  4 Experiments

We evaluate our method on two state-of-the-art datasets used for future
trajectory prediction of pedestrians.

Data UCY [28] has 3 videos, recorded from a bird’s eye point of view in 2
different outdoor locations. It contains annotated positions for pedestrians
only. The average length of a video is 361 seconds and it has 250
approximately people per video. Following[6] we use the leave-one-out approach
during training. The Stanford Drone Dataset [29] contains annotated positions
of various types of agents that navigate inside a university campus. On
average, there exist 104 agents per video. We use the standard train and test
splits provided by TrajNet [15].

Evaluation metrics We use the standard metrics reported in the literature [6,
1, 2], namely the minimum Average Displacement Error,
$\texttt{mADE}=\min_{K}\frac{1}{T}\sum_{t=1}^{T}\left\lVert{Y_{i}^{t}-\hat{Y}_{i}^{t}}\right\rVert_{2}$,
and minimum Final Displacement Error,
$\texttt{mFDE}=\min_{K}\left\lVert{Y_{i}^{T_{pred}}-\hat{Y}_{i}^{T_{pred}}}\right\rVert_{2}$
over best sample $K$.

Last, to evaluate the capacity of the methods to avoid collisions, we measure
the number of collisions with dynamic and static elements. The number of
collisions is computed as
$\texttt{NC}=\sum_{i,j=1}^{k}\sum_{t=1}^{T}\delta_{ij}\cdot\mathbf{1}(\hat{Y}_{i}^{t}-\hat{Y}_{j}^{t}<\epsilon)$
where $\delta_{ij}$ is the kronecker delta function, $\mathbf{1}$ is the
indicator function, and $\epsilon=0.10$ (here $i$ and $j$ are two different
agents).

Implementation details Following standard protocol [15] we observe
trajectories for $T_{\text{obs}}=8$ times steps (3.2 seconds). We then predict
the next $T_{\text{pred}}=12$ (4.8 seconds) time steps. The encoder and
decoder LSTM have an embedding size of 32 each. We use batch normalization
[30] and Adam optimization [31] with learning rates $10^{-}3$ for the
generator and $10^{-}4$ for the Critic and the discriminator. Like [6] we
bound and normalize distances between agents and the boundary points within a
neighborhood.

Baselines. We report comparisons with state-of-the-art methods on the UCY and
SDD: SocialGAN [1], SocialLSTM [6], Car-Net [32], SoPhie [2], DESIRE [3], as
well as a linear regressor that estimates parameters by minimizing the least
square error and Social Force implemented in [2].

###  4.1 Evaluating attentive scene representation features

First, we evaluate whether the proposed attentive scene representations allow
for better generation of future trajectories. In Table 1 we show results for
trajectory generation when we turn on and off the proposed static and dynamic
pooling, while also turning off the Critic module to minimize external
influences. Following dataset conventions, errors are measured in meters for
UCY and pixels for SDD. We observe that when switching on the proposed
attentive scene representations we obtain lower errors in terms of mADE and
mFDE for both UCY and SDD. We conclude that the attentive scene
representations allows for more accurate future trajectory generation.

Table 1:  Minimum Average (ADE) and Final Displacement Error (FDE) for
SafeCritic with and without Attentive Scene Representation (ASR). We conclude
that the proposed interaction poolings lead to better modelling of the
interactions for predicting future trajectorieslight increase of the speed.
Table 2: Number of collisions between agents on UCY. Our model scores 0.01 vs
0.06 and 0.05 compared to LSTM, SocialLSTM and SocialGAN on the UCYSDD.

Data Without ASR With ASR mADE / mFDE mADE / mFDE UCY 0.40 / 0.93 0.39 / 0.84
SDD 16.18 / 29.11 15.06 / 28.45

Data GT Linear SocialGAN SoPhie SafeCritic UCY Students3 0.12 1.24 0.55 0.62
0.48 Zara1 0.00 3.77 1.74 1.02 0.18 Zara2 0.73 3.63 2.02 1.46 1.22 Average
0.28 2.88 1.43 1.03 0.62

###  4.2 Evaluating the safety of SafeCritic

Next, we examine how safe the generated trajectories are in terms of the
number of collisions caused. We report results in Table 2. A collision is
defined with a threshold of 0.10m between two locations. The videos in SDD are
less crowded and there are no people within a proximity of 0.10 m. We make two
observations. Under all settings and data sets SafeCritic generates three
times fewer collisions compared to SocialGAN, Sophie and a linear regressor.
What is more, SafeCritic is quite close to the number of the collisions that
are actually observed in the ground truth trajectories. This is important
because collisions may actually happen in real life and thus a model should
also generate collisions when these are likely to occur.

###  4.3 Evaluating the accuracy of SafeCritic

Next, we examine how accurate the generated future trajectories are. We report
results in Table 3 for UCY and Table 4 for SDD. We observe that SafeCritic
generates more accurate future trajectories than all other baselines in UCY.

Table 3: Given 8 past time steps the proposed SafeCritic is more accurate in
predicting the future 12 timesteps on UCY (3 videos, 2492 frames).

Data SocialLSTM SocialGAN SoPhie Linear LSTM SafeCritic mADE / mFDE mADE /
mFDE mADE / mFDE mADE / mFDE mADE / mFDE mADE / mFDE UCY Students3 0.67 / 1.40
0.76 / 1.52 0.54 / 1.24 0.82 / 1.59 0.61 / 1.31 0.52 / 1.23 Zara1 0.47 / 1.00
0.35 / 0.68 0.30 / 0.63 0.62 / 1.21 0.41 / 0.88 0.29 / 0.60 Zara2 0.56 / 1.17
0.42 / 1.01 0.38 / 0.78 0.77 / 1.48 0.52 / 1.11 0.30 / 0.69 Average 0.56 /
1.19 0.51 / 0.96 0.41 / 0.88 0.74 / 1.43 0.51 / 1.10 0.37 / 0.84

Importantly, SafeCritic comfortably outperforms all methods on the much larger
and complex SDD dataset. We also evaluate the diversity of the generated
trajectories in Fig. 2. We conclude that SafeCritic predicts future
trajectories that are not only accurate, but also cause much fewer collisions
and that are more diverse.

Table 4: Given 8 past time steps the proposed SafeCritic is more accurate in
predicting the future 12 timesteps on SDD (48 videos, 14558 frames).

Data SocialLSTM SocialGAN SoPhie Desire Linear SF CAR-Net SafeCritic mADE /
mFDE mADE / mFDE mADE / mFDE mADE / mFDE mADE / mFDE mADE / mFDE mADE / mFDE
mADE / mFDE SDD 31.19 / 56.97 27.25 / 41.44 16.27 / 29.38 19.25 / 34.05 37.11
/ 63.51 36.48 / 58.14 25.72 / 51.80 14.39 / 27.39

###  4.4 Qualitative results

Last, we present some qualitative results for SafeCritic. In Fig. 2 we plot
various generated trajectories at a specific moment in time, examining whether
mode collapse is happening. We observe that SafeCritic does not suffer from
mode collapse [24], as also quantified earlier.

Figure 2: As expected SafeCritic generates diverse trajectories. Five
trajectories of two persons walking towards a goal (star) in the observed
(left) and predicted (right) trajectories.

Further, in Fig. 3 we show that the generated future trajectories are nicely
sampled such that the attended agents and static elements are avoided and no
collisions are generated.

Figure 3: The attention (left) of a person (red) show the importance of other
agents (blue and green) (left columns) and of nearby trees and buildings
(right columns), resulting in alternating his trajectory (right).

Last, in Fig. 4 we show how the predicted locations of SafeCritic do not
collide. The second row shows how multiple sampled trajectories of SafeCritic
move around a lawn.

Observed

SocialGAN

SafeCritic

Observed

SocialGAN

SafeCritic

Figure 4: As expected SafeCritic generates safe trajectories. Trajectories of
two walking persons with final ground truth location (star). Collisions
between people for Coupa 0 (left columns) with the locations at t=1s and 2s
(circles) and going onto lawn for Bookstore 2 (right columns).

##  5 Conclusion

In this work we introduced SafeCritic, a model that is based on training in an
adversarial way along with a reinforcement learning based reward network to
predict both accurate and safe trajectories. A novel way of handling static
and dynamic features in both the generator and reward network makes it
possible to discard unsafe trajectories. Our proposed collision metric enables
the evaluation of safety aspect of predicted trajectories. We show that
SafeCritic outperforms baseline models both on accuracy as well as on the
number of collisions. Furthermore, our experiments show that a convex
combination of trajectory’s ‘safe’ and ‘real’ aspects achieves a perfect
balance for the final prediction.

## References

  * [1] Agrim Gupta, Justin Johnson, Li Fei-Fei, Silvio Savarese, and Alexandre Alahi.  Social gan: Socially acceptable trajectories with generative adversarial networks.  In IEEE Conference on Computer Vision and Pattern Recognition (CVPR), number CONF, 2018.
  * [2] Amir Sadeghian, Vineet Kosaraju, Ali Sadeghian, Noriaki Hirose, and Silvio Savarese.  Sophie: An attentive gan for predicting paths compliant to social and physical constraints.  arXiv preprint arXiv:1806.01482, 2018.
  * [3] Namhoon Lee, Wongun Choi, Paul Vernaza, Christopher B Choy, Philip HS Torr, and Manmohan Chandraker.  Desire: Distant future prediction in dynamic scenes with interacting agents.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 336–345, 2017.
  * [4] Ian Goodfellow, Jean Pouget-Abadie, Mehdi Mirza, Bing Xu, David Warde-Farley, Sherjil Ozair, Aaron Courville, and Yoshua Bengio.  Generative adversarial nets.  In Advances in neural information processing systems, pages 2672–2680, 2014.
  * [5] Andrew Y. Ng and Stuart J. Russell.  Algorithms for inverse reinforcement learning.  In ICML, 2000.
  * [6] Alexandre Alahi, Kratarth Goel, Vignesh Ramanathan, Alexandre Robicquet, Li Fei-Fei, and Silvio Savarese.  Social lstm: Human trajectory prediction in crowded spaces.  In Proceedings of the IEEE International Conference on Computer Vision and Pattern Recognition, 2016.
  * [7] Jon Gauthier.  Conditional generative adversarial nets for convolutional face generation.  Class Project for Stanford CS231N: Convolutional Neural Networks for Visual Recognition, Winter semester, 2014(5):2, 2014.
  * [8] Mehdi Mirza and Simon Osindero.  Conditional generative adversarial nets.  arXiv preprint arXiv:1411.1784, 2014.
  * [9] Nicholas Rhinehart, Kris M Kitani, and Paul Vernaza.  R2p2: A reparameterized pushforward policy for diverse, precise generative path forecasting.  In Proceedings of the European Conference on Computer Vision (ECCV), pages 772–788, 2018.
  * [10] Anirudh Vemula, Katharina Muelling, and Jean Oh.  Social attention: Modeling attention in human crowds.  In 2018 IEEE International Conference on Robotics and Automation (ICRA), pages 1–7. IEEE, 2018.
  * [11] Dirk Helbing and Peter Molnar.  Social force model for pedestrian dynamics.  Physical review E, 51(5):4282, 1995.
  * [12] Avneesh Sud, Russell Gayle, Erik Andersen, Stephen Guy, Ming Lin, and Dinesh Manocha.  Real-time navigation of independent agents using adaptive roadmaps.  In ACM SIGGRAPH 2008 classes, page 56. ACM, 2008.
  * [13] Peter Trautman and Andreas Krause.  Unfreezing the robot: Navigation in dense, interacting crowds.  In Intelligent Robots and Systems (IROS), 2010 IEEE/RSJ International Conference on, pages 797–803. IEEE, 2010.
  * [14] Stefan Becker, Ronny Hug, Wolfgang Hübner, and Michael Arens.  An evaluation of trajectory prediction approaches and notes on the trajnet benchmark.  arXiv preprint arXiv:1805.07663, 2018.
  * [15] Amir Sadeghian, Vineet Kosaraju, Agrim Gupta, Silvio Savarese, and Alexandre Alahi.  Trajnet: Towards a benchmark for human trajectory prediction.  arXiv preprint, 2018.
  * [16] Chelsea Finn, Paul F. Christiano, Pieter Abbeel, and Sergey Levine.  A connection between generative adversarial networks, inverse reinforcement learning, and energy-based models.  CoRR, abs/1611.03852, 2016.
  * [17] Nicola De Cao and Thomas Kipf.  Molgan: An implicit generative model for small molecular graphs.  arXiv preprint arXiv:1805.11973, 2018.
  * [18] Gabriel Lima Guimaraes, Benjamin Sanchez-Lengeling, Carlos Outeiral, Pedro Luis Cunha Farias, and Alán Aspuru-Guzik.  Objective-reinforced generative adversarial networks (organ) for sequence generation models.  arXiv preprint arXiv:1705.10843, 2017.
  * [19] Vijay R Konda and John N Tsitsiklis.  Actor-critic algorithms.  In Advances in neural information processing systems, pages 1008–1014, 2000.
  * [20] Richard S Sutton, David A McAllester, Satinder P Singh, and Yishay Mansour.  Policy gradient methods for reinforcement learning with function approximation.  In Advances in neural information processing systems, pages 1057–1063, 2000.
  * [21] Richard Stuart Sutton.  Temporal credit assignment in reinforcement learning.  1984\.
  * [22] AG Barto.  Neuron-like adaptive elements that can solve difficult learning control-problems.  Behavioural Processes, 9(1), 1984.
  * [23] Richard S Sutton and Andrew G Barto.  Introduction to reinforcement learning, volume 135.  MIT press Cambridge, 1998.
  * [24] Tim Salimans, Ian Goodfellow, Wojciech Zaremba, Vicki Cheung, Alec Radford, and Xi Chen.  Improved techniques for training gans.  In Advances in Neural Information Processing Systems, pages 2234–2242, 2016.
  * [25] Timothy P Lillicrap, Jonathan J Hunt, Alexander Pritzel, Nicolas Heess, Tom Erez, Yuval Tassa, David Silver, and Daan Wierstra.  Continuous control with deep reinforcement learning.  arXiv preprint arXiv:1509.02971, 2015.
  * [26] Nicolas Heess, Jonathan J Hunt, Timothy P Lillicrap, and David Silver.  Memory-based control with recurrent neural networks.  arXiv preprint arXiv:1512.04455, 2015.
  * [27] Christer Ericson.  Real-time collision detection.  CRC Press, 2004.
  * [28] Laura Leal-Taixé, Michele Fenzi, Alina Kuznetsova, Bodo Rosenhahn, and Silvio Savarese.  Learning an image-based motion context for multiple people tracking.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 3542–3549, 2014.
  * [29] Alexandre Robicquet, Amir Sadeghian, Alexandre Alahi, and Silvio Savarese.  Learning social etiquette: Human trajectory understanding in crowded scenes.  In European conference on computer vision, pages 549–565. Springer, 2016.
  * [30] Sergey Ioffe and Christian Szegedy.  Batch normalization: Accelerating deep network training by reducing internal covariate shift.  arXiv preprint arXiv:1502.03167, 2015.
  * [31] Diederik P Kingma and Jimmy Ba.  Adam: A method for stochastic optimization.  arXiv preprint arXiv:1412.6980, 2014.
  * [32] Amir Sadeghian, Ferdinand Legros, Maxime Voisin, Ricky Vesel, Alexandre Alahi, and Silvio Savarese.  Car-net: Clairvoyant attentive recurrent network.  arXiv preprint arXiv:1711.10061, 2017.


versus age (Figure 14 of Agertz et al., 2021) is nearly identical to ours for
m12i, both in value and shape: both of our simulations show near exponential
increases in $\sigma_{\rm Z,now}$ with age for stars younger than 6 Gyr and a
jump at 6 Gyr. Our results also agree well with the NIHAO-UHD suite (Buck et
al., 2020): their 5 galaxies, which span $M_{\text{star}}=1.5-16\times
10^{10}\,\text{M}_{\odot}$, have $\sigma_{\rm Z,now}\approx
20-35~{}\text{km}\,\text{s}^{-1}$ for young stars, and $\approx
60-120~{}\text{km}\,\text{s}^{-1}$ for stars 12 Gyr old, similar to ours. Any
difference in the exact age dependence of each galaxy’s dispersion follows
from its unique formation history.

Grand et al. (2016) analyzed 16 MW-mass galaxies from the Auriga simulations,
also finding broadly similar disk-wide $\sigma_{\rm Z,now}$ and $\sigma_{\rm
Z,form}$ as we do. Similarly to FIRE, NIHAO-UHD, and VINTERGARTEN, Auriga
simulations have typical disk-wide $\sigma_{\rm Z,now}\approx
20-25~{}\text{km}\,\text{s}^{-1}$. They also concluded that cosmological disk
settling was generally the primary effect that set the relation between the
velocity dispersion and the age of stars today. However, Grand et al. (2016)
found that bars were the strongest contributor to vertical dynamical heating.
This seemingly contradicts our results, because stars in all of our galaxies
had significant post-formation increases to their disk-wide $\sigma_{\rm Z}$,
including stars in unbarred or weakly-barred galaxies. Furthermore, even our
galaxies that do house significant bars near $z=0$ (such as m12m) did not at
early times, and our results indicate that heating is most significant for
stars that formed in the Pre-Disk Era. Although we do not study the effect of
bars in this work, Ansar et al., (in prep.) quantifies the incidence of bars
in the FIRE-2 simulations, and in future work we will address the drivers of
dynamical heating in our simulations, including bars. While it is possible
that the bar heated these old stars more recently, bars most likely do not
heat only old stars.

Our results may differ from those of Grand et al. (2016) because of numerical
differences. First, Auriga uses the subgrid model from Springel & Hernquist
(2003), treating the ISM as a two-phase medium with an effective equation of
state that inhibits the formation of dense gas, so that star formation occurs
at much lower density, $n>0.13$ cm-3, compared with $n>1000$ cm-3 in FIRE-2.
House et al. (2011) found that low thresholds in density impose a dispersion
‘floor’ in $\sigma_{\rm form}$, even in high-resolution simulations (50 pc),
because stars form from higher-temperature gas. House et al. (2011) find that
this floor effectively goes away if one increases the star-formation threshold
from $n=0.1$ to 1000 cm-3 (see also Martig et al., 2014; Kumamoto et al.,
2017). That said, given that Grand et al. (2016) find that $\sigma_{\rm form}$
decreases over time, this dispersion ‘floor’ may have little effect on their
results.

Second, certain mechanisms that induce dynamical heating in our simulations
may be different (and less significant) in the Auriga simulations. The
simulations in Grand et al. (2016) had baryonic mass resolution of $4\times
10^{4}\,\text{M}_{\odot}$, $8-10$ times larger than ours, and they used lower
spatial resolution, with gravitational force softenings of 375 pc physical for
stars and gas (which is comparable to the thin-disk scale height of the MW),
while our FIRE-2 simulations use $2.7-4$ pc for star particles and $\approx 1$
pc minimum and $\sim 4$ pc ISM-averaged for gas cells. Thus, Auriga
simulations do not resolve the dense ISM, including important dynamical
structures like GMCs that likely heat stars particles after their formation in
FIRE. Furthermore, the modeling of stellar feedback in Auriga is markedly
different: supernovae launch ‘wind particles’ that are temporarily decoupled
from the hydrodynamics until the particle encounters sufficiently low-density
gas, which does not couple that feedback directly to dense gas in a star-
forming region. In turn, our simulations may exhibit more gas turbulence. That
said, the Auriga simulations do model certain physics that these FIRE-2
simulations do not include, most importantly, AGN feedback.

Grand et al. (2016) also found that spiral arms did not considerably heat the
disk, while we proposed that spiral arms are likely responsible for the
significant increase in $\sigma_{\rm tot,now}$ of young stars, in part because
$\sigma_{\rm R}$ dominates the overall dispersion (see also Orr et al., 2022).
However, Grand et al. (2016) only considered vertical heating. Because the
Auriga simulations do not resolve GMCs, spiral arms likely have negligible
effects on the vertical dispersion, but may have more significant impacts on
the in-plane dispersions. For example, while the FIRE-2 simulations resolve
GMCs, young stars in FIRE-2 have vertical effective heating rates that are
$\approx 4\times$ lower than their radial rates, as Figure 6 showed. Many
previous works – both observational and theoretical – only examined the
vertical velocity dispersion. However, vertical post-formation dynamical
heating operates differently than total and in-plane heating. For example,
Figure 3.4 (top) shows that the amount of vertical dispersion added post-
formation increases approximately linearly across all ages, while in-plane
heating plateaus for Early-Disk stars. Similarly, Figure 3.4 (bottom) shows
that the present-day vertical dispersion has the largest relative impact from
post-formation heating: the dispersion added from heating typically begins to
equal (or surpass) the dispersion from formation at $\approx$ 4 Gyr for the
vertical component, but does not rival the in-plane dispersions from formation
until $\approx$ 11 Gyr. Furthermore, the vertical dispersion typically
exhibits the smallest values. For these reasons, we do not recommend using
vertical heating to draw inferences regarding total heating.

As we discussed above, nearly all current cosmological zoom-in simulations
form MW-mass galaxies with typical disk-wide $\sigma_{\rm Z,now}\gtrsim
20~{}\text{km}\,\text{s}^{-1}$, not as dynamically cold as the (solar-
neighborhood) measurements of MW today. However, Bird et al. (2021) analyzed a
single MW-mass galaxy, h277, selected from the g14 suite of simulations
(Christensen et al., 2012), generated using the parallel $N$–body+SPH code
GASOLINE (Wadsley et al., 2004). That simulation had disk-wide $\sigma_{\rm
Z,now}<10~{}\text{km}\,\text{s}^{-1}$, measured within an annulus of
$R=8-10~{}\text{kpc}$, for stars younger than 1 Gyr, similar to (primarily
locally-measured) values for the MW (for example Nordström et al., 2004;
Mackereth et al., 2019b).

The cold disk-wide kinematics of h277 likely reflect its formation history:
its disk began to settle early, $\approx 9~{}\text{Gyr}$ ago, and it had no
major mergers since $z\approx 2$. That said, one of our galaxies, Romeo, has a
similar time of disk onset and merger history as h277 but has a larger (disk-
wide) $\sigma_{\rm Z}\approx 16~{}\text{km}\,\text{s}^{-1}$ for stars 1 Gyr
old today. However, we note that Romeo has a local $\sigma_{\rm Z}\approx
12~{}\text{km}\,\text{s}^{-1}$, which is closer, albeit still higher, than
measurements of the MW. Therefore, differences in our simulations’ physical
and numerical modeling likely contribute to our unequal dispersions. Although
the baryonic mass resolution, modeling of dense ISM, and molecular-based star-
formation criteria are broadly similar between our FIRE-2 simulations and
h277, key numerical differences exist between our simulations: h277 used SPH
for the hydrodynamics solver (instead of MFM); used larger gravitational force
softening of 173 pc for all particle species; used the ‘blastwave’ feedback
model for core-collapse supernovae in which cooling is turned off for gas
within the blast wave radius to mimic the adiabatic expansion of a remnant
(Stinson et al., 2006); used only thermal energy injection for white-dwarf
(Type Ia) supernovae; and did not include stellar winds or radiative feedback.
While it is unclear how any of these individual aspects do or do not affect
their results, differences in feedback modeling may contribute to their lower
dispersions to some degree. FIRE-2 accounts for energy and momentum injection
for both core-collapse and white-dwarf supernovae, and includes stellar winds,
radiation pressure, and photo-ionization. This feedback launches starburst-
induced, star-forming outflows at early times (see the top-right panel of our
Figure 1 and El-Badry et al. (2018a)), which Bird et al. (2021) notes are
absent in their simulation. As such, it is possible that our additional
feedback channels contribute to larger gas velocity dispersions and
subsequently larger $\sigma_{\rm form}$ values in our simulations, though it
is unclear if this persists to later cosmic times. Similarly, Hopkins et al.
(2012) and Agertz et al. (2013) both show that the structure of the ISM can
non-trivially vary depending on the feedback model used. Despite any
differences, we again emphasize that both the disk-wide $\sigma_{\rm Z,form}$
and $\sigma_{\rm Z,now}$ exhibited similar trends with age in h277 and our
FIRE-2 simulations.

####  4.2.4 Connections to the early formation of the MW

As we discussed in Sections 3.6 and 3.8, the cold kinematics of MW stars may
indicate that the MW’s disk began to settle unusually early in its history.
For example, using our predicted result that the disk-wide kinematics of young
stars today correlates with when the disk began to settle, and placing the
MW’s value of $v_{\rm\phi}/\sigma_{\rm tot}\approx 8$ along the relation from
our simulations in Figure 10, this implies that the MW’s time of disk onset
was $\gtrsim 11~{}\text{Gyr}$ ago, which agrees with recent estimates from
observations (for example Belokurov & Kravtsov, 2022; Conroy et al., 2022;
Xiang & Rix, 2022). Thus, our analysis of disk settling times in cosmological
simulations and its relation to present-day kinematics complements and
augments these works.

That said, measurements of the MW’s kinematics are largely restricted to one
local patch: the solar neighborhood. As we show in Figure 8, locally-measured
velocity dispersions in our simulations (at formation) are typically $\approx
2$ times smaller than disk-wide measurements. In turn, while all of our disk-
wide dispersions are larger than the MW’s, 3 of our galaxies have present-day
local dispersions lower than the MW’s. Furthermore, the galaxy-to-galaxy
scatter in the relation between disk-wide and local kinematics is significant.
For example, stars today younger than 100 Myr in one of our most MW-like
galaxies, Romeo, have disk-wide $\sigma$ that is $1.1-1.2$ times higher than
its local value, while m12i, which is otherwise similar to Romeo, has disk-
wide values that are $1.5-1.7$ times higher than local values. More work is
needed to determine how the MW’s disk-wide dispersion compares to its local
dispersion, and if these differences affect its placement on our disk onset -
young stellar kinematics relation. In McCluskey et al. (in prep.) we will more
rigorously compare our simulations to the MW as well as other observed
galaxies.

We now discuss connections with some of these recent observational works.
Conroy et al. (2022) combined Gaia astrometry and H3 Survey spectroscopy to
posit that the MW disk formed in three dynamical eras, broadly consistent with
our Pre-Disk, Early-Disk, and Late-Disk Eras. Conroy et al. (2022) proposed
that the MW began in a ‘simmering phase’, during which the star formation
efficiency was low and stars formed kinematically hot, but that $\approx
12-13~{}\text{Gyr}$ ago the MW transitioned into a ‘boiling phase’, during
which the SFE strongly increased and stars formed with increasingly ‘disky’
dynamics, marking the ‘birth of the Galactic disk’. This thick disk continued
to grow and settle in a ‘boiling phase’ for $\approx 3-4~{}\text{Gyr}$ (from
$z\approx 4$ to $z\approx 1$) until the Gaia-Sausage-Enceladus (GSE) merger,
after which the star formation efficiency decreased and stars formed in a
dynamically-colder, thin disk until today.

Similarly, Belokurov & Kravtsov (2022) combined Gaia astrometry with
spectroscopy from the APOGEE Data Release 17. Although Belokurov & Kravtsov
(2022) lacked empirical age estimates, their general picture agrees with ours:
old, metal-poor, in-situ stars were born in a ‘messy’ phase, then the MW
‘spun-up’ as stars rapidly became metal-rich and rotational earlier than
$\approx 8$ Gyr ago. Belokurov & Kravtsov (2022) posited that this rapid
formation of the MW’s disk took place over 1 - 2 Gyr, with the GSE merger
occurring after the initial formation of a thick-disk, such that the merger
heated disk stars onto halo-like orbits.

On the other hand, Xiang & Rix (2022) conclude that the formation of the in-
situ halo and thick-disk overlapped. In this picture, the GSE merger then
occurred $\approx 11~{}\text{Gyr}$ ago, $1-2~{}\text{Gyr}$ earlier than
previous estimates, and did not strongly heat pre-existing disk stars, or
cause a transition to a thin-disk, but instead enhanced thick-disk formation.

Our results help connect these present-day observations to understanding the
MW’s formation history. In particular, the presence of coherent rotation in a
stellar population today does not necessarily indicate the presence of
coherent rotation at the time of formation. Observations show that older, more
metal-poor stars in the MW rotate slower than younger stars, but they still
display some coherent rotation. In our simulations, this is not primarily the
result of these Pre-Disk stars being born with coherent rotation, but rather,
this tends to arise from post-formation torquing (see Figure 1), for example,
caused by mergers (see Bellardini et al. in preparation). m12m provides an
illustrative example: it experienced several major mergers at early times,
broadly similar to estimates of when the GSE merger occurred. As Figure 3
shows, m12m has old stars with $v_{\rm\phi,now}\approx
200~{}\text{km}\,\text{s}^{-1}$, despite having $v_{\rm\phi,form}\approx
0-50~{}\text{km}\,\text{s}^{-1}$. Thus, mergers like the GSE merger could have
torqued old stars in the MW onto more coherently-rotating orbits.

Ultimately, the MW’s kinematics reflect its individual formation history.
While the MW provides a useful benchmark, the fact that most current
cosmological simulations form galaxies with hotter kinematics than the MW may
speak to the MW’s early settling time rather than to limitations of the
simulations. Indeed, the comparison of the MW to M31 and M33 in Figure 7
affirms that the MW may be an outlier among galaxies at its mass. Therefore,
future works that study the kinematics of stars and their relation to stellar
ages in external galaxies will provide excellent opportunities to place the MW
in a statistical context, to study the combination of cosmological disk
settling and dynamical heating across a range of formation histories (Dorman
et al., 2015; Leaman et al., 2017; Quirk et al., 2022).

## Acknowledgements

We thank Jonathan Stern for helpful comments and the participants of the Disk
Formation Workshop at UC Irvine for their enlightening discussions. FM and AW
received support from: NSF via CAREER award AST-2045928 and grant AST-2107772;
NASA ATP grant 80NSSC20K0513; HST grants AR-15809, GO-15902, GO-16273 from
STScI; and a Scialog Award from the Heising-Simons Foundation. CAFG was
supported by NSF through grants AST-2108230 and CAREER award AST-1652522; by
NASA through grants 17-ATP17-0067 and 21-ATP21-0036; by STScI through grant
HST-GO-16730.016-A; by CXO through grant TM2-23005X; and by the Research
Corporation for Science Advancement through a Cottrell Scholar Award. PFH
recieved supported by NSF Research Grants 1911233 and 20009234, NSF CAREER
grant 1455342, and NASA grants 80NSSC18K0562 and HST-AR-15800.001-A. We
performed some of this work at the Aspen Center for Physics, supported by NSF
grant PHY-1607611. We ran simulations using: XSEDE, supported by NSF grant
ACI-1548562; Blue Waters, supported by the NSF; Frontera allocations AST21010
and AST20016, supported by the NSF and TACC; Pleiades, via the NASA HEC
program through the NAS Division at Ames Research Center.

## Data Availability

All of the Python code that we used to generate these figures are available at
https://fmccluskey.github.io, which uses the publicly available package
https://bitbucket.org/awetzel/gizmo_analysis (Wetzel & Garrison-Kimmel, 2020).
The FIRE-2 simulations are publicly available (Wetzel et al., 2023) at
http://flathub.flatironinstitute.org/fire. Additional FIRE simulation data is
available at https://fire.northwestern.edu/data. A public version of the Gizmo
code is available at http://www.tapir.caltech.edu/~phopkins/Site/GIZMO.html.

## References

  * Agertz et al. (2013) Agertz O., Kravtsov A. V., Leitner S. N., Gnedin N. Y., 2013, ApJ, 770, 25
  * Agertz et al. (2021) Agertz O., et al., 2021, MNRAS, 503, 5826
  * Anders et al. (2023) Anders F., et al., 2023, arXiv e-prints, p. arXiv:2304.08276
  * Anglés-Alcázar et al. (2017) Anglés-Alcázar D., Faucher-Giguère C.-A., Kereš D., Hopkins P. F., Quataert E., Murray N., 2017, MNRAS, 470, 4698
  * Aumer et al. (2016) Aumer M., Binney J., Schönrich R., 2016, MNRAS, 462, 1697
  * Barbanis & Woltjer (1967) Barbanis B., Woltjer L., 1967, ApJ, 150, 461
  * Bellardini et al. (2022) Bellardini M. A., Wetzel A., Loebman S. R., Bailin J., 2022, MNRAS, 514, 4270
  * Belokurov & Kravtsov (2022) Belokurov V., Kravtsov A., 2022, MNRAS, 514, 689
  * Belokurov et al. (2018) Belokurov V., Erkal D., Evans N. W., Koposov S. E., Deason A. J., 2018, MNRAS, 478, 611
  * Benincasa et al. (2020) Benincasa S. M., et al., 2020, MNRAS, 497, 3993
  * Binney & Tremaine (2008) Binney J., Tremaine S., 2008, Galactic Dynamics: Second Edition
  * Bird et al. (2012) Bird J. C., Kazantzidis S., Weinberg D. H., 2012, MNRAS, 420, 913
  * Bird et al. (2021) Bird J. C., Loebman S. R., Weinberg D. H., Brooks A. M., Quinn T. R., Christensen C. R., 2021, MNRAS, 503, 1815
  * Birkin et al. (2023) Birkin J. E., et al., 2023, arXiv e-prints, p. arXiv:2301.05720
  * Bonaca et al. (2017) Bonaca A., Conroy C., Wetzel A., Hopkins P. F., Kereš D., 2017, ApJ, 845, 101
  * Bonaca et al. (2020) Bonaca A., et al., 2020, ApJ, 897, L18
  * Boulares & Cox (1990) Boulares A., Cox D. P., 1990, ApJ, 365, 544
  * Brook et al. (2004) Brook C. B., Kawata D., Gibson B. K., Freeman K. C., 2004, ApJ, 612, 894
  * Brook et al. (2012) Brook C. B., et al., 2012, MNRAS, 426, 690
  * Buck et al. (2020) Buck T., Obreja A., Macciò A. V., Minchev I., Dutton A. A., Ostriker J. P., 2020, MNRAS, 491, 3461
  * Buder et al. (2018) Buder S., et al., 2018, MNRAS, 478, 4513
  * Carlberg & Innanen (1987) Carlberg R. G., Innanen K. A., 1987, AJ, 94, 666
  * Carlberg & Sellwood (1985) Carlberg R. G., Sellwood J. A., 1985, ApJ, 292, 79
  * Casagrande et al. (2011) Casagrande L., Schönrich R., Asplund M., Cassisi S., Ramírez I., Meléndez J., Bensby T., Feltzing S., 2011, A&A, 530, A138
  * Ceverino et al. (2017) Ceverino D., Primack J., Dekel A., Kassin S. A., 2017, MNRAS, 467, 2664
  * Chan et al. (2019) Chan T. K., Kereš D., Hopkins P. F., Quataert E., Su K. Y., Hayward C. C., Faucher-Giguère C. A., 2019, MNRAS, 488, 3716
  * Chan et al. (2021) Chan T. K., Keres D., Gurvich A. B., Hopkins P., Trapp C., Ji S., Faucher-Giguere C.-A., 2021, arXiv e-prints, p. arXiv:2110.06231
  * Christensen et al. (2012) Christensen C., Quinn T., Governato F., Stilp A., Shen S., Wadsley J., 2012, MNRAS, 425, 3058
  * Conroy et al. (2022) Conroy C., et al., 2022, arXiv e-prints, p. arXiv:2204.02989
  * Conselice et al. (2005) Conselice C. J., Blackburne J. A., Papovich C., 2005, ApJ, 620, 564
  * Cui et al. (2012) Cui X.-Q., et al., 2012, Research in Astronomy and Astrophysics, 12, 1197
  * Das et al. (2020) Das P., Hawkins K., Jofré P., 2020, MNRAS, 493, 5195
  * Dashyan et al. (2018) Dashyan G., Silk J., Mamon G. A., Dubois Y., Hartwig T., 2018, MNRAS, 473, 5698
  * De Silva et al. (2015) De Silva G. M., et al., 2015, MNRAS, 449, 2604
  * Debattista et al. (2019) Debattista V. P., Gonzalez O. A., Sanderson R. E., El-Badry K., Garrison-Kimmel S., Wetzel A., Faucher-Giguère C.-A., Hopkins P. F., 2019, MNRAS, 485, 5073
  * Dehnen & Binney (1998) Dehnen W., Binney J. J., 1998, MNRAS, 298, 387
  * Dekel et al. (2009) Dekel A., et al., 2009, Nature, 457, 451
  * Dorman et al. (2015) Dorman C. E., et al., 2015, ApJ, 803, 24
  * Dubois et al. (2013) Dubois Y., Gavazzi R., Peirani S., Silk J., 2013, MNRAS, 433, 3297
  * Edvardsson et al. (1993) Edvardsson B., Andersen J., Gustafsson B., Lambert D. L., Nissen P. E., Tomkin J., 1993, A&A, 275, 101
  * El-Badry et al. (2016) El-Badry K., Wetzel A., Geha M., Hopkins P. F., Kereš D., Chan T. K., Faucher-Giguère C.-A., 2016, ApJ, 820, 131
  * El-Badry et al. (2018a) El-Badry K., et al., 2018a, MNRAS, 473, 1930
  * El-Badry et al. (2018b) El-Badry K., et al., 2018b, MNRAS, 480, 652
  * Elmegreen et al. (2007) Elmegreen D. M., Elmegreen B. G., Ravindranath S., Coe D. A., 2007, ApJ, 658, 763
  * Escala et al. (2018) Escala I., et al., 2018, MNRAS, 474, 2194
  * Faucher-Giguère (2018) Faucher-Giguère C.-A., 2018, MNRAS, 473, 3717
  * Faucher-Giguère et al. (2009) Faucher-Giguère C.-A., Lidz A., Zaldarriaga M., Hernquist L., 2009, ApJ, 703, 1416
  * Faucher-Giguère et al. (2013) Faucher-Giguère C.-A., Quataert E., Hopkins P. F., 2013, MNRAS, 433, 1970
  * Feldmann et al. (2017) Feldmann R., Quataert E., Hopkins P. F., Faucher-Giguère C.-A., Kereš D., 2017, MNRAS, 470, 1050
  * Flores Velázquez et al. (2021) Flores Velázquez J. A., et al., 2021, MNRAS, 501, 4812
  * Förster Schreiber & Wuyts (2020) Förster Schreiber N. M., Wuyts S., 2020, ARA&A, 58, 661
  * Förster Schreiber et al. (2009) Förster Schreiber N. M., et al., 2009, ApJ, 706, 1364
  * Gaia Collaboration et al. (2016) Gaia Collaboration et al., 2016, A&A, 595, A1
  * Gallart et al. (2019) Gallart C., Bernard E. J., Brook C. B., Ruiz-Lara T., Cassisi S., Hill V., Monelli M., 2019, Nature Astronomy, 3, 932
  * Garrison-Kimmel et al. (2018) Garrison-Kimmel S., et al., 2018, MNRAS, 481, 4133
  * Garrison-Kimmel et al. (2019) Garrison-Kimmel S., et al., 2019, MNRAS, 489, 4574
  * Genzel et al. (2015) Genzel R., et al., 2015, ApJ, 800, 20
  * Giocoli et al. (2007) Giocoli C., Moreno J., Sheth R. K., Tormen G., 2007, MNRAS, 376, 977
  * Girard et al. (2021) Girard M., et al., 2021, ApJ, 909, 12
  * Grand et al. (2016) Grand R. J. J., Springel V., Gómez F. A., Marinacci F., Pakmor R., Campbell D. J. R., Jenkins A., 2016, MNRAS, 459, 199
  * Grand et al. (2017) Grand R. J. J., et al., 2017, MNRAS, 467, 179
  * Gurvich et al. (2020) Gurvich A. B., et al., 2020, MNRAS, 498, 3664
  * Gurvich et al. (2023) Gurvich A. B., et al., 2023, MNRAS, 519, 2598
  * Guszejnov et al. (2020) Guszejnov D., Grudić M. Y., Offner S. S. R., Boylan-Kolchin M., Faucher-Giguère C.-A., Wetzel A., Benincasa S. M., Loebman S., 2020, MNRAS, 492, 488
  * Hafen et al. (2022) Hafen Z., et al., 2022, MNRAS, 514, 5056
  * Hahn & Abel (2011) Hahn O., Abel T., 2011, MNRAS, 415, 2101
  * Haywood et al. (2018) Haywood M., Di Matteo P., Lehnert M. D., Snaith O., Khoperskov S., Gómez A., 2018, ApJ, 863, 113
  * He et al. (2023) He H., Bottrell C., Wilson C., Moreno J., Burkhart B., Hayward C. C., Hernquist L., Twum A., 2023, arXiv e-prints, p. arXiv:2301.13250
  * Helmi et al. (2018) Helmi A., Babusiaux C., Koppelman H. H., Massari D., Veljanoski J., Brown A. G. A., 2018, Nature, 563, 85
  * Holmberg et al. (2009) Holmberg J., Nordström B., Andersen J., 2009, A&A, 501, 941
  * Hopkins (2015) Hopkins P. F., 2015, MNRAS, 450, 53
  * Hopkins (2016) Hopkins P. F., 2016, MNRAS, 455, 89
  * Hopkins et al. (2009) Hopkins P. F., Cox T. J., Younger J. D., Hernquist L., 2009, ApJ, 691, 1168
  * Hopkins et al. (2012) Hopkins P. F., Quataert E., Murray N., 2012, MNRAS, 421, 3488
  * Hopkins et al. (2018) Hopkins P. F., et al., 2018, MNRAS, 480, 800
  * Hopkins et al. (2020) Hopkins P. F., et al., 2020, MNRAS, 492, 3465
  * Hopkins et al. (2023) Hopkins P. F., et al., 2023, arXiv e-prints, p. arXiv:2301.08263
  * House et al. (2011) House E. L., et al., 2011, MNRAS, 415, 2652
  * Hung et al. (2019) Hung C.-L., et al., 2019, MNRAS, 482, 5125
  * Irodotou et al. (2022) Irodotou D., et al., 2022, MNRAS, 513, 3768
  * Jenkins & Binney (1990) Jenkins A., Binney J., 1990, MNRAS, 245, 305
  * Johnson et al. (2018) Johnson H. L., et al., 2018, MNRAS, 474, 5076
  * Kassin et al. (2012) Kassin S. A., et al., 2012, ApJ, 758, 106
  * Kim et al. (2023) Kim J.-G., Gong M., Kim C.-G., Ostriker E. C., 2023, ApJS, 264, 10
  * Koudmani et al. (2021) Koudmani S., Henden N. A., Sijacki D., 2021, MNRAS, 503, 3568
  * Kroupa (2001) Kroupa P., 2001, MNRAS, 322, 231
  * Krumholz & Gnedin (2011) Krumholz M. R., Gnedin N. Y., 2011, ApJ, 729, 36
  * Krumholz et al. (2018) Krumholz M. R., Burkhart B., Forbes J. C., Crocker R. M., 2018, MNRAS, 477, 2716
  * Kumamoto et al. (2017) Kumamoto J., Baba J., Saitoh T. R., 2017, PASJ, 69, 32
  * Lacey (1984) Lacey C. G., 1984, MNRAS, 208, 687
  * Lang et al. (2020) Lang P., et al., 2020, ApJ, 897, 122
  * Leaman et al. (2017) Leaman R., et al., 2017, MNRAS, 472, 1879
  * Lehnert et al. (2009) Lehnert M. D., Nesvadba N. P. H., Le Tiran L., Di Matteo P., van Driel W., Douglas L. S., Chemin L., Bournaud F., 2009, ApJ, 699, 1660
  * Leitherer et al. (1999) Leitherer C., et al., 1999, ApJS, 123, 3
  * Leroy et al. (2008) Leroy A. K., Walter F., Brinks E., Bigiel F., de Blok W. J. G., Madore B., Thornley M. D., 2008, AJ, 136, 2782
  * Ma et al. (2017) Ma X., Hopkins P. F., Wetzel A. R., Kirby E. N., Anglés-Alcázar D., Faucher-Giguère C.-A., Kereš D., Quataert E., 2017, MNRAS, 467, 2430
  * Mackereth et al. (2019a) Mackereth J. T., et al., 2019a, MNRAS, 482, 3426
  * Mackereth et al. (2019b) Mackereth J. T., et al., 2019b, MNRAS, 489, 176
  * Madau & Dickinson (2014) Madau P., Dickinson M., 2014, ARA&A, 52, 415
  * Majewski et al. (2017) Majewski S. R., et al., 2017, AJ, 154, 94
  * Margalef-Bentabol et al. (2022) Margalef-Bentabol B., Conselice C. J., Haeussler B., Casteels K., Lintott C., Masters K., Simmons B., 2022, MNRAS, 511, 1502
  * Martig et al. (2014) Martig M., Minchev I., Flynn C., 2014, MNRAS, 443, 2452
  * McElroy et al. (2022) McElroy R., et al., 2022, MNRAS, 515, 3406
  * Mercedes-Feliz et al. (2023) Mercedes-Feliz J., et al., 2023, arXiv e-prints, p. arXiv:2301.01784
  * Minchev & Quillen (2006) Minchev I., Quillen A. C., 2006, MNRAS, 368, 623
  * Moreno et al. (2021) Moreno J., et al., 2021, MNRAS, 503, 3113
  * Muratov et al. (2015) Muratov A. L., Kereš D., Faucher-Giguère C.-A., Hopkins P. F., Quataert E., Murray N., 2015, MNRAS, 454, 2691
  * Naab & Ostriker (2017) Naab T., Ostriker J. P., 2017, ARA&A, 55, 59
  * Nordström et al. (2004) Nordström B., et al., 2004, A&A, 418, 989
  * Orr et al. (2018) Orr M. E., et al., 2018, MNRAS, 478, 3653
  * Orr et al. (2020) Orr M. E., et al., 2020, MNRAS, 496, 1620
  * Orr et al. (2022) Orr M. E., et al., 2022, arXiv e-prints, p. arXiv:2209.14159
  * Pandya et al. (2021) Pandya V., et al., 2021, MNRAS, 508, 2979
  * Pessa et al. (2023) Pessa I., et al., 2023, A&A, 673, A147
  * Planck Collaboration et al. (2020) Planck Collaboration et al., 2020, A&A, 641, A6
  * Quinn et al. (1993) Quinn P. J., Hernquist L., Fullagar D. P., 1993, ApJ, 403, 74
  * Quirk et al. (2022) Quirk A. C. N., et al., 2022, AJ, 163, 166
  * Renaud et al. (2013) Renaud F., et al., 2013, MNRAS, 436, 1836
  * Rizzo et al. (2020) Rizzo F., Vegetti S., Powell D., Fraternali F., McKean J. P., Stacey H. R., White S. D. M., 2020, Nature, 584, 201
  * Roškar et al. (2010) Roškar R., Debattista V. P., Brooks A. M., Quinn T. R., Brook C. B., Governato F., Dalcanton J. J., Wadsley J., 2010, MNRAS, 408, 783
  * Saha et al. (2010) Saha K., Tseng Y.-H., Taam R. E., 2010, ApJ, 721, 1878
  * Sanderson et al. (2020) Sanderson R. E., et al., 2020, ApJS, 246, 6
  * Santistevan et al. (2020) Santistevan I. B., Wetzel A., El-Badry K., Bland-Hawthorn J., Boylan-Kolchin M., Bailin J., Faucher-Giguère C.-A., Benincasa S., 2020, MNRAS, 497, 747
  * Santistevan et al. (2021) Santistevan I. B., Wetzel A., Sanderson R. E., El-Badry K., Samuel J., Faucher-Giguère C.-A., 2021, MNRAS, 505, 921
  * Sellwood (2013) Sellwood J. A., 2013, ApJ, 769, L24
  * Shapley (2011) Shapley A. E., 2011, ARA&A, 49, 525
  * Silva Aguirre et al. (2018) Silva Aguirre V., et al., 2018, MNRAS, 475, 5487
  * Soderblom (2010) Soderblom D. R., 2010, ARA&A, 48, 581
  * Somerville & Davé (2015) Somerville R. S., Davé R., 2015, ARA&A, 53, 51
  * Sparre et al. (2017) Sparre M., Hayward C. C., Feldmann R., Faucher-Giguère C.-A., Muratov A. L., Kereš D., Hopkins P. F., 2017, MNRAS, 466, 88
  * Spitzer & Schwarzschild (1951) Spitzer Lyman J., Schwarzschild M., 1951, ApJ, 114, 385
  * Springel & Hernquist (2003) Springel V., Hernquist L., 2003, MNRAS, 339, 289
  * Stern et al. (2021) Stern J., et al., 2021, ApJ, 911, 88
  * Stinson et al. (2006) Stinson G., Seth A., Katz N., Wadsley J., Governato F., Quinn T., 2006, MNRAS, 373, 1074
  * Stott et al. (2016) Stott J. P., et al., 2016, MNRAS, 457, 1888
  * Strömberg (1946) Strömberg G., 1946, ApJ, 104, 12
  * Su et al. (2017) Su K.-Y., Hopkins P. F., Hayward C. C., Faucher-Giguère C.-A., Kereš D., Ma X., Robles V. H., 2017, MNRAS, 471, 144
  * Sun et al. (2020) Sun J., et al., 2020, ApJ, 901, L8
  * Swinbank et al. (2012) Swinbank A. M., Smail I., Sobral D., Theuns T., Best P. N., Geach J. E., 2012, ApJ, 760, 130
  * Tacconi et al. (2020) Tacconi L. J., Genzel R., Sternberg A., 2020, ARA&A, 58, 157
  * Teyssier (2002) Teyssier R., 2002, A&A, 385, 337
  * Ting & Rix (2019) Ting Y.-S., Rix H.-W., 2019, ApJ, 878, 21
  * Trapp et al. (2022) Trapp C. W., et al., 2022, MNRAS, 509, 4149
  * Übler et al. (2019) Übler H., et al., 2019, ApJ, 880, 48
  * Villalobos & Helmi (2008) Villalobos Á., Helmi A., 2008, MNRAS, 391, 1806
  * Vogelsberger et al. (2020) Vogelsberger M., Marinacci F., Torrey P., Puchwein E., 2020, Nature Reviews Physics, 2, 42
  * Wadsley et al. (2004) Wadsley J. W., Stadel J., Quinn T., 2004, New A, 9, 137
  * Wechsler et al. (2002) Wechsler R. H., Bullock J. S., Primack J. R., Kravtsov A. V., Dekel A., 2002, ApJ, 568, 52
  * Weinberger et al. (2017) Weinberger R., et al., 2017, MNRAS, 465, 3291
  * Wellons et al. (2023) Wellons S., et al., 2023, MNRAS, 520, 5394
  * Wetzel & Garrison-Kimmel (2020) Wetzel A., Garrison-Kimmel S., 2020, GizmoAnalysis: Read and analyze Gizmo simulations (ascl:2002.015)
  * Wetzel et al. (2023) Wetzel A., et al., 2023, ApJS, 265, 44
  * Wibking & Krumholz (2023) Wibking B. D., Krumholz M. R., 2023, MNRAS, 521, 5972
  * Wielen (1977) Wielen R., 1977, A&A, 60, 263
  * Wisnioski et al. (2015) Wisnioski E., et al., 2015, ApJ, 799, 209
  * Xiang & Rix (2022) Xiang M., Rix H.-W., 2022, Nature, 603, 599
  * Yu et al. (2021) Yu S., et al., 2021, MNRAS, 505, 889
  * Yu et al. (2022) Yu S., et al., 2022, arXiv e-prints, p. arXiv:2210.03845

##  Appendix A Effect of Vertical Selection

By default, we select stars currently at $|Z|<3~{}\text{kpc}$, which includes
stars in both the thin and thick disk, to sample disk stars at all ages in an
unbiased manner and avoid assumptions about thin- versus thick-disk dynamics.
Figure 11 shows the effect of the vertical selection on $\sigma_{\rm Z,now}$
and $\sigma_{\rm tot,now}$ versus stellar age, specifically, the fractional
change relative to our fiducial selection of $|Z|<3~{}\text{kpc}$. We also
examined the effect on $\sigma_{\rm R,now}$ and $\sigma_{\rm\phi,now}$ and
they show trends nearly identical to $\sigma_{\rm tot,now}$ so we do not show
them separately.

At most ages, our vertical selection has minimal impact on any component of
the dispersion, $\lesssim 3\%$. This is because most disk stars lie well
within $|Z|<1~{}\text{kpc}$, so they dominate regardless of vertical
selection. Thinner vertical selection does yield slightly higher $\sigma_{\rm
tot,now}$, although with significant scatter. Interestingly, thinner selection
only yields higher $\sigma_{\rm Z,now}$ for stars older than $\approx
8~{}\text{Gyr}$, and instead yield slightly lower values of $\sigma_{\rm
Z,now}$ for younger stars. The youngest population formed near the mid-plane
and experienced $\approx 125~{}\text{Myr}$ of post-formation heating, on
average, so most are still near the mid-plane and have low dispersions. Any
that are above $|Z|=0.25-1~{}\text{kpc}$ likely experienced an atypical amount
of heating. Again, all of these effects are at the level of $\lesssim 3\%$,
where non-equilibrium effects like disk warping and asymmetries may become
important.

If instead we measure the dispersion for the galaxy overall (not in age bins)
as in Section 7, larger vertical selections non-trivially increase the
dispersion. Vertical dispersions within $|Z|<3~{}\text{kpc}$ are 10 - 35%
larger than within $|Z|<1~{}\text{kpc}$, because larger vertical selections
include a higher fraction of older, higher-dispersion stars.

Figure 11:  Effect of vertical selection on $\sigma_{\rm Z}$ (top) and
$\sigma_{\rm tot}$ (bottom): present-day velocity dispersion for stars
selected within 4 vertical ranges, divided by the dispersion for our fiducial
range of $|Z|<3~{}\text{kpc}$, versus stellar age. We select all stars at
$R=6-10~{}\text{kpc}$. Thinner vertical selections lead to slightly higher
$\sigma_{\rm Z}$ and $\sigma_{\rm tot}$ for stars older than $\approx
8~{}\text{Gyr}$, but slightly lower $\sigma_{\rm Z}$ and higher $\sigma_{\rm
tot}$ for stars younger. However, different vertical selections change the
dispersion by $\lesssim 3\%$, with the exception of the oldest and youngest
stars, which have $\approx 5\%$ lower dispersions for thinner vertical ranges.

##  Appendix B Effect of In-Situ Selection

By default, we select only stars that are at $R=6-10~{}\text{kpc}$ and
$|Z|<3~{}\text{kpc}$ at present. We include all of these stars when measuring
present-day kinematics, but when measuring kinematics at formation, we select
only in-situ stars that formed within 30 kpc comoving of the galactic center.
Figure 12 shows the effect of also including this in-situ selection on the
current velocity dispersion versus stellar age. As above, we only show
$\sigma_{\rm Z,now}$ and $\sigma_{\rm tot,now}$, because $\sigma_{\rm R,now}$
and $\sigma_{\rm\phi,now}$ exhibit the same trends as $\sigma_{\rm tot,now}$.
Overall, our in-situ selection has no effect on stars younger than $\approx
6~{}\text{Gyr}$, because almost all of them formed in situ. However, the
effects of in-situ selection become increasingly significant for older stars,
when a more significant fraction formed ex situ. Not surprisingly, including
those older ex-situ stars increases the dispersion. That said, this effect is
$\lesssim 10\%$ for all but the oldest stars (see El-Badry et al., 2018b, for
analysis of these ancient stars).

Figure 12:  Effect of excluding ex-situ stars on the present-day velocity
dispersion: $\sigma_{\rm Z,now}$ (top) and $\sigma_{\rm tot,now}$ (bottom) of
only stars that formed in situ divided by that of all (in-situ + ex-situ)
stars, versus age, averaged across 11 galaxies. Stars with ages $\lesssim
6~{}\text{Gyr}$ show no significant difference, because almost all young disk
stars formed in situ. Older than 6 Gyr, the dispersion of in-situ stars
relative to all stars decreases with age, because the fraction of ex-situ
stars increases with age, and ex-situ stars have larger present-day
dispersions. However, at almost all ages that we investigate, the difference
in velocity dispersion is $\lesssim 10\%$.

##  Appendix C Effect of Cosmic-Ray Feedback

Figure 13:  Effects of cosmic-ray feedback on disk dynamical history: current
stellar kinematics versus age, averaged across 8 isolated MW-mass galaxies
simulated with fiducial FIRE-2 physics (solid) and additionally including
cosmic ray (CR) injection and transport from stars (dotted), assuming a
constant effective diffusion coefficient $\kappa=3\times 10^{29}\rm
cm^{2}s^{-1}$. The shaded regions show the galaxy-to-galaxy standard
deviation. Top: median azimuthal velocity, $v_{\rm\phi}$. The addition of CR
feedback reduces overall star formation, leading to lower-mass galaxies with
lower velocities at all ages. Middle: vertical velocity dispersion,
$\sigma_{\rm Z}$. Given their lower stellar mass, the CR simulations have
lower velocity dispersions, though the CR simulations converge to the fiducial
simulations for the youngest stars. Bottom: median $v_{\rm\phi}/\sigma_{\rm
Z}$, as a measure of rotational support and ‘diskiness’. Given the reduction
of both $v_{\rm\phi}$ and $\sigma_{\rm Z}$ in the CR simulations,
$v_{\rm\phi}/\sigma_{\rm Z}$ is similar at nearly all ages for the fiducial
and CR simulations. For the youngest stars, the CR simulations are
fractionally dynamically hotter and less rotationally supported. We find
similar trends for the other velocity components and the total dispersion.

The fiducial FIRE-2 simulations that we examine do not include self-consistent
cosmic ray (CR) injection and transport. However, CRs may be important in
regulating the dynamical state of the ISM, and therefore the resultant stellar
populations, because in the solar neighborhood the CR energy density is
comparable to the gas thermal and turbulent energy densities (for example
Boulares & Cox, 1990).

Here, we explore the effects of CR feedback from stars on the dynamical state
of our MW-mass disks. One proposed way to form kinematically colder disks in
cosmological simulations is the inclusion of previously neglected feedback
channels, including CRs, because they can provide a spatially and temporally
smoother form of feedback and potentially launch cooler galactic winds.
Previous work using FIRE-2 simulations found that the inclusion of magnetic
fields, anisotropic conduction, and viscosity has no appreciable effects on
any bulk galaxy properties, but that the inclusion of CRs, within reasonable
uncertainty for the treatment of the CR diffusion coefficient, significantly
can alter galaxy properties (Su et al., 2017; Hopkins et al., 2020). FIRE-2
simulations with CRs also include anisotropic CR transport with stream and
advection / diffusion (with constant parallel diffusivity, $\kappa=3\times
10^{29}\rm cm^{2}s^{-1}$), CR cooling (hadronic and Compton, adiabatic, and
streaming losses), CR injection in supernova shocks, and CR–gas coupling (see
Chan et al., 2019, for details of the numerical implementation of CRs in these
simulations).

Hopkins et al. (2020) showed in the FIRE-2 simulations that the primary effect
of CR feedback is to prevent the accretion of halo gas into the galaxy,
leading to lower star formation and stellar mass by $\approx 2-4\times$ in MW-
mass galaxies at $z\lesssim 2$. Furthermore, Chan et al. (2021) found that the
inclusion of CRs in FIRE-2 leads to lower velocity dispersions in gas, which
might suggest dynamically cooler disks, although they did not examine
$v_{\rm\phi}/\sigma$. We show here that the opposite is true for young stars
when considering $v_{\rm\phi}/\sigma$, that is, the inclusion of CRs leads to
dynamically hotter stellar disks.

Figure 13 shows the effect of CRs on present-day stellar kinematics versus
age. We include only the isolated galaxies from the Latte suite, which have
versions with CR feedback. We also include the lower-mass galaxies m12r and
m12z, because our goal is only to examine the differential effects of
including CRs. Figure 13 (top) shows the current median azimuthal velocity,
$v_{\rm\phi,now}$, versus stellar age, averaged across the same galaxies in
the CR simulations and the fiducial simulations. At all ages, galaxies with CR
feedback have lower $v_{\rm\phi}$, with the largest absolute difference for
the youngest stars, which arises primarily because of their lower stellar
masses.

Figure 13 (middle) shows the present vertical velocity dispersion,
$\sigma_{\rm Z}$, versus age. The CR simulations have lower $\sigma_{\rm Z}$
at all ages, except for the youngest stars, which may arise from the effects
of CRs on the ISM and halo gas occurring over long timescales. We also
examined $\sigma_{\rm R}$, $\sigma_{\rm\phi}$, and $\sigma_{\rm tot}$ and
found qualitatively similar results. Overall, the inclusion of CR feedback
leads to lower velocity dispersions in stars, consistent with analyses of gas
in these simulations (Hopkins et al., 2020; Chan et al., 2021).

Figure 13 (bottom) shows the sample-averaged $v_{\rm\phi}/\sigma_{\rm Z}$
versus age. For most ages, stars in the CR simulations have slightly less
rotational support than in the simulations without CR feedback; the youngest
stars show the most significant reduction at $\approx 30\%$. Thus, including
CR feedback leads to fractionally dynamically hotter galaxies, despite their
lower absolute velocity dispersions. However, Figure 14 shows that CR feedback
leads to galaxies with $\approx 3\times$ lower stellar mass at present, as
Hopkins et al. (2020) showed.

Figure 14 (top row) shows $v_{\rm\phi}/\sigma_{\rm Z}$ for all stars versus
stellar mass, while Figure 14 (bottom row) shows the same but for young stars
(age < 250 Myr). We show our 6 fiducial LG-like galaxies (as black stars) to
clarify our fiducial sample’s mass-relation across this entire mass range,
given that none of our fiducial isolated simulations have masses within
$2-4\times 10^{1}0\,\text{M}_{\odot}$, the range containing half of their CR
counterparts. That said, we advise against making direct comparisons between
any individual fiducial LG-like simulations and CR isolated simulations,
because each galaxy has its own unique history and our fiducial LG-like
simulations generally have slightly colder kinematics than our isolated
simulations. We emphasize that Figure 14 shows the same trends with stellar
mass as for our fiducial simulations without CR feedback in Figure 7. Thus,
the inclusion of CR feedback primarily acts to reduce a galaxy’s stellar mass
and correspondingly shifts its $v_{\rm\phi}$ and $v_{\rm\phi}/\sigma_{\rm Z}$
along the existing relation. This leads to galaxies that are slightly less
disky, with slightly lower $v_{\rm\phi}/\sigma_{\rm Z}$. We conclude that the
inclusion of CR feedback, at least as implemented in FIRE-2 and assuming a
constant diffusion coefficient, does not lead to dynamically colder disks.

Figure 14:  Effects of cosmic-ray feedback on current stellar
$v_{\rm\phi}/\sigma_{\rm Z}$ versus stellar mass, across 8 isolated MW-mass
galaxies simulated with fiducial FIRE-2 physics (circles) and additionally
with CR injection and transport (triangles). Black stars show LG-like galaxies
simulated with fiducial physics; although these galaxies lack CR counterparts,
we include them in order to illustrate our fiducial sample’s
$v_{\phi}/\sigma_{\rm Z}$ \- stellar mass relation. CR feedback leads to
galaxies with $\approx 3$ times lower stellar mass. Top row:
$v_{\phi}/\sigma_{\rm Z}$ for stars of all ages. The CR simulations have lower
$\sigma_{\rm Z}$ but also lower stellar masses, leading them to land on a
similar relation in $v_{\phi}/\sigma_{\rm Z}$. Bottom row: same but for young
stars with ages $<250~{}\text{Myr}$. CR simulations have slightly lower
$v_{\phi}/\sigma_{\rm Z}$ for 7 of the 8 galaxies, though again, they remain
on a similar scaling relation. Adding CR feedback does not lead to dynamically
colder disks; it primarily reduces the stellar mass formed and moves galaxies
along a similar relation.

  *[MW]: