††thanks: Corresponding author

# A First Principles Approach to Trust-Based Recommendation Systems

Paras Stefanopoulos<EMAIL_ADDRESS>School of Computing, The
Australia National University, Canberra, Australia  Ahad N. Zehmakan
<EMAIL_ADDRESS>School of Computing, The Australia National
University, Canberra, Australia  Sourin Chatterjee<EMAIL_ADDRESS>Department of Mathematics and Statistics, Indian Institute of Science
Education and Research, Kolkata, West Bengal 741246, India

###### Abstract

This paper explores recommender systems in social networks which leverage
information such as item rating, intra-item similarities, and trust graph. We
demonstrate that item-rating information is more influential than other
information types in a collaborative filtering approach. The trust graph-based
approaches were found to be more robust to network adversarial attacks due to
hard-to-manipulate trust structures. Intra-item information, although sub-
optimal in isolation, enhances the consistency of predictions and lower-end
performance when fused with other information forms. Additionally, the
Weighted Average framework is introduced, enabling the construction of
recommendation systems around any user-to-user similarity metric.

Trust graph; Cold-start recommendation; Adversarial attacks; Similarity
measures

##  I Introduction

Recommendation systems have become an essential component of online platforms,
helping users navigate vast amounts of information and make decisions based on
personalized suggestions. In social networks, these systems play a key role in
suggesting products, content, friends, or pages to users based on their
preferences, behaviors, and social connections. Social networks provide a rich
source of information on user interactions and preferences, which can be
leveraged to create more accurate and personalized recommendations. This type
of recommendation system considers not only the past behaviors and preferences
of users but also the behaviors and preferences of their connections or
friends in the network.

Trust is a crucial concept in social networks and plays a vital role in
recommendation systems [1, 2, 3, 4, 5]. The abstraction of social connections
as ’trust’ is essential because it allows the system to consider implicit
behaviors and preferences of users and their connections. For example, a user
liking or sharing content, or having many mutual connections, can indicate a
level of trust or shared interest, even if the users are not explicitly
connected or have not directly interacted. This abstraction helps in capturing
the nuances of social interactions and enables the recommendation system to
provide more personalized and accurate suggestions. Additionally, the concept
of trust incorporates a user’s social circle’s opinions and preferences,
making recommendations more robust and comprehensive [6]. In this paper, we
investigate trust connections, and the implications they have on a user’s
preferences and the robustness to adversarial attacks.

There are several types of recommendation systems, including content-based,
collaborative filtering, hybrid, and others. Content-based recommendation
systems suggest items based on the characteristics of the items and the
preferences of the users. On the other hand, collaborative filtering
recommends items based on the preferences and behaviors of similar users [7].
Hybrid systems combine both approaches. Trust-based collaborative filtering is
particularly interesting because it not only considers the preferences of
similar users but also incorporates the trust relationships between users.
This approach helps in mitigating some of the limitations of traditional
collaborative filtering, such as the cold start problem and data sparsity, by
leveraging the trust network to generate recommendations even when there is
limited interaction data available [8]. This is a technique we will explore in
this paper, but with a first principles approach to understand how to utilize
these different methods and information forms best.

Outline. The paper is organized in the following manner. Section II gives an
overview of prior work, and section III outlines our contribution. After that,
in section IV we introduce the methodology and in section V we discuss the
single information-based algorithms and subsequently multiple information-
based algorithms in section VI. Validation of the findings on more datasets is
given in section VII and the impact of adversarial attacks is explored in
section VIII. A discussion of our main findings is given in section IX.

##  II Literature Review

###  II.1 Traditional Approaches

Traditional trust-based recommendation systems mainly focus on incorporating
explicit and implicit trust relationships into the recommendation process. For
example, some approaches consider the ratings or preferences of a user’s
trusted friends when generating recommendations. [9] delved into a trust-aware
collaborative filtering approach. Their method distinctively pivots on
explicit trust relationships. This means that the system directly harnesses
expressed trust values (like friendship connections) to modulate the influence
of a user’s peers on their recommendations.

Other methods incorporate trust propagation mechanisms to infer trust
relationships between users who are not directly connected. The idea is
grounded in the notion that trust can be transitive to some degree [10]. For
example, if A trusts B and B trusts C, even if A hasn’t interacted with C,
there might be a certain level of implicit trust A has for C due to the mutual
connection. Numerous methods have been developed based on this trust
transitivity principle. Some models utilize matrix factorization techniques to
propagate trust across the user-item matrix, while others employ graph-based
algorithms to trace the spread of trust across a network of users. Such
propagation methods attempt to fill the gaps in the trust matrix, offering a
more comprehensive and interconnected trust landscape that can be leveraged
for more informed and nuanced recommendations [11].

###  II.2 Cold-Start Approaches

Addressing the “cold start” problem remains a pivotal challenge for trust-
based recommendation systems. The “cold start” issue surfaces when a new user
joins the platform or a new item is added, and there’s minimal to no
interaction data available for meaningful recommendation generation. In such
scenarios, traditional collaborative filtering methods often fall short, as
they rely predominantly on historical interaction data [12]. Various solutions
have been proposed, including adaptive clustering and incremental learning
[13], and a demographic approach based on finding similarities between old and
new users [14].

Another porpular approach is leveraging the trust network [15]. The assumption
here is that even if a user hasn’t interacted with many items, their trusted
connections might have. By leveraging this network, systems can extrapolate
the potential interests of the new user based on the preferences of their
trusted peers. One particularly noteworthy approach to leverage the trust
network against the “cold start” problem was conceived by [16]. They
introduced TrustWalker, a random walk model grounded in the trust network. It
commences its journey from the target user and propagates through their
trusted connections, and with some probability, swapping to other similar
items that may have more information. In this way, TrustWalker harnesses both
explicit and implicit trust connections along with skipping to different
items, which lets “cold start” users with little information still receive
“acceptable” recommendations.

###  II.3 Machine Learning Approaches

Machine learning, over the past decade, has revolutionized numerous domains
with its ability to extract patterns from vast data sets. Traditional methods
of collaborative filtering were primarily linear. However, the dynamics of
user preferences, trust relationships, and item attributes often exhibit
nonlinear patterns that require more sophisticated models for accurate
representation [17]. Hence, machine learning is a widely used tool in
recommendation systems [18]. Bayesian [19, 20] and decision tree [21, 22]
approaches are mostly used due to their simplicity [23]. Moreover, gradient-
decent-based [24], matrix factorization-based [25], and hidden Markov model-
based [26] algorithms are also in use.

Neural collaborative filtering models underscore the potential of deep neural
networks in handling the collaborative filtering challenge on implicit
feedback data sets [27]. They propose three different ML-based models for
developing recommendation systems and show the models compete with other
state-of-the-art approaches. Moreover, these approaches have also been used to
tackle the cold-start problem [28, 29].

In parallel, the application of graph neural networks (GNNs) in recommender
systems has witnessed a significant surge, establishing GNNs as a state-of-
the-art approach in this area. [30] offer a comprehensive review of the
literature centered around GNN-based recommender systems. The motivation to
incorporate GNNs into recommender systems primarily stems from their ability
to capture the higher-order structure of the graph.

###  II.4 Adversarial Defense Approaches

Recommendation systems are vulnerable to adversarial attacks, subverting their
performance. These adversarial attacks can manifest in various forms: fake
user profiles, strategically crafted reviews, or even swarming behaviors
intended to boost or suppress particular items. A robust recommendation system
isn’t just about accurate recommendations; it’s also about ensuring its
resilience against these adversarial intrusions [31]. Even the most advanced
models are vulnerable to certain adversarial tactics [32]. However, the silver
lining is the evolution of defense mechanisms. Techniques such as adversarial
training, where the model is trained with perturbed data, have shown promise
in enhancing model robustness.

Low-rank defenses [33] against adversarial attacks are one type of adversarial
defense strategy used in recommendation systems. Low-rank defenses entail
transforming the user-item matrix into a low-rank matrix to decrease the
impact of adversarial attacks. This method has been demonstrated to be
effective in fighting against several sorts of adversarial attacks, including
poisoning attacks and model extraction attacks. In recommendation systems,
generative adversarial networks (GANs) can detect and fight against
adversarial attacks. GANs can be used to generate adversarial instances and
teach the recommendation model to defend against them. This technique
effectively defends against adversarial assaults such as poisoning and model
extraction [31]. Another study by [34] focused on how to accurately portray
the representative rating for an item to mitigate the risk of an adversary
adding fake ratings. The paper found that placing weights on network-related
factors leads to increased robustness to adversaries.

##  III Our Contribution

We aim to leverage trust-based social network data to develop recommendation
systems that not only provide accurate recommendations but also can handle the
aforementioned challenges around cold-start users and adversarial attacks. We
discover that incorporating the trust graph data can facilitate the
establishment of recommendation systems that are robust against adversarial
attacks and provide high-quality recommendations for cold starters.

This research takes a first principles approach to building a recommendation
system. We will develop novel approaches to combining different data forms in
recommendation systems, such as a simple, but effective approach for
determining item similarity and an adaptable framework for designing
recommenders based on user similarities. We build several recommendation
systems, starting by only accessing one type of data at a time. This will let
us then evaluate these recommendation systems and understand how derivatives
of different types of data act in various contexts.

To structure the development of these recommendation systems, the data can be
split into three sub-sets, which are explained and then visualized below. The
idea is, that if we can build recommendation systems that are limited to only
these fundamental elements of “information”, we can then take what we learn to
produce more effective recommenders that operate on combinations of these
subsets.

  * •

Trust Graph. The trust graph is simply the nodes and the trust connections
between one another. Recommendation systems built upon this sub-set of the
data only concern themselves with the structure of the social network and not
with the ratings on items.

  * •

Intra-Item Information. The intra-item information is a subset of the data
that contains the intra-item similarity, essentially just a table, providing
similarity scores between pairs of items along with each node’s neighbors.

  * •

Item-Rating Information. The Item-Rating information is a subset of the data
that only has nodes, items, and their ratings. Recommendation systems built
upon this subset draw similarities between nodes based on how they’ve jointly
rated items, not considering the actual social connection between them.

Figure 1: An example of three types of information. 4 individuals whose social
structure is shown via edges. Rating corresponding to each individual and
overall similarity metric among items provided using tables.

As visualized in Figure 2, we study algorithms in three levels. In level 1, we
investigate the algorithms which use one type of information. Building on our
findings in level 1, we consider all three possibilities of combining this
information in pairs in level 2. Finally, in level 3, we leverage all three
types of information to devise algorithms. This systematic study of the
problem will allow us to evaluate the strengths of various novel methods, such
as similarity indices-based approaches [35], random walks techniques [2], and
opinion formation models [36]. More precisely, we characterize the approaches
and types of information that facilitate more accurate performance and
robustness against adversarial attacks.

Figure 2: Schematic diagram of building a recommendation system based on
different subsets of data and combining them. In particular, cold-start
compatible approaches are outlined in red.

##  IV Methodology

###  IV.1 Primary Training Data Set

For the sake of iteration in our experiments, we used a single data set for
developing and testing the algorithms and then validated the results with
other datasets as well.

Epinions data set is one of the most commonly used datasets in building
recommendation systems. Epinions was a consumer review site, where users could
place reviews on items, and “trust” other users. Ratings were based on stars
and are hence in the range of 1-5. The trust network contains 49k users, with
480k “trust” edges between them. There are also 150k items, with 660k ratings
between them. This data set contains nodes that have at least 1 “in-edge”,
implying at least one other user trusts them. This network is distributed
similarly to other social networks, where node popularity exponentially
decays, and where there are a few trusted authorities who provide powerful and
centralized opinions.

Figure 3: Distribution of ratings for all items

By observing the ratings data (see Figure 3), we found that people’s ratings
are negatively skewed, i.e. users had a positive bias on their ratings
relative to the 5-star scale. Upon focusing on the top 10 items rated by the
largest number of users, we can observe negative skewness. Hence, being able
to accurately provide a recommendation will require interpolation from the
graph data and is not easily solved by simply assuming high ratings.

###  IV.2 Similarity Metrics

Here, we discuss the similarity metrics used in the design of the
recommendation systems.

Note that the trust relations between the users can be modeled using a graph
$G=(V,E)$, where the node set $V$ corresponds to the users, and the edge set
$E$ models the relations.

####  IV.2.1 Jaccard Index

One of the most popular similarity measures between nodes (users) in a graph
is the Jaccard index. Let $N(u)$ for a node $u$ denote the set of neighbors of
$u$. The Jaccard similarity index between two nodes $u_{i}$ and $u_{j}$ is
defined as:

$J(u_{i},u_{j})=\frac{|N(u_{i})\cap N(u_{j})|}{|N(u_{i})\cup N(u_{j})|}$

One drawback of this index is that if two nodes are adjacent but share no
other neighbors, then the value is 0. To mitigate this issue, we redefine the
index in the following way:

$J^{\prime}(u_{i},u_{j})=J(u_{i},u_{j})+\begin{cases}\frac{1}{|N(u_{i})\cup
N(u_{j})|}&u_{j}\in N(u_{i})\\\ 0&\text{otherwise}\end{cases}$

####  IV.2.2 Jaccard-Item Index

Jaccard-Item Index is a form of the Jaccard index, but instead of relating
users using their neighbors, we relate them based on the items they have
rated.

$J_{I}(u_{i},u_{j})=\frac{|I(u_{i})\cap I(u_{j})|}{|I(u_{i})\cup I(u_{j})|}$

where $J_{I}$ is the Jaccard-Item Index and $I(u)$ is the items rated by node
$u$.

####  IV.2.3 Item-Rating Difference

The Item-Rating difference compares two nodes as follows:

$R_{D}(u_{i},u_{j})=1-\frac{\sum_{k\in I(u_{i})\cap
I(u_{j})}{|R(u_{i},k)-R(u_{j},k)|}}{4\cdot|I(u_{i})\cap I(u_{j})|}$

$R(u_{i},k)$ denotes the rating received by the k-th item from user $u_{i}$.
we are effectively enumerating the items that two users have both rated. For
each item, we take the absolute difference in their rating. We then divide by
4 to normalize (maximum possible difference in ratings) and divide by the
magnitude of the intersection to obtain the average.

####  IV.2.4 Pearson Correlation

One similarity metric used to determine the relationship between two items is
the Pearson correlation of the ratings given to them. This was proposed by
[16].

The similarity function between two items, i and j, works as follows:

$S(i,j)=\text{set of users who have rated both items }i\text{ and }j.$

$\rho(i,j)=\frac{\sum_{u\in
S(i,j)}{(R(u,i)-\overline{R}(u))(R(u,j)-\overline{R}(u))}}{\sqrt{(\sum_{u\in
S(i,j)}{(R(u,i)-\overline{R}(u))^{2}})(\sum_{u\in
S(i,i)}{(R(u,i)-\overline{R}(u))^{2})}}}$

We then want to skew this correlation, such that if $|S(i,j)|$ is large, we
are more confident. So, the final similarity metric, which is equivalent to
that given in [16], is:

$sim(i,j)=\frac{1}{1+e^{-\frac{|S(i,j)|}{2}}}\times\rho(i,j)$

The sigmoid function is used to balance how much we favor the magnitude of $S$
and to ensure our similarity metric is in the range $[0,1]$.

Thus, for any pair of items, using just the Intra-Item graph, we can determine
a quantified value that indicates how “related” the pair of items is.

####  IV.2.5 Intra-Item Jaccard

In this similarity metric, to relate items we use Jaccard. i.e. the number of
users who have rated $i$ and $j$ divided by the union of the users who rated
both $i,j$. So the similarity between the two items is:

$J_{II}(i,j)=\frac{|S(i,j)|}{|S(i)\cup S(j)|}$

where $S(i)$ and $S(j)$ are the set of nodes rated $i$ and $j$, respectively.

###  IV.3 Evaluating Recommendation Systems

The information needed by recommendation systems is quite sparse. Considering
the number of items and graph size and often ratings concentrated on only
popular items, computing recommendations for every item, is very costly.
Hence, we will be focusing on an approach that will give results relatively
quickly. To make recommendations quickly, and to build an environment where
there is substantive data, we reduce the data sets to a subset where we keep
the top $k$ items and nodes that have rated these items. Now the subset
contains active users, who all have at least one item in common (the most
popular item), and it is likely that they also share a number of the other top
items due to their popularity.

To quantify the performance, we take $15\%$ of the users at random, who are
called “test users” and now for each item $i$, we remove its true rating and
then compute the recommendation for the item. After that, we calculate the
binary accuracy (correct discrete rating or not) and mean absolute error (MAE)
[37] across all the test users. This process is repeated 5 times and for each
time, we record the MAE (this is important since some of our algorithms are
randomized). For simplicity, we assume the recommender’s performance follows a
normal distribution. We compute the standard deviation and mean of MAE to gain
an idea of how reliable the recommender is.

As a reference point for designed recommenders, we use two simple strategies
as control recommenders.

  1. 1.

Random: The random recommender simply provides a random rating of 1-5 for each
item.

  2. 2.

Universal Random: The universal random recommender takes the distribution of
ratings into account by generating a rating proportional to the ratings given
to the item.

Both recommenders on the Epinions data set, give moderate results and will be
good reference points for future evaluation.

Table 1: Mean and Standard Deviation of MAE for Baseline Algorithms Algorithm | $\mu_{MAE}$ | $\sigma_{MAE}$
---|---|---
Random | 1.67 | 0.21
Universal Random | 1.30 | 0.38
Figure 4: Distribution of MAE for Control Recommenders

##  V Understanding single level information

In this section, we focus on the algorithms which use only one of the three
types of information available.

###  V.1 Trust Graph

The most simple approach to take with the trust graph is, given a node, to
investigate its neighbors. We can infer the current node’s rating by measuring
its neighbors’ central tendencies like mean, median, and mode.

####  V.1.1 Centrality of Neighbours

Let’s say we are rating item i, we go through a node’s neighbors and calculate
the mean, median, or mode of their ratings for item i (if it exists). If no
ratings for i exist, we return a random rating. Note, this algorithm mutates
the underlying data source, whereby a recommended value for a node, can be
used by another node in its rating inference. The performance of these
recommenders can be seen in Table 2. We observe that the mean and median of
neighbors were the best algorithms in this range. The mode version did not
perform as well. This is because mode cannot produce decimal-precision ratings
it also defaults to randomness when a mode can’t be determined.

Table 2: Mean and Standard Deviation of MAE for Neighbourhood Algorithms Algorithm | $\mu_{MAE}$ | $\sigma_{MAE}$
---|---|---
Mean of Neighbours | 1.27 | 0.25
Median of Neighbours | 1.27 | 0.19
Mode of Neighbours | 1.32 | 0.23
Universal Random | 1.30 | 0.38
Random | 1.67 | 0.21

####  V.1.2 Jaccard-Weighted Model

The Jaccard-Weighted Model extends the neighbor-based models from before.
Effectively, given the neighbors of node $u$, their rating of an item is
considered based on their modified Jaccard Index ($J^{\prime}$) with node $u$.
Finally, we identify the rating that received the maximum ‘score’. The
performance of this algorithm, along with the previous Median of Neighbours
model can be seen in Table 3. We observe that the Jaccard-Weighted model
performs quite similarly to the Median of Neighbours model when it comes to
mean accuracy.

Table 3: Mean and Standard Deviation of MAE for Jaccard Weighted Neighbours Algorithm | $\mu_{MAE}$ | $\sigma_{MAE}$
---|---|---
Median of Neighbours | 1.27 | 0.19
Jaccard Weighted Neighbours | 1.31 | 0.32
Universal Random | 1.30 | 0.38
Random | 1.67 | 0.21

####  V.1.3 Monte-Carlo Random Walk

We now move to investigate recommendation systems that consider the entire
trust graph. The Monte-Carlo (MC) Random Walk model works by starting at the
node we are trying to produce a rating for and then performing a random walk.
As the distance from the original node, $d$, increases, so does the
probability that a random rating is returned. This random walk continues until
a node with a rating is reached, or this random rating is produced. We run
this experiment $k$ times, for some fixed $k$, for each node and take the mean
rating from the nodes reached via random walk.

We then extended this algorithm by choosing the neighbor, at each step of the
random walk based on a weighted random incorporating the Jaccard index. The
performance of both MC algorithms, as compared with the others can be seen in
Table 4. We observe that the Monte Carlo algorithm outperforms the previous
models. Interestingly, the Jaccard weighted version performed better. This
implies that Jaccard indexes seem to hold more valuable information than mere
trust connections.

Table 4: Mean and Standard Deviation of MAE for Random Walk Algorithms Algorithm | $\mu_{MAE}$ | $\sigma_{MAE}$
---|---|---
Jaccard Monte-Carlo Random Walk | 1.14 | 0.22
Monte-Carlo Random Walk | 1.16 | 0.23
Median of Neighbours | 1.27 | 0.19
Universal Random | 1.30 | 0.38
Random | 1.67 | 0.21

####  V.1.4 Jaccard Majority of Majorities (JMoM) and Weighted Average
(Jaccard WA)

These recommenders work by considering every other node, which has rated item
i, in the graph. We then perform a summation of the Jaccard values for every
node corresponding to each rating value. In the Majority of Majorities
version, we then take the rating that has received the highest contribution of
Jaccard scores. However, in the Weighted Average version, we instead take the
weighted average of these ratings, weighted by the contribution of Jaccard
scores. We will call these techniques JMoM and WA moving forward. (These
methods are inspired by the popular majority based opinion diffusion models,
cf. [38, 39, 40].)

As can be observed from Table 5, the JMoM performed reasonably well. Its mean
accuracy is better than the random walk. However, it was less robust. The WA
version performed very well significantly beating both the JMoM version and
the random walk. By overcoming the constraint of only considering direct
neighbors, the JMoM and Jaccard WA models can produce decent results, with a
very simple algorithm. Note that these approaches can be used with any user-
to-user similarity metric.

Table 5: Mean and Standard Deviation of MAE for Full-graph Jaccard Algorithms Algorithm | $\mu_{MAE}$ | $\sigma_{MAE}$
---|---|---
Jaccard WA | 1.05 | 0.25
Jaccard Monte-Carlo Random Walk | 1.14 | 0.22
Jaccard MoM | 1.13 | 0.30
Universal Random | 1.30 | 0.38
Random | 1.67 | 0.21

###  V.2 Item Rating Information

In the Item Jaccard model, we used the Jaccard-Item Index $J_{I}$ which is
based on items rated. We also used the weighted average framework approach but
with the Item-Rating Difference. The performance of these alternative WA
models can be seen in Table 6.

Table 6: Mean and Standard Deviation of MAE for Item Rating Algorithms Algorithm | $\mu_{MAE}$ | $\sigma_{MAE}$
---|---|---
Item-Jaccard WA | 1.02 | 0.23
Jaccard WA | 1.05 | 0.25
Item-Rating Difference WA | 1.17 | 0.49
Universal Random | 1.30 | 0.38
Random | 1.67 | 0.21

We can see that the Item-Jaccard form of WA outperforms the previous Jaccard
WA. Item-Rating Difference doesn’t seem to perform very well and is quite
inconsistent in its performance. A potential reason for the Item-Jaccard WA
model performing so well is that a user’s social circle may not be the best
predictor for who is similar to them, rather relating people by their items
themselves does a better job.

###  V.3 Intra-Item Information

The Intra-Item Information concerns itself with the relationships between
items themselves. In a similar fashion to the previous algorithms, we
manipulated the WA approach to work with Intra-Item similarity data. This is
slightly different from the WA framework as we iterate through each user’s
items rather than using a user-user similarity metric. The performance of the
algorithms, one based on $J_{II}$ and one based on the Pearson $sim$, can be
seen in Table 7.

Table 7: Mean and Standard Deviation of MAE for Intra-Item Algorithms Algorithm | $\mu_{MAE}$ | $\sigma_{MAE}$
---|---|---
Item-Jaccard WA | 1.02 | 0.23
Intra-Item WA | 1.24 | 0.14
Universal Random | 1.30 | 0.38
Intra-Item WA (Pearson) | 1.34 | 0.13
Random | 1.67 | 0.21

It is evident both of the intra-item-based models underperform our previous
models. The intra-item similarity using the Jaccard approach outperforms the
Pearson correlation approach. This is a very interesting result, potentially
this could be a better metric for determining the similarity of items and
could be applied in [16]. It is evident the intra-item similarity conveys some
information as we perform better than random, hence this could be useful in
combination models in the coming sections.

###  V.4 Summary

By splitting the data up in this way, we now understand how algorithms focused
on the different types of information perform. In a way, we know the value of
the different types of information. Not very surprisingly, our item-based
similarity metrics resulted in better-performing algorithms. However, as we
will see in the rest of the paper, leveraging the other types of data
facilitate the design of algorithms which are more robust against adversarial
attacks and support cold start users.

##  VI Combining Multiple Information

###  VI.1 Jaccard Intra-Item WA

As previously discussed, both Trust and Intra-Item algorithms are effective on
cold start users, and so their combination. In an attempt to combine both
these modes of data, we developed an intersection between the Intra-Item WA
($J_{II}$) which operates on neighbors, and the Jaccard WA which operates on
the entire graph. The idea of this algorithm is for every (user, item) pair,
we contribute the Jaccard index of the user multiplied by the item’s
similarity to our rating bins. The performance can be seen in Table 8.

Table 8: Mean and Standard Deviation of MAE for Jaccard Intra-Item WA Algorithm Algorithm | $\mu_{MAE}$ | $\sigma_{MAE}$
---|---|---
Item-Jaccard WA | 1.02 | 0.23
Jaccard Intra-Item WA | 1.20 | 0.13
Universal Random | 1.30 | 0.38
Random | 1.67 | 0.21

Combining these two forms of information seem to produce a very consistently
performing recommender. However, it did not perform very well relative to our
current best, Item-Jaccard WA in terms of mean MAE.

###  VI.2 Weighted Item-Rating Difference WA (WIRD WA)

This model combines the Intra-Item graph with the Item graph. Thus combing
information about how nodes rated individual items and how the items are
interrelated. The general idea was to slightly alter the item rating
similarity metric, $R_{D}$, to provide a bias taking into account the
similarity of items to the item we are recommending.

$WR_{D}(u_{i},u_{j},k)=1-\frac{\sum_{j\in
I(u_{j})}{J_{II}(k,j)\times|R(u_{i},k)-R(u_{j},j)|}}{4\cdot|I(u_{j})|}$

Note for intra-item similarity we use $J_{II}$ rather than $sim$ as we found
it was more performant in the previous section. The results can be seen in
Table 9:

Table 9: Mean and Standard Deviation of MAE for WIRD WA Algorithm Algorithm | $\mu_{MAE}$ | $\sigma_{MAE}$
---|---|---
WIRD WA | 1.05 | 0.28
Item-Jaccard WA | 1.02 | 0.23
Intra-Item WA | 1.34 | 0.13
Item-Rating Difference WA | 1.17 | 0.49
Universal Random | 1.30 | 0.38
Random | 1.67 | 0.21

Though the WIRD model performed quite well still not better than the Item-
Jaccard WA model. However, it performs better than just the intra-item
information or item-rating difference WA alone. Hence, we have successfully
combined the intra-item and item-rating information to achieve a better
result.

###  VI.3 Jaccard Item-Jaccard WA

In this model, we combine the Trust information with the Item-Rating
information. To do so, we explored different ways of combining the Jaccard WA
and Item-Jaccard WA models. Both similarity metrics are in the range $[0,1]$,
so we attempted linear combinations, multiplying them and using the $max$
function.

Thus, the metrics would be combined as:

$I_{S}(u_{i},u_{j})=\alpha\times J^{\prime}(u_{i},u_{j})+(1-\alpha)\times
J_{I}(u_{i},u_{j}),\alpha\in[0,1]$

Table 10: Mean and Standard Deviation of MAE for Jaccard Item-Jaccard WA Algorithm Algorithm | $\mu_{MAE}$ | $\sigma_{MAE}$
---|---|---
Jaccard Item-Jaccard WA (Addition, $\alpha=0.5$) | 1.00 | 0.26
Jaccard Item-Jaccard WA (Maximum) | 1.01 | 0.26
Jaccard Item-Jaccard WA (Multiplication) | 1.04 | 0.22
Item-Jaccard WA | 1.02 | 0.23
Universal Random | 1.30 | 0.38
Random | 1.67 | 0.21

Looking at the multiplication-based model, it produced an interesting result
which was mirrored in the Jaccard and $J_{II}$ combined from the Trust and
Intra-Item section. The accuracy was essentially a blend between the Jaccard
model and the Item-Difference one. The $max$-based model performed well but
not as well as addition. This implies that even small values of Jaccard or
Item-Jaccard contribute information to the recommender that the addition model
captures. So, we find that the Jaccard Item-Jaccard WA combination model with
equal contribution is the best performer.

###  VI.4 All Combined Algorithm

Now, we will combine all the data we have in an attempt to build our best-
performing model. From previous experiments, we found that the Jaccard cross
Item-Jaccard Weighted Average model performed the best. This incorporates the
trust information and the rating information. How can we extend this to
incorporate intra-item information? To do this, we defined the similarity
between two given users as $U_{sim}=J(u_{1},u_{2})+I_{j}(u_{1},u_{2})$. Let’s
say we are predicting a rating for $u_{1}$, item $I$. We would iterate through
all other nodes in the graph, and for each of the items they’ve rated,
contribute $U_{sim}\times J_{II}(I,j)$ where $j$ is the other item.
Essentially, other users influence our user’s ratings by how similar they are
based on the social network structure, how similar they are based on the items
they’ve rated, and how similar each of their items is to the one we are
predicting. Interestingly, the combination model under-performed our best
model in terms of mean MAE.

##  VII Validation With Other Data Sets

Throughout this paper, we have designed our models and tested them using the
Epinions data set. This is a data set popular in literature and gives us a
good playground. However, to test the universality of the algorithms, we have
tested them on the following two extra datasets.

  1. 1.

FilmTrust [41]

It was a movie review platform, where users could ‘trust’ one another’s
reviews. This data set was scraped from the site in 2011 and is of common use
for recommendation system research.

  2. 2.

CiaoDVD [42]

This is again another movie review website, with a similar structure to
FilmTrust. This data set was scraped in 2013.

Table 11: Mean and Standard Deviation of MAE for All Algorithms on Different Data Sets Algorithm | Epinions | FilmTrust | CiaoDVD
---|---|---|---
$\mu_{MAE}$ | $\sigma_{MAE}$ |  | $\mu_{MAE}$ | $\sigma_{MAE}$ |  | $\mu_{MAE}$ | $\sigma_{MAE}$
Jaccard Item-Jaccard WA | 1.00 | 0.26 |  | 0.66 | 0.08 |  | 0.53 | 0.28
Item-Jaccard WA | 1.02 | 0.23 |  | 0.67 | 0.08 |  | 0.58 | 0.32
Jaccard WA | 1.05 | 0.25 |  | 1.14 | 0.08 |  | 1.73 | 0.36
WIRD WA | 1.05 | 0.28 |  | 0.67 | 0.04 |  | 0.72 | 0.37
Jaccard Item-Jaccard JII Combination WA | 1.07 | 0.22 |  | 0.69 | 0.06 |  | 0.63 | 0.30
JWIRD WA | 1.09 | 0.27 |  | 0.67 | 0.07 |  | 0.72 | 0.40
Jaccard MoM | 1.13 | 0.30 |  | 1.19 | 0.09 |  | 1.77 | 0.37
Jaccard Monte-Carlo Random Walk | 1.14 | 0.22 |  | 1.20 | 0.08 |  | 1.81 | 0.49
Monte-Carlo Random Walk | 1.16 | 0.23 |  | 1.20 | 0.08 |  | 1.82 | 0.35
Item-Rating Difference WA | 1.17 | 0.49 |  | 0.67 | 0.08 |  | 0.79 | 0.41
Jaccard Intra-Item WA | 1.20 | 0.13 |  | 1.07 | 0.06 |  | 1.75 | 0.32
Intra-Item WA | 1.24 | 0.14 |  | 1.18 | 0.08 |  | 1.70 | 0.45
Median of Neighbours | 1.27 | 0.19 |  | 1.26 | 0.08 |  | 1.76 | 0.28
Mean of Neighbours | 1.27 | 0.25 |  | 1.25 | 0.07 |  | 1.67 | 0.43
Universal Random | 1.30 | 0.38 |  | 0.89 | 0.11 |  | 0.72 | 0.58
Jaccard Weighted Neighbours | 1.31 | 0.32 |  | 1.22 | 0.10 |  | 1.66 | 0.44
Mode of Neighbours | 1.32 | 0.23 |  | 1.23 | 0.09 |  | 1.67 | 0.43
Intra-Item WA (Pearson) | 1.34 | 0.13 |  | 1.21 | 0.07 |  | 1.67 | 0.50
Random | 1.67 | 0.21 |  | 1.34 | 0.07 |  | 1.82 | 0.31

It appears the Jaccard Item-Jaccard WA recommendation algorithm performed best
across all data sets with the lowest MAE value. The algorithm consistently
outperformed other algorithms in each of the data sets - Epinions, FilmTrust,
and CiaoDVD - with MAEs of 1.00, 0.66, and 0.53 respectively. This algorithm’s
Mean MAE ($\mu_{MAE}$) was significantly lower than other algorithms.

Several observations can be made from the results in Table 11:

  * •

As expected, the Random algorithm had the highest MAE in every experiment.

  * •

The performance of algorithms varied across data sets. For instance, the
Jaccard WA algorithm had a relatively low MAE for Epinions but performed
poorly on the FilmTrust and CiaoDVD data sets.

  * •

Algorithms that combined different recommendation techniques, like Jaccard
Item-Jaccard JII Combination WA, did not necessarily guarantee better
performance than their simpler counterparts.

##  VIII Robustness Against Adversarial Attacks

We will now investigate how our recommendation systems respond to adversarial
attacks on the network. To conduct this experiment we must construct an attack
method and mutate the graph to represent it. The attack we will simulate is
that of an actor who wants to artificially reduce the ratings of certain items
on the system. The actor can create fake “bot” accounts and befriend an
adversarial celebrity. We choose the 10th most popular node by in-edges as a
celebrity node and then create some number of fake nodes, who rate every item
a $1$ out of $5$. These nodes will trust the celebrity node and will be
trusted by the celebrity too.

To test the robustness of our algorithms, we will use the “Epinions” data set
and use 150 fake accounts (around 5% of the graph by node size). We deem this
to be quite a significant attack on the integrity of the network. We only
tested the algorithms that performed reasonably well. Algorithms tested here
are as follows: Jaccard Item-Jaccard WA, Item-Jaccard WA, Jaccard WA, WIRD WA,
Jaccard Item-Jaccard JII Combination WA, and Jaccard Monte-Carlo Random Walk.

Table 12: Algorithm Performance Metrics Under Adversarial Conditions Algorithm | Normal | Adversarial | Reduction (%)
---|---|---|---
| $\mu_{MAE}$ | $\sigma_{MAE}$ | $\mu_{MAE}$ | $\sigma_{MAE}$ | $\mu_{MAE}$
Jaccard Item-Jaccard WA | 1.00 | 0.26 | 1.25 | 0.19 | 24.93
Item-Jaccard WA | 1.02 | 0.23 | 1.24 | 0.20 | 21.77
Jaccard WA | 1.05 | 0.25 | 1.17 | 0.22 | 12.12
WIRD WA | 1.05 | 0.28 | 1.32 | 0.23 | 25.53
Jaccard Item-Jaccard JII Combination WA | 1.07 | 0.22 | 1.23 | 0.20 | 15.58
Jaccard Monte-Carlo Random Walk | 1.14 | 0.22 | 1.21 | 0.19 | 6.39
Figure 5:  MAE distributions for different recommender algorithms. The dotted
lines depict the performance with an adversarial network and the solid lines
show the same recommender under normal conditions.

In these results, we can observe that the minimally impacted algorithms were
those that relied on the trust information. These are the Jaccard WA model and
the Jaccard Monte-Carlo Random Walk. The most impacted were those that relied
on Item-Rating and Intra-Item information, i.e. Item-Jaccard WA and WIRD WA.
Interestingly, the combination model, which made use of all forms of data,
seemed to not suffer as big a reduction in performance.

##  IX Discussion

We successfully broke down the concept of a recommender system to its core
principles, starting with the information it ingests. Through breaking down
the different information forms, we noticed that the most important feature
when assigning the similarity between two users for their preferences is the
way they’ve rated other items and which items they have engaged with. The
trust and intra-item-based algorithms led to recommenders with higher MAE as
seen in Table 11. All other models included the item-rating information in
some way, which had significantly lower MAE. This result makes logical sense
as it’s more likely we are similar in opinion to a person who we share
interests with than a person who we share friends with. However, the problem
with item-rating information is the information that we lack when cold-
starting a user.

In the challenging case of cold start users, the algorithms that only use the
Trust and Intra-Item information would be very beneficial, but these
algorithms, unsuperisingly, perform not as well. Among these, we found that
the Jaccard WA approach provided the best performance, outperforming random
walk and intra-item-based approaches. A further idea that could be explored is
simply using an approach like the Universal Random, which was one of our
baseline recommenders. This recommender performed reasonably when averaged
across our data sets and predicts preferences based on the distribution of
current ratings for said item.

While combining information, we found it wasn’t easy to outperform the Item-
Jaccard-based recommender. The one case where we managed to improve upon the
base Item-Jaccard model was when we performed addition with the Jaccard Index.
So the similarity between the two users is based on the sum of the Jaccard
Index and Item-Jaccard Index. We think that the performance of the recommender
is increased as the Jaccard provides some marginal information which can
impact the ratings in the case when a user trusts a celebrity. That way the
user has a non-zero Jaccard index with many other nodes who also trust said
celebrity, hence, slightly shifting the rating prediction in the direction of
these users.

We attempted combination models of various forms, combining the similarity
scores with multiplication, summation, and a maximum function. However, we
discovered that scaling the distributions had minimal impact on the resulting
performance. We also concluded that simply using an equally weighted summation
of the two scores leads to the simplest and most effective model. Hence, the
highest-performing recommendation model combined trust information and item
rating information. Though, it is not suitable for cold-start users.

We also found that recommenders based on the WA framework outperformed random-
walk methods. Furthermore, the WA framework outperforms a similar approach
called Majority of Majorities. We believe this approach to collaborative
filtering is effective because it is not constrained by the structure of the
graph (as is a random walk) and can produce fine-grained recommendations that
can be non-integer values. Whereas, the MoM approach forces a specific integer
rating to be assigned and hence is not as representative of the user’s
preference.

When it came to the intra-item information, we found it performed quite poorly
from the perspective of $\mu_{MAE}$. However, an interesting observation, as
per Table 11, was that the algorithms that included intra-item information
were the most consistent in their performance—featuring the lowest
$\sigma_{MAE}$ across all data sets. Thus, it can be deduced that the intra-
item information is additive from a stability perspective, making a
recommender perform with similar accuracy for all users. This idea was
reinforced when we created the fully combined model, which utilized all
information types. This was the “Jaccard Item-Jaccard JII Combination WA”
model. We found this model performed worse than the Jaccard Item-Jaccard WA
model in $\mu_{MAE}$, however, it performed better in mean $\sigma_{MAE}$.
Hence, in a scenario where providing consistently good predictions for all
users is of importance, the introduction of intra-item information could
facilitate this.

Through our experimentation with opposing intra-item similarity metrics, we
determined that the Intra-Item Jaccard approach outperformed the Pearson
Correlation metric on downstream tasks. This was determined by building models
upon these metrics and comparing the resulting accuracy. When further testing
with the Intra-Item Jaccard similarity metric was undertaken, the resulting
models performed better than random, and in combination led to stable models
as described earlier. A recurring theme in these experiments was the very
impressive performance of the Jaccard Index when used in a variety of
applications. It was the go-to for drawing similarity scores using all
information forms and is a simple but logical way to reason about the
similarity of sets.

Testing the performance of our algorithms on an adversarial network provided
some interesting insights. The most obvious of which is that item rating
information is highly susceptible to an attack with fake accounts. The trust
graph is relatively robust to adversaries, as one must influence individuals
to shift the dynamics of this information. This is represented in the data by
the recommendation systems involving the item rating information suffering
from the most severe reductions in their MAE. We found that the least impacted
algorithm was the MC Random Walk. As the algorithm propagates from a user,
through the network to derive the rating, it will arrive at a legitimate node
with a rating before venturing into the adversaries. The only case where the
adversaries have an impact is when no ratings are found and the randomness
means a rating is taken uniformly from the network. In these cases, the volume
of ratings from adversaries impacts results. The next best algorithm for being
robust in the face of adversaries was Jaccard WA. This is a Weighted Average
based entirely on the Jaccard Index. This algorithm performs better due to the
value it places on the trust graph, very rarely do nodes have a non-zero
Jaccard similarity with the bad actors, hence the high performance.

Future Work: In investigating intra-item information, we found that the
”Intra-Item Jaccard” metric was more effective than Pearson Correlation in
determining item similarity. This result suggests further exploration of the
Trust-Walker [16], which may lead to improved performance.

Further research is needed to improve the computational efficiency of the WA
approach. Unlike, random-walk methods, the WA-based recommenders require
iterating over all nodes in the graph. Optimizations such as vectorization or
caching similarities between users could be explored.

Furthermore, this work focused on fake accounts and fake ratings, but other
attacks such as bribing popular nodes or adversarial censorship should be
investigated in future studies.

## Authors’ Contribution

Paras Stefanopoulos: Conceptualization, Formal Analysis, Methodology,
Software, Validation, Writing - Original draft preparation. Ahad N. Zehmakan:
Conceptualization, Methodology, Supervision, Validation, Writing - Review.
Sourin Chatterjee: Software, Validation, Writing - Original draft preparation.

## Declaration of Competing Interest

The authors declare that they have no known competing financial interests or
personal relationships that could have appeared to influence the work reported
in this paper.

## References

  * O’Donovan and Smyth [2005] J. O’Donovan and B. Smyth, Proceedings of the 10th international conference on Intelligent user interfaces - IUI ’05 (2005), https://doi.org/10.1145/1040830.1040870.
  * Andersen  _et al._ [2008] R. Andersen, C. Borgs, J. Chayes, U. Feige, A. Flaxman, A. Kalai, V. Mirrokni, and M. Tennenholtz, in  _Proceedings of the 17th international conference on World Wide Web_ (2008) pp. 199–208.
  * Zehmakan  _et al._ [2023] A. N. Zehmakan, C. Out, and S. H. Khelejan, arXiv preprint arXiv:2305.08558 (2023).
  * Htun and Tar [2013] Z. Htun and P. P. Tar, Int. J. Innov. Eng. Technol.(IJIET) Technol.(IJIET) 2, 271 (2013).
  * Liu  _et al._ [2023] C. Liu, X. Zhou, A. N. Zehmakan, and Z. Zhang, IEEE Transactions on Knowledge and Data Engineering (2023).
  * Sherchan  _et al._ [2013] W. Sherchan, S. Nepal, and C. Paris, ACM Computing Surveys 45, 1–33 (2013).
  * Su and Khoshgoftaar [2009] X. Su and T. M. Khoshgoftaar, Advances in artificial intelligence 2009 (2009).
  * Yang  _et al._ [2017] B. Yang, Y. Lei, J. Liu, and W. Li, IEEE Transactions on Pattern Analysis and Machine Intelligence 39, 1633–1647 (2017).
  * Massa and Avesani [2007] P. Massa and P. Avesani, Proceedings of the 2007 ACM conference on Recommender systems - RecSys ’07 (2007), https://doi.org/10.1145/1297231.1297235.
  * Richters and Peixoto [2011] O. Richters and T. P. Peixoto, PloS one 6, e18384 (2011).
  * Jamali and Ester [2010] M. Jamali and M. Ester, Proceedings of the fourth ACM conference on Recommender systems - RecSys ’10 (2010), https://doi.org/10.1145/1864708.1864736.
  * Moghaddam and Elahi [2019] F. B. Moghaddam and M. Elahi, in  _Big data recommender systems: Recent trends and advances_ (IET, 2019).
  * Sanjeevi  _et al._ [2015] S. G. Sanjeevi, S. Grandhe, and H. Shrivastav, Int J Comput Sci Electron Eng (IJCSEE) 3 (2015).
  * Pandey and Rajpoot [2016] A. K. Pandey and D. S. Rajpoot, in  _2016 International Conference on Signal Processing and Communication (ICSC)_ (IEEE, 2016) pp. 213–218.
  * Lam  _et al._ [2008] X. N. Lam, T. Vu, T. D. Le, and A. D. Duong, Proceedings of the 2nd international conference on Ubiquitous information management and communication - ICUIMC ’08 (2008), https://doi.org/10.1145/1352793.1352837.
  * Jamali and Ester [2009] M. Jamali and M. Ester, in  _Proceedings of the 15th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining_, KDD ’09 (Association for Computing Machinery, New York, NY, USA, 2009) p. 397–406.
  * Marlin [2004] B. Marlin,  _Collaborative filtering: A machine learning perspective_ (University of Toronto Toronto, 2004).
  * Khanal  _et al._ [2020] S. S. Khanal, P. Prasad, A. Alsadoon, and A. Maag, Education and Information Technologies 25, 2635 (2020).
  * Ericson and Pallickara [2013] K. Ericson and S. Pallickara, Future Generation Computer Systems 29, 1024 (2013).
  * Felden and Chamoni [2007] C. Felden and P. Chamoni, in  _2007 40th Annual Hawaii International Conference on System Sciences (HICSS’07)_ (IEEE, 2007) pp. 168a–168a.
  * Shinde  _et al._ [2013] A. Shinde, M. Haghnevis, M. A. Janssen, G. C. Runger, and M. Janakiram, International Journal of Innovation and Technology Management 10, 1340019 (2013).
  * Lucas  _et al._ [2012] J. P. Lucas, S. Segrera, and M. N. Moreno, Expert Systems with Applications 39, 1273 (2012).
  * Portugal  _et al._ [2018] I. Portugal, P. Alencar, and D. Cowan, Expert Systems with Applications 97, 205 (2018).
  * Chen  _et al._ [2013] H. Chen, Y. Tang, L. Li, Y. Yuan, X. Li, and Y. Tang, IEEE transactions on cybernetics 43, 898 (2013).
  * Schelter  _et al._ [2013] S. Schelter, C. Boden, M. Schenck, A. Alexandrov, and V. Markl, in  _Proceedings of the 7th ACM Conference on Recommender Systems_ (2013) pp. 281–284.
  * Fang [2022] Y. Fang, in  _2022 IEEE International Conference on Electrical Engineering, Big Data and Algorithms (EEBDA)_ (IEEE, 2022) pp. 1209–1213.
  * He  _et al._ [2017] X. He, L. Liao, H. Zhang, L. Nie, X. Hu, and T.-S. Chua, in  _Proceedings of the 26th international conference on world wide web_ (2017) pp. 173–182.
  * Wei  _et al._ [2017] J. Wei, J. He, K. Chen, Y. Zhou, and Z. Tang, Expert Systems with Applications 69, 29 (2017).
  * Volkovs  _et al._ [2017] M. Volkovs, G. Yu, and T. Poutanen, Advances in neural information processing systems 30 (2017).
  * Gao  _et al._ [2023] C. Gao, Y. Zheng, N. Li, Y. Li, Y. Qin, J. Piao, Y. Quan, J. Chang, D. Jin, X. He, and Y. Li, ACM Transactions on Recommender Systems (2023), https://doi.org/10.1145/3568022.
  * Deldjoo  _et al._ [2021] Y. Deldjoo, T. D. Noia, and F. A. Merra, ACM Computing Surveys 54, 1–38 (2021).
  * Anelli  _et al._ [2021] V. W. Anelli, Y. Deldjoo, T. DiNoia, and F. A. Merra, Recommender Systems Handbook , 335–379 (2021).
  * Entezari and Papalexakis [2022] N. Entezari and E. E. Papalexakis, in  _2022 IEEE International Conference on Big Data (Big Data)_ (IEEE, 2022) pp. 5708–5714.
  * Lian and Zhao [2022] X. Lian and D. Zhao, Springer eBooks , 216–227 (2022).
  * Bag  _et al._ [2019] S. Bag, S. K. Kumar, and M. K. Tiwari, Information Sciences 483, 53 (2019).
  * Zehmakan [2021] A. N. Zehmakan, in  _Proceedings of the AAAI Conference on Artificial Intelligence_ , Vol. 35 (2021) pp. 5611–5619.
  * Bobadilla  _et al._ [2013] J. Bobadilla, F. Ortega, A. Hernando, and A. Gutiérrez, Knowledge-based systems 46, 109 (2013).
  * Out and Zehmakan [2021] C. Out and A. N. Zehmakan, arXiv preprint arXiv:2109.14265 (2021).
  * Auletta  _et al._ [2020] V. Auletta, D. Ferraioli, and G. Greco, Artificial Intelligence 284, 103288 (2020).
  * Zehmakan [2020] A. N. Zehmakan, Discrete Applied Mathematics 277, 280 (2020).
  * Guo  _et al._ [2013] G. Guo, J. Zhang, and N. Yorke-Smith, in  _Proceedings of the 23rd International Joint Conference on Artificial Intelligence (IJCAI)_ (2013) pp. 2619–2625.
  * Guo  _et al._ [2014] G. Guo, J. Zhang, D. Thalmann, and N. Yorke-Smith, in  _Proceedings of the 2014 International Conference on Advances in Social Networks Analysis and Mining (ASONAM)_ (2014).


# Cosmic ray density variations in nearby giant molecular clouds

Jiahao Liu  Bing Liu  and Ruizhi Yang

###### Abstract

In this paper, we analyzed 12 years of Fermi LAT gamma-ray data towards three
nearby giant molecular clouds (GMCs), i.e., R CrA, Chamaeleon, and Lupus. We
calibrated the gas column density of these regions by using the Planck dust
opacity map as well as the Gaia extinction map. With both the gamma-ray
observations and gas column density maps, we derived the cosmic ray densities
in the three GMCs. We found the derived CR spectra have almost the same shape
but significantly different normalizations, which may reflect that the
distributions of CRs in the vicinity of solar systems are inhomogeneous.

##  1 Introduction

Cosmic rays (CRs), as one of the most important ingredients in the
interstellar medium (ISM), dominate the heating and ionization of gas inside
molecular clouds, regulate the star-forming rate, and play a leading role in
the astrochemistry process. Composed of mainly charged particles, CRs are well
mixed in the Galactic magnetic field and lost the information of their
acceleration sites. Furthermore, it is debatable whether or not the CR density
measured in the solar system is a good representative of that in the whole
Galaxy. In this regard, gamma rays, the secondary products of interactions
between CRs and the ISM, provide us with crucial information about the
distribution of CRs in the Galaxy. Giant molecular clouds (GMCs), as the CR
barometers [1], are considered the ideal sites to study CRs. The gamma rays in
GMCs have already been extensively studied after the Launch of Fermi
satellite, and similar CR densities and spectra inside the GMCs in the Gould
Belt are derived [2, 3].

In this paper, we analyzed 12 years of Fermi Large Area Telescope (LAT)
observations towards three GMCs, including R CrA, Chamaeleon, and Lupus.
Compared with previous studies[3, 2, 4], in addition to the updated Fermi-LAT
analysis with more exposure, we also applied the method described above for a
new estimation of the gas mass, which may lead us to new results. Indeed, in
GMCs the gamma rays are mainly produced by the interaction of CRs with the
ambient gas. Thus, the gamma-ray spectra can be used to determine the spectral
shape of the parent CRs, meanwhile, the absolute normalization of CR density
depends also on the gas mass. In the former work, both CO and dust opacity
maps are used to trace the gas [3, 2, 4]. However, the CO emission may
underestimate the total gas due to the existence of the well-known CO dark gas
[5]. The dust opacity map is free of the ”dark gas” problem, but the ratio
between dust opacity and gas column may vary significantly in different
positions of the sky [6]. In this regard, the extinction map may be a more
reliable tracer of gas, since the dust absorption depends less on the
environment compared with the dust emission, which was used to derive the dust
opacity. However, in very dense regions such as the core of GMCs, extinction
may not be measured. Thus in this work, we derived the linear relation between
the Gaia extinction map [7] and the Planck dust opacity in different GMCs.
After that we derived the ratio between Planck dust opacity and gas column
density, assuming the ratio between extinction and gas column density is
constant. In this way, we derived more reliable gas column density maps, which
were used to study the CR density in different clouds.

##  2 Fermi-LAT Data Analysis

In order to study the gamma-ray emission of the GMCs, we collected 12 years of
Fermi-LAT Pass 8 data, which is the latest data release. We used the
Fermitools from Conda distribution and applied the latest version of the
instrument response functions ( $P8R3\\_SOURCE\\_V3$). We chose
$15^{\circ}\times 15^{\circ}$ regions centered at the positions of R CrA
(R.A.=285.89∘, Dec.=-36.50∘), Chamaeleon (R.A.=236.19∘, Dec.=-34.15∘), and
Lupus (R.A.=178.75∘, Dec.=-78.00∘) as the regions of interest (ROI). For Lupus
and R CrA, we performed the analysis under the equatorial coordinate system,
but for Chamaeleon, we chose the galactic coordinate due to the aberration
caused by the high declination.

We used gtselect to select the ”source” class events within the energy range
of 0.1-500 GeV to conduct our research, and applied the data filtrating
equation $(\text{DATA}\\_\text{QUAL}>0)\&\&(\text{LAT}\\_\text{CONFIG}==1)$ to
choose the good time intervals. Additionally, we only included the events
whose zenith angles are less than 90∘ to filter out the background
contamination from the Earth’s limb. The gamma-ray counts maps of the GMCs are
presented in Fig.1. We generated the source models using
make4FGLxml.py111https://fermi.gsfc.nasa.gov/ssc/data/analysis/user/make4FGLxml.py,
which consists of the sources in 4FGL-DR2 within the ROI enlarged by 7∘, the
galactic diffuse background emission (gll_iem_v07.fits), and the isotropic
emission background (iso_P8R3_SOURCE_V3_v1.txt).

(a) R CrA

(b) Chamaeleon

(c) Lupus

Figure 1: Observed gamma-ray counts maps of the GMCs, smoothed with Gaussian
filter of 0.3∘.

###  2.1 Spatial analysis and results

Because the galactic diffuse model mentioned above already included the gamma
rays produced in GMCs, so at first, we generated our own spatial templates to
separate the gamma rays associated with the GMCs from the galactic diffuse
background. Assuming that gamma rays trace the distribution of the molecular
gas, we derived the gas column density distribution using data from Planck
all-sky observation of thermal dust emission [6], and the details are
described in Sec.3. Here, we retained the data with gas column density higher
than $10^{21}$ cm-2 around the GMCs as the spacial templates of the GMCs (see
Fig.2). Meanwhile, the gas column density maps within the ROI enlarged by 10∘
in which the GMCs are subtracted are set as the background templates. Besides
the gamma rays produced from the pion-decay process induced by the inelastic
collision between CR protons and ambient gas, the diffuse emission also
includes the inverse Compton (IC) scattering of CR electrons in the
interstellar radiation fields (ISRFs). Thus, We generated the diffuse IC
emission template using GALPROP222http://galprop.stanford.edu/webrun/ [8] with
the parameter set ${}^{S}Y^{Z}4^{R}30^{T}150^{C}2$ in Ref. [9]. Moreover, we
used the log-parabola function ($dN/dE=N_{0}(E/E_{\rm b})^{-\alpha-\beta
log(E/E_{\rm b})}$) as the spectral model of the gamma-ray emission associated
with the gas templates.

(a) R CrA

(b) Chamaeleon

(c) Lupus

Figure 2: Gas column density maps of the GMCs (in the unit of ${\rm
cm}^{-2}$). Details are described in Sec.3.

Next, using gtlike, we performed likelihood analysis with the modified source
models, in which the galactic diffuse background emission is replaced by the
emission from the GMCs, the gas background, and the IC component. During the
fitting process, the normalizations and spectral parameters for diffuse
sources mentioned above as well as the sources with significance $\leq
5\,\sigma$ and angular distance $\geq 10^{\circ}$ to the center of each ROI
are set free. Based on the best-fit results, we generated the model maps using
gtmodel, and by subtracting the model maps from the observed gamma-ray counts
maps we obtained the residuals for the ROI. To check the fitting results
statistically, we then divided the residual counts by the square root of the
observed gamma-ray counts, i.e., the residual significance ((signal-to-noise,
in $\sigma$), of each pixel. As can be seen from Fig.3, for each ROI, the
significance of almost all the residual emissions is within $\pm 3\sigma$. In
other words, our gas templates can approximate the observed gamma-ray
distribution of the GMCs well, which allows for further analysis.

(a) R CrA

(b) Chamaeleon

(c) Lupus

Figure 3: Residual significance (signal-to-noise, S/N) maps of the GMCs,
smoothed with Gaussian filter of 0.3∘.

###  2.2 Spectral analysis and results

Figure 4: Gamma-ray spectral energy distributions of the GMCs. For each data
point, the error bar indicates the uncertainty caused by the statistics and
the systematic error. The dashed lines show the corresponding best-fit
spectra.

##  3 Gas column density estimation

To investigate the property of CRs within the clouds, we also need to estimate
the gas column density of the GMCs. The main components of the clouds, i.e.,
atomic and molecular hydrogen gases are usually traced by 21 cm H i line and
2.6 mm CO line emission, respectively. However, in many cases, only CO and H i
observation data are not enough to derive the entire mass of the neutral gas
[5]. Therefore, we chose the infrared emission from cold interstellar dust to
derive the gas column density, which is an independent measurement and will
take the untraceable ”dark” gas account.

To establish the relation between dust opacity and column density, we adopted
the Eq.(4) from [10], in which the dust opacity $\tau_{\rm M}$ is represented
as a function of wavelength $\lambda$,

$\tau_{\rm M}(\lambda)=\left(\frac{\tau_{\rm D}(\lambda)}{N_{\rm
H}}\right)^{\rm dust}[N_{{\rm HI}}+2X_{\rm CO}W_{\rm CO}].$ (3.1)

In the above equation, $(\tau_{\rm D}/N_{\rm H})^{\rm dust}$ is the reference
dust emissivity measured in low-$N_{\rm H}$ regions, and $X_{\rm CO}=N_{\rm
H_{2}}/W_{\rm CO}$ is the conversion factor between the integrated brightness
temperature of the CO ($W_{\rm CO}$) and the density of molecular hydrogen
($N_{\rm H_{2}}$). According to Eq.3.1, the total gas density ($N_{\rm H}$)
can be estimated as

$N_{\rm H}=N_{\rm{HI}}+2N_{\rm H_{2}}=\tau_{\rm
M}(\lambda)\left[\left(\frac{\tau_{\rm D}(\lambda)}{N_{\rm H}}\right)^{\rm
dust}\right]^{-1}.$ (3.2)

By substituting $(\tau_{\rm D}/N_{\rm H})^{\rm dust}_{353{\rm~{}GHz}}=(8.4\pm
3.0)\times 10^{-27}$ cm2 for the dust emissivity at 353 GHz (which is adopted
from the Table 4 of [10]) into Eq.3.2, we derived the gas column density maps
for the GMC regions. Then we set the data points as zero while the column
density $<10^{21}\,{\rm cm}^{-2}$ to generate spatial templates of the GMCs
for Fermi-LAT data analysis (see Sec.2.1). The gas column density distribution
maps of the GMCs are shown in Fig.2.

In previous studies (see, e.g., Ref. [3, 11]), the reference dust emissivity
was considered to be a constant value. However, such a ratio is not uniform in
our Galaxy, and the difference can be as high as 100% as shown in Fig.17 of
Ref. [6]. To calibrate the ratio between gas column density and dust opacity,
we adopted the linear relation between the optical depth and the extinction as
in Ref. [12], which is

$A_{V}=\gamma\tau_{850}+\delta$ (3.3)

where $A_{V}$ stands for the extinction in $V$-band, $\tau_{850}$ stands for
the optical depth at frequency $\nu=353$ GHz ($\lambda=850$ $\upmu$m) ,while
$\gamma$ and $\delta$ are constants. Assuming the extinction is proportional
to the gas column density, then we can get

$N_{\rm H}=N_{\rm HI}+2N_{\rm H_{2}}=\beta_{V}{A_{V}},$ (3.4)

in which $\beta_{V}=5.8\times 10^{21}\mathrm{cm}^{-2}\mathrm{mag}^{-1}$ is the
mean ratio of total neutral hydrogen to color excess [13].

Combining Eq.3.3 and Eq.3.4, now we can use the extinction to calibrate the
gas-to-dust ratio for each GMC. For the extinction, we used the Gaia total
galactic extinction map [14, 15], where $A_{0}$ is nearly equivalent to the
$V$-band extinction $A_{\rm V}$ within an accuracy of better than 3% for the
extinction range considered in this paper [16]. For optical depth, as
mentioned in the article [12], we only used the data points whose $\tau_{850}$
was in the range of $0-2\times 10^{-4}$. Both the extinction maps and optical
depth maps are generated with the same regions and division of pixels as that
of the spatial templates of the GMCs for the gamma-ray data analysis, i.e.,
Fig.2. Next, we extracted the data points of each pixel from the extinction
maps and optical depth maps, and applied maximum likelihood fitting to get the
value of $\gamma$ and $\delta$ in Eq.3.3 for the GMC regions. In addition, we
calculated $\bar{\tau}$, which is the average $\tau_{850}$ of each GMC region.
The best-fit $\gamma$, $\delta$ as well as $\bar{\tau}$ are listed in Table.1.
As can be seen from Table.1, although $\gamma$ is much larger than $\delta$,
the value of $\gamma\bar{\tau}$ and $\delta$ are of a similar order of
magnitude. Therefore, $\delta$ can not be ignored. The best-fit curves
representing Eq.3.3 as well as the data points of the extinction maps and the
dust opacity maps are plotted in Fig.5. Using the newly derived $\gamma$ and
$\delta$, we recalculated $A_{K}$, then generated the recalibrated gas column
density map by applying Eq.3.4 to the same dust opacity data.

Finally, using the two kinds of gas column density maps, which are derived
from a constant gas-to-dust ratio or recalibrated gas-to-dust ratio, we
calculated the corresponding total gas masses of the GMCs and listed the
results in Table1.

Table 1: Distance and the derived parameters of the GMCs

GMC | distance | $\gamma$ | $\delta$ | $\left(\tau_{\mathrm{D}}/N_{\mathrm{H}}\right)_{353\mathrm{GHz}}^{\text{dust }}$ | $\left(\tau_{\mathrm{D}}/N_{\mathrm{H}}\right)_{353\mathrm{GHz}}^{\text{dust }}$
---|---|---|---|---|---
| (pc) | ($\times 10^{4}$) | ($\times 10^{-2}$) | ($\times 10^{-27}$)(C) | ($\times 10^{-27}$)(R)
R CrA | 150 | $0.933\pm 0.013$ | $14.78\pm 0.38$ | $8.4\pm 3.0$ | $12.9\pm 0.2$
Chamaeleon | 215 | $1.213\pm 0.005$ | $8.38\pm 0.16$ | $8.4\pm 3.0$ | $12.5\pm 0.1$
Lupus | 190 | $1.358\pm 0.008$ | $7.58\pm 0.20$ | $8.4\pm 3.0$ | $11.3\pm 0.1$

The meaning of the fitted parameters, the constant (marked with ”C” in table)
and recalibrated (marked with ”R” in table) gas-to-dust
ratio$\left(\tau_{\mathrm{D}}/N_{\mathrm{H}}\right)_{353\mathrm{GHz}}^{\text{dust
}}$ are described in Sec.3.

(a) R CrA

(b) Chamaeleon

(c) Lupus

Figure 5: Relationship between the optical depth of given wavelength and the
extinction of the GMCs. The color bar shows the counts of data. The red line
shows the best-fit result of the linear relation.

##  4 CR properties in the GMCs

Assuming the diffuse gamma rays are mainly generated from the pp interaction
between CR protons and the gases, then we could derive the spectra of CR in
the GMCs, based on their gamma-ray spectra and the derived total gas masses.
Here, we adopted the gamma-ray production cross section from Ref. [17] for
calculation, assuming a simple power-law spectrum $F=f_{i}(E+E_{\rm
b})^{\Gamma}$ for the CR protons of each GMC. Here $F$ is the flux density of
protons in unit of $\mathrm{cm}^{-2}\,\mathrm{s}^{-1}\,\mathrm{GeV}^{-1}$, $E$
is the proton energy in unit of GeV, and $\Gamma$ stands for the spectral
index. During the fitting process, we fixed $E_{\rm b}$ to 1 GeV since the
current $\gamma$-ray data is not very constraining for the low energy part of
CR spectra. Applying the masses estimated from the two kinds of methods as
described in Sec.3, we derived two different normalizations respectively for
the CR spectra of the GMCs. For $f_{i}$(constant), we assumed a constant
$\left(\frac{\tau_{\rm D}(\lambda)}{N_{\rm H}}\right)^{\rm dust}$ while for
$f_{i}$(recalibrated) we used the values recalibrated for each GMC,
respectively. The fitted parameters of CR spectra are shown in Table.2.

Table 2: Fitting results of CR spectra for the GMCs

GMC | Mass | Mass | $f_{i}(\times 10^{-10})$ | $f_{i}(\times 10^{-10})$ | $\Gamma$
---|---|---|---|---|---
| (constant) | (recalibrated) | (constant) | (recalibrated) |
R CrA | 2231.3 | 1450.5 | $12.20\pm 1.39$ | $18.77\pm 2.10$ | $2.716\pm 0.051$
Chamaeleon | 14279.5 | 9621.2 | $19.48\pm 1.14$ | $28.89\pm 1.69$ | $2.817\pm 0.021$
Lupus | 3964.6 | 2957.6 | $11.31\pm 0.80$ | $15.19\pm 1.10$ | $2.868\pm 0.049$

The different methods of deriving the mass (in the unit of $M_{\odot}$) of
each GMC lead to different $f_{i}$, corresponding to the two panels in Fig.6.

(a) CR spectrum (constant)

(b) CR spectrum (recalibrated)

Figure 6: Derived CR spectra of the GMCs. The dashed lines show the best-fit
CR spectra, the shaded regions represent the 1 $\sigma$ errors including both
the statistics and the systematic uncertainty, and the red dots are the AMS-02
observation data acquired from Ref. [18]. The left panel shows the results
assuming a constant $\left(\frac{\tau_{\rm D}(\lambda)}{N_{\rm H}}\right)^{\rm
dust}$, and the right panel shows the results using the recalibrated values
(see Sec.3 for details) for each GMC.

Next, we plotted the derived CR spectra together with the local measurement by
AMS-02 [18] in Fig.6. We showed the results with both a constant ratio between
the gas column and dust opacity $\left(\frac{\tau_{\rm D}(\lambda)}{N_{\rm
H}}\right)^{\rm dust}$ (Fig.6(a)) and the ratios recalibrated in Sec.3
(Fig.6(b)). We found that for both cases the shape of the derived CR spectra
are in good consistency with the AMS-02 observation. But the total
normalizations show significant differences, even though we have taken into
account the $30\%$ uncertainty of the ratio between gas column and dust
opacity and the $20\%$ uncertainty in the $\gamma$-ray production cross
sections as systematic errors (shown as the shaded area in Fig.6).

##  5 Discussion

In this paper, we analyzed 12 years of LAT observation data towards three
nearby GMCs, and used both the dust opacity map and the extinction map to
derive the gas column density map in GMCs. We found the ratio between gas
column density and dust opacity can be quite different from cloud to cloud,
which may be due to the differences in dust composition and grain distribution
[12]. With both the gas distribution and gamma-ray measurement, we derived the
CR spectra in all three clouds. We found in all three clouds the spectral
shapes are consistent with the measurement of local CRs in the solar system,
but the absolute normalizations are significantly different.

The differences in the ratios between gas column density and dust opacity in
different GMCs have strong implications for studying CRs via gamma-ray
observations. Since the derived CR density from the gamma-ray flux would be
inversely proportional to the derived total mass, which is dependent on such
ratios. In the former studies [e.g., 3, 11], a uniform ratio between gas
column density and dust opacity is used in the Galaxy. In this work, we found
that the ratio for the nearby GMCs are similar to each other, but can have a
difference of about 50% compared with the value derived from thew whole sky.
Such difference can be even higher when we consider further GMCs, in which not
only dust composition and grain distribution but also the dust-to-gas ratio
can be quite different. Thus to derive the exact CR distribution in our
Galaxy, using only the dust opacity as the gas tracer may bring significant
systematic errors in the absolute normalizations of the CR spectra, and more
efforts are required to determine the gas distributions. We note that the
current results rely on the assumption that the dust-to-gas ratios are
identical in the three GMCs, which is reasonable considering the proximity of
all three clouds. But if such an assumption breaks down, one can not even use
the extinction map to calibrate the gas mass.

The different CR densities derived in this paper may reflect the inhomogeneous
distribution of CRs near the solar system. The height above the Galactic plane
for the three clouds is comparable, thus the difference cannot be attributed
to the CR gradient above the Galactic plane. On the other hand, we note that
the CR flux in Chamaeleon, which locates in the anti-Galactic Center (GC)
direction, is significantly higher than that in Lupus, which locates towards
the GC. Such a difference is similar to the CR anisotropy, whose phase points
to the anti-GC direction below 100 TeV[19]. Ref. [19] have proposed a nearby
CR source in the anti-GC direction may be responsible for such anisotropy. And
Ref. [20] further proposed that the Geminga SNR may be a promising candidate.
In this scenario, the higher CR flux in Chameleon can be explained by its
vicinity to the local CR source. Further modeling of the CR propagation may be
needed to confirm such a scenario.

## Acknowledgments

Bing Liu acknowledges the support of the NSFC under grant 12103049. Rui-zhi
Yang is supported by the NSFC under grants 12041305, 11421303 and the national
youth thousand talents program in China.

This work has made use of data from the European Space Agency (ESA) mission
Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing
and Analysis Consortium (DPAC,
https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has
been provided by national institutions, in particular, the institutions
participating in the Gaia Multilateral Agreement.

## References

  * [1] F.A. Aharonian, _Gamma Rays From Molecular Clouds_ , _Space Science Reviews_ 99 (2001) 187 [astro-ph/0012290].
  * [2] A. Neronov, D.V. Semikoz and A.M. Taylor, _Low-Energy Break in the Spectrum of Galactic Cosmic Rays_ , _Physical Review Letters_ 108 (2012) 051105 [1112.5541].
  * [3] Yang, Rui-zhi, de Oña Wilhelmi, Emma and Aharonian, Felix, _Probing cosmic rays in nearby giant molecular clouds with the fermi large area telescope_ , _A &A_ 566 (2014) A142.
  * [4] G. Peron, F. Aharonian, S. Casanova, R. Yang and R. Zanin, _Probing the Cosmic-Ray Density in the Inner Galaxy_ , _Astrophysical Journal Letters_ 907 (2021) L11 [2101.09510].
  * [5] I.A. Grenier, J.-M. Casandjian and R. Terrier, _Unveiling Extensive Clouds of Dark Gas in the Solar Neighborhood_ , _Science_ 307 (2005) 1292.
  * [6] Planck Collaboration, A. Abergel, P.A.R. Ade, N. Aghanim, M.I.R. Alves, G. Aniano et al., _Planck 2013 results. XI. All-sky model of thermal dust emission_ , _Astronomy & Astrophysics_ 571 (2014) A11 [1312.1300].
  * [7] L. Delchambre, C.A.L. Bailer-Jones, I. Bellas-Velidis, R. Drimmel, D. Garabato, R. Carballo et al., _Gaia DR3: Apsis III – Non-stellar content and source classification_ , _arXiv e-prints_ (2022) arXiv:2206.06710 [2206.06710].
  * [8] A.E. Vladimirov, S.W. Digel, G. Jóhannesson, P.F. Michelson, I.V. Moskalenko, P.L. Nolan et al., _GALPROP WebRun: An internet-based service for calculating galactic cosmic ray propagation and associated photon emissions_ , _Computer Physics Communications_ 182 (2011) 1156 [1008.3642].
  * [9] M. Ackermann, M. Ajello, W.B. Atwood, L. Baldini, J. Ballet, G. Barbiellini et al., _Fermi-LAT Observations of the Diffuse $\gamma$-Ray Emission: Implications for Cosmic Rays and the Interstellar Medium_, _ApJ_ 750 (2012) 3 [1202.4039].
  * [10] Planck Collaboration, A. Abergel, P.A.R. Ade, N. Aghanim, M.I.R. Alves, G. Aniano et al., _Planck 2013 results. XI. All-sky model of thermal dust emission_ , _Astron. Astrophys._ 571 (2014) A11 [1312.1300].
  * [11] R. Yang, F. Aharonian and C. Evoli, _Radial distribution of the diffuse $\gamma$ -ray emissivity in the Galactic disk_, _Phys. Rev. D_ 93 (2016) 123007 [1602.04710].
  * [12] E. Zari, M. Lombardi, J. Alves, C.J. Lada and H. Bouy, _Herschel-Planck dust optical depth and column density maps. II. Perseus_ , _Astronomy & Astrophysics_ 587 (2016) A106 [1511.08503].
  * [13] R.C. Bohlin, B.D. Savage and J.F. Drake, _A survey of interstellar H I from Lalpha absorption measurements. II._ , _ApJ_ 224 (1978) 132.
  * [14] Gaia Collaboration, T. Prusti, J.H.J. de Bruijne, A.G.A. Brown, A. Vallenari, C. Babusiaux et al., _The Gaia mission_ , _Astron. Astrophys._ 595 (2016) A1 [1609.04153].
  * [15] Gaia Collaboration, A. Vallenari, A.G.A. Brown, T. Prusti, J.H.J. de Bruijne, F. Arenou et al., _Gaia Data Release 3: Summary of the content and survey properties_ , _arXiv e-prints_ (2022) arXiv:2208.00211 [2208.00211].
  * [16] A. Ulla, O.L. Creevey, M.A. Álvarez, R. Andrae, C.A.L. Bailer-Jones, I. Bellas-Velidis et al., “Gaia DR3 documentation Chapter 11: Astrophysical parameters.” Gaia DR3 documentation, European Space Agency; Gaia Data Processing and Analysis Consortium., June, 2022.
  * [17] E. Kafexhiu, F. Aharonian, A.M. Taylor and G.S. Vila, _Parametrization of gamma-ray production cross sections for p p interactions in a broad proton energy range from the kinematic threshold to PeV energies_ , _Physical Review D_ 90 (2014) 123014 [1406.7369].
  * [18] M. Aguilar, D. Aisa, B. Alpat, A. Alvino, G. Ambrosi, K. Andeen et al., _Precision Measurement of the Proton Flux in Primary Cosmic Rays from Rigidity 1 GV to 1.8 TV with the Alpha Magnetic Spectrometer on the International Space Station_ , _Physical Review Letters_ 114 (2015) 171103.
  * [19] W. Liu, Y.-Q. Guo and Q. Yuan, _Indication of nearby source signatures of cosmic rays from energy spectra and anisotropies_ , _J. Cosmol. Astropart. Phys._ 2019 (2019) 010 [1812.09673].
  * [20] B. Zhao, W. Liu, Q. Yuan, H.-B. Hu, X.-J. Bi, H.-R. Wu et al., _Geminga SNR: Possible Candidate of the Local Cosmic-Ray Factory_ , _ApJ_ 926 (2022) 41 [2104.07321].