\section{Introduction}

Data analysis tools are essential for data science. Robust movement data analysis tools are therefore key to advancing mobility data science. However, the development of movement data analysis tools is hampered by a lack of shared understanding and standardization. There are numerous implementations, including dozens of R libraries, as well as Python libraries and moving object databases. For example, \cite{joo2020navigating} review 57 R libraries related to movement in ecology. The development of Python libraries for movement analysis is picking up as well. For example, \cite{graser_movingpandas_2019} introduces the Python library MovingPandas and compares it to the R trajectories library \cite{moradi2018trajectories} and the movement database Hermes \cite{pelekis2015hermes}. Other Python libraries for movement analysis include scikit-mobility (focusing on human movement) \cite{pappalardo2019scikitmobility} and Traja (for animal movement) \cite{justin_shenk_2019_3237827}. However, even though they are all built for movement data analysis, they still vary considerably in their underlying concepts and provided functionality.

The ISO Standard 19141 Geographic information -- Schema for moving features was first published in 2008 \cite{ISO19141}. It defines a data model for representing moving points and moving rigid regions. Based on this standard, the Open Geospatial Consortium (OGC) has worked on various encodings for representing moving features as well as standardized operations for manipulating moving geometries. However, so far, these standards have failed to reach significant adoption.

This discussion paper summarizes the current status as well as open issues towards a universal API for mobility data science. We start with a short introduction to OGC Moving Features standard and provide conceptual and technical commentary on these standardization efforts. These findings should serve as a point of departure for a community-wide effort to develop a shared understanding of what would be the right API for mobility data science -– if such a thing exists. Based on these findings, we then propose essential building blocks to advance movement data science.

\begin{figure}[tb]
\includegraphics[width=\columnwidth]{example}
\caption{Data model of the Moving Features standard illustrated with two moving points A and B. Stars mark changes in attribute values.}
\label{fig:example}
\end{figure}

\section{Moving Features}

The ISO and OGC Simple Features standard \cite{ISO19125} has gained widespread adoption in Geographic Information Systems (GIS) and related systems dealing with spatial data. Following this example, the new set of OGC Moving Features standards (which are based on the ISO standard \cite{ISO19141}) define how moving features should be encoded and how they should be accessed. A moving feature contains a temporal geometry, whose location changes over time, as well as dynamic non-spatial attributes whose values vary with time. The standards supports 0-dimensional (points), 1-dimensional (lines), 2-dimensional (polygons), and 3-dimensional (polyhedrons) geometries that vary over time. The standard allows the representation of the following phenomena:

\begin{enumerate}
\item Discrete phenomena, which exist only on a set of instants, such as road accidents,
\item Step phenomena, where the changes of locations are abrupt at an instant, such as the location of mobile speed cameras for monitoring traffic, and
\item Continuous phenomena, whose locations move continuously for a period in time, such as vehicles, typhoons, or floods.
\end{enumerate}

\begin{table*}[h]
  \centering
  \begin{tabularx}{\textwidth}{lXXX}
    \hline
    & \multicolumn{1}{c}{\textbf{CSV}}
    & \multicolumn{1}{c}{\textbf{XML}}
    & \multicolumn{1}{c}{\textbf{JSON}} \\
    \hline
    \textbf{Concept}
      & Segments with start and end time and attributes, no static object properties
      & Segments with start and end time and attributes, with static properties
      & Points with timestamps and attributes with (independent) timestamps
      \\
    \hline
    \textbf{Advantages}
      & \begin{enumerate}
            \item Supports temporal gaps in observations
        \end{enumerate}
      & \begin{enumerate}
            \item Supports temporal gaps in observations
            \item Can handle complex geometries
            \item Support for static / non-temporal descriptive object properties
        \end{enumerate}
      & \begin{enumerate}
            \item Most compact representation
            \item Can handle complex geometries
            \item Time stamps of location and attribute changes are modelled independently -- no synchronization necessary
            \item Interpolation modes can be specified individually for each attribute
            \item Support for non-temporal descriptive attributes
        \end{enumerate}
        \\
    \hline
    \textbf{Limitations}
      & \begin{enumerate}
            \item Verbose: redundant information (start and end time and location for each segment)
            \item Temporal geometry and temporal attributes have to be synced
            \item Only linear interpolation
            \item Only moving points
            \item Not readily usable in GIS (missed the opportunity to use WKT to represent the geometry)
        \end{enumerate}
      & \begin{enumerate}
            \item Very verbose: XML \& redundant information (start and end time and location for each segment)
            \item Temporal geometry and temporal attributes have to be synced
            \item Same interpolation for geometry and all attributes
        \end{enumerate}
      & \begin{enumerate}
            \item Multiple options for encoding the same situation (e.g. unclear bounds of time periods)
            \item No support for temporal gaps in observations
        \end{enumerate}
        \\
    \hline
  \end{tabularx}
  \caption{Overview of OGC Moving Feature encodings}
  \label{tab:1}
\end{table*}

Fig.~\ref{fig:example} illustrates the Moving Features data model. It shows the movement of two vehicles denoted by $A$ and $B$. The goal is to model the vehicles' changing location and gear settings. The vehicle locations are modelled as moving points while the selected gear is modelled as a temporal property. Both vehicles start their movement at time $t=0$. While vehicle $A$ records its location at time $t=1$ and $t=2$, vehicle $B$ records its location only at time $t=2$. Furthermore, timestamps when $A$ changes gears are marked by a star symbol.

This example follows the default assumption of linear interpolation between locations.
In general, a moving feature can be modelled as positions recorded at discrete timestamps, where the position between two recorded timestamps is computed by interpolation. If the interpolation is discrete, then no assumption is made about the location of the object between two observations. On the other hand, stepwise interpolation assumes that the location of the object is constant between two observations. Finally, linear or polynomial interpolations specify how to compute the location of the object between two observations.

Based on this data model, OGC Simple Features defines XML/GML \cite{OGC14-083r2} and JSON \cite{OGC16-140r1} data encodings, as well as simpler CSV \cite{OGC14-084r2} and binary encodings (which are limited to moving points). However, even though these encodings are part of the same standard and are based on the same general data model, there are still considerable differences that significantly affect which kind of information can and can not be modelled (as shown in Table~\ref{tab:1}). To dive deeper and illustrate our point, the following sections show the example introduced in  Fig.~\ref{fig:example} encoded using the CSV, XML, and JSON encodings.

\subsection{CSV Encoding}

The CSV encoding is the simplest option in the Moving Features standard.
Fig.~\ref{fig:csv} shows the CSV representation of the example introduced in Fig.~\ref{fig:example}. The CSV structure is divided into two parts: First, the header lines (starting with the ``\verb+@+" character) provide the meta-information. Then the trajectory lines describe the movement. The CSV encoding is segment based \cite{OGC14-084r2}, that is, each line describes a trajectory segment with two or more coordinate pairs (or triplets for 3D trajectories).

\begin{figure}[htb]
\begin{lstlisting}
@stboundedby,urn:x-ogc:def:crs:EPSG:6.6:4326,
	2D,10.0 10.0,10.6 12.2,2011-07-14T22:00:00Z,
	2011-07-14T22:00:20Z,sec
@columns,mfidref,trajectory,gear,xsd:integer
A,0,5,10.0 10.0 10.2 10.6,1
A,5,10,10.2 10.6 10.4 11.2,2
A,10,15,10.4 11.2 10.5 11.7,2
A,15,20,10.5 11.7 10.6 12.2,3
B,0,20,2.0 2.0 2.1 2.1,1
\end{lstlisting}
\caption{CSV encoding of vehicles $A$ and $B$}
\label{fig:csv}
\end{figure}

\begin{figure*}[htb]
\begin{lstlisting}[language=XML]
<?xml version="1.0" encoding="UTF-8"?>
	<mf:MovingFeatures xmlns:mf="http://schemas.opengis.net/mf-core/1.0" ...>
	<mf:STBoundedBy offset="sec">
		<gml:EnvelopeWithTimePeriod srsName="urn:x-ogc:def:crs:EPSG:6.6:4326">
			<gml:lowerCorner>10.0 10.0</gml:lowerCorner>
			<gml:upperCorner>10.6, 12.2</gml:upperCorner>
			<gml:beginPosition>2011-07-14T22:00:00Z</gml:beginPosition>
			<gml:endPosition>2011-07-14T22:00:20Z</gml:endPosition>
		</gml:EnvelopeWithTimePeriod>
	</mf:STBoundedBy>
	<mf:Member>
		<mf:MovingFeature gml:id="A">
			<gml:name>NissanA</gml:name>
			<gml:description>Nissan Sentra ...</gml:description>
		</mf:MovingFeature>
	</mf:Member>
	<mf:Header>
		<mf:VaryingAttrDefs>
			<mf:attrDef name="gear" type="xsd:integer">
			<mf:AttrAnnotation>The gear number used... </mf:AttrAnnotation>
			</mf:attrDef>
		</mf:VaryingAttrDefs>
	</mf:Header>
	<mf:Foliation>
		<mf:LinearTrajectory gml:id="LT0001" mfIdRef="A" start="0" end="5">
			<gml:posList>10.0 10.0 10.2 10.6</gml:posList>
			<mf:Attr>1</mf:Attr>
		</mf:LinearTrajectory>
		<mf:LinearTrajectory gml:id="LT0003" mfIdRef="A" start="5" end="10">
			<gml:posList>10.2 10.6 10.4 11.2</gml:posList>
			<mf:Attr>2</mf:Attr>
		</mf:LinearTrajectory>
		<mf:LinearTrajectory gml:id="LT0003" mfIdRef="A" start="10" end="15">
			<gml:posList>10.4 11.2 10.5 11.7</gml:posList>
			<mf:Attr>2</mf:Attr>
		</mf:LinearTrajectory>
		<mf:LinearTrajectory gml:id="LT0003" mfIdRef="A" start="15" end="20">
			<gml:posList>10.5 11.7 10.6 12.2</gml:posList>
			<mf:Attr>3</mf:Attr>
		</mf:LinearTrajectory>
	</mf:Foliation>
</mf:MovingFeatures>
\end{lstlisting}
\caption{XML encoding of vehicle $A$}
\label{fig:xml}
\end{figure*}

The first header line specifies the spatio-temporal boundary of the moving features. It lists the spatial reference system used, the number of dimensions of the geometries, the coordinates of the corners of the spatial boundary, the start and end times, and the time encoding, that is, the units of time used in the trajectory lines for encoding the offset from the start time.

The second header line specifies the columns in the trajectory lines. The first column \verb+mfidref+ is used to identify the moving object. The second column \verb+trajectory+ defines the spatio-temporal geometry, followed by the definition of the time-varying attributes (name and its type).

The following trajectory lines contain the spatio-temporal geometries for moving features. Each line specifies the \verb+mfidref+, the start and end time, which can be specified as absolute or relative (offset) values, the points of the line string and the attribute values.
Since temporal geometry and temporal attributes (in this case the gear) are specified together, it is necessary to encode intermediate locations each time a temporal attribute changes its value.

The CSV encoding assumes linear interpolation between locations. There is no way to store timestamps for intermediate positions along trajectory segments. The interpolation assumes constant speed along a segment.

A notable design choice is that that the segment geometries are not encoded as well-known text (WKT), which would have made the CSV more easily readable by existing GIS tools.

\subsection{XML Encoding}

Like the CSV encoding, the XML encoding is also segment based.
Fig.~\ref{fig:xml} shows the XML representation of the example introduced in Fig.~\ref{fig:example}. It is obvious that this encoding is considerably more verbose than the CSV encoding. The most essential XML tags are:

\begin{itemize}
\item \verb+mf:STBoundedBy+ specifies the spatiotemporal bounding box and time units
\item \verb+mf:Member+  specifies properties of the moving features
\item \verb+mf:Header+  contains meta-information about the moving features, for example, attribute definition
\item \verb+mf:Foliation+ specifies the moving geometry and the dynamic attributes
\end{itemize}

While the XML encoding does support different interpolation types, the same interpolation applies to all attributes (linear by default). Therefore, we cannot specify that the geometry evolves linearly while the gear changes in a stepwise manner.

\begin{figure*}[htb]
\begin{lstlisting}
{
  "type": "MovingFeature",
  "temporalGeometry": {
      "type": "MovingPoint",
      "coordinates": [ [10.0, 10.0], [10.4, 11.2], [10.6, 12.2] ],
      "datetimes": ["2011-07-14T22:00:00Z", "2011-07-14T22:00:10Z", "2011-07-14T22:00:20Z"],
      "interpolations": "Linear"
  },
  "temporalProperties":  [ {
      "name": "gear",
      "values": [1, 2, 3, 3],
      "datetimes": ["2011-07-14T22:00:00Z", "2011-07-14T22:00:05Z", "2011-07-14T22:00:15Z",
                         "2011-07-14T22:00:20Z"],
      "interpolations": "Stepwise"
  }, ],
  "stBoundedBy": {
      "bbox": [10.0, 10.0, 10.6, 12.2],
      "period": { "begin": "2011-07-14T22:00:00Z", "end" : "2011-07-15T22:00:20Z" }
  },
  "properties": {
      "name": "NissanA", "description": "Nissan Sentra ..."
  }
}
\end{lstlisting}
\caption{JSON encoding of vehicle $A$ }
\label{fig:json}
\end{figure*}

\subsection{JSON Encoding}

Contrary to the previous two encodings, the JSON encoding is point based rather than segment based. Fig.~\ref{fig:json} shows the JSON representation of the example introduced in Fig.~\ref{fig:example}. The essential elements are:

\begin{itemize}
\item \verb+temporalGeometry+ represents the coordinates and timestamps in parallel arrays and states that linear interpolation is used
\item \verb+temporalProperties+ allows us to represent several such properties, where the values and the timestamps of each of them are represented as parallel arrays, and the interpolation is also stated, which is stepwise in our example for the attribute gear
\item \verb+stBoundedBy+ states the spatio-temporal bounding box
\item \verb+properties+ represents static / non-temporal descriptive attributes
\end{itemize}

The JSON encoding is considerably more compact than the CSV and XML encodings since the temporal geometry and the temporal attributes and their respective timestamps are encoded independently. In a real-world vehicle tracking use case, for example, there would considerably fewer gear changes than location changes. Therefore, it saves space to only record gear change events when necessary instead of reporting the gear setting at every location.

However, there are still ambiguous situations: for example, for attributes with stepwise interpolation, it is not possible to specify whether the change of value is exactly at the timestamp or just after it. Therefore, when trying to determine the period of time when the first gear is used, we cannot specify whether it was used during [2011-07-14 22:00:00, 2011-07-14 22:00:05] or [2011-07-14 22:00:00, 2011-07-14 22:00:05) (notice the right-open end of the period).
This has implications when computing topological predicates for moving geometries, because -- depending on the predicate considered -- the fact that a moving geometry is inside or crosses a boundary needs to take into account whether the bounds of the period are left- and right-inclusive or not.

\subsection{Moving Feature Encodings Summary}

As the examples in the previous sections show, the three Moving Feature encoding schemes implement different concepts for modeling movement data. This makes it hard to understand and use the standard since each encoding choice has different consequences.

To summarize the situation, Table~\ref{tab:1} lists the advantages and limitations of the three Moving Feature encoding options.
Overall, the JSON encoding (which is the most recent addition to the standard) is the most compact and has the fewest limitations.
However, it is worth mentioning that the segment based CSV and XML encodings (unlike the point based JSON encoding) allow us to represent temporal gaps due to the fact that every segment has both a start and an end timestamps. This is important from a modeling perspective since it makes it possible to represent, for example, that the GPS signal was lost while a moving car was inside a tunnel.

\subsection{Moving Features Access}

Besides encodings, Moving Features also standardizes functions for working with movement data. The corresponding standard is called OGC Moving Features Access \cite{OGC16-120r3} and it covers a wide range of functions for the retrieval of trajectory feature attributes, as well as operations between trajectories and geometry objects. Functions to retrieve feature attributes include basic functions to access the location, speed, or acceleration at a given time, or the subtrajectory between two timestamps. Conversely, there is also a function to extract the time at a given point. These and similar functions can be found in many existing tools, such as the ones analyzed by \cite{graser_movingpandas_2019}.

Beyond these basic functions, however, the standard also defines more obscure functions, such as, for example, timeToDistance which “shall return a graph of the time to distance function as a set of curves in the Euclidean space consisting of coordinate pairs of time and distance”. The standard remains unclear as to how exactly this curved graph should be implemented.

\subsection{Discussion}

The limited success of OGC Moving Features so far may be due to its failure to involve a wide user base during its development stage. The lack of wide-spread engagement certainly limits public awareness of the standard's existence. Some of the limitations (particularly of the CSV and XML encodings) as well as ambiguities in the standard may be due to limited diversity of perspectives represented during the development of the standard.

The decision to define multiple encoding standards that implement different concepts for modelling movement data further complicates the situation and potentially hinders the standard's adoption.

Another limiting factor is that, currently, there is no official reference implementation of OGC Moving Features. However, MobilityDB \cite{zimanyi2019mobilitydb, MobilityDB} builds upon the OGC Moving Features specification to add support for spatiotemporal objects to PostGIS databases. This effort may be considered a reality check to gauge the real-world applicability of the standard.

\section{Towards an API for Mobility Data Science}

Multiple notable recent publications, such as \cite{dodge_analysis_2016, demsar_analysis_2015}, have emphasized the need for a movement data science that bridges the boundaries of established domains, such as ecology, human mobility analytics for health and planning, or transportation research. Why then has there been such limited progress towards a common understanding and implementations of the key elements of an API for movement data analysis or even just a common data exchange format? The following conjectures summarize some of the main challenges that we have identified:

\textbf{Conjecture 1: The curse of movement data heterogeneity.} From high-resolution continuous tracking data of cooperative moving objects to sparse checkpoint-based trajectories with large spatiotemporal uncertainties, the flavours of movement datasets vary vastly. Furthermore, there is no consistent terminology \cite{graser_movingpandas_2019}: for example, terms such as trajectory, track, path, move, trip or travel may used to describe the same or different concepts related to movement.
Moreover, while some domains still struggle to collect sufficient data, others are amassing huge amounts thanks to recent advances in location and communication technologies. But even if data is plentiful, there is a lack of commonly available data management tools for such volumes of data, and a pressing need for them. However, this heterogeneity means that it is unlikely that there will be one solution that fits all requirements. Therefore, our chances to build a flourishing mobility data science environment would be vastly improved by a shared terminology, as well as a set of standard data formats and well-defined analysis functions.

\textbf{Conjecture 2: Lack of pragmatic solutions.} Many recently popularized data formats followed a somewhat opposite approach to established standardization processes. For example, GeoJSON, Vector tiles, or GTFS (General Transit Feed Specification) and the software tools around them quickly gained popularity as useful and pragmatic solutions for common problems. Individual mobility data projects, such as, for example, \cite{TraceJSON} try to follow a similar approach but have yet to reach a critical mass of users. The challenge lies in addressing the common problems faced by movement data analysts in different domains.

\textbf{Conjecture 3: Wrong priorities and incentives.} Spatiotemporal data analysis is a hard nut to crack. However, there is no shortage of publications and sophisticated concepts, reaching back to, works such as \cite{sistla1997modeling}. Similarly, there are numerous research software prototypes but most of them are not usable since they were either never made publicly available or not engineered for use by others other than the original author. The current scientific environment provides very limited incentives for the development of well-engineered scientific software. It is risky to spend time on development when academic hiring committees may only look at the conventional publication record.

\section{Conclusion}

To solve the above mentioned issues and advance the mobility data science domain, we have to take steps to overcome these challenges. Looking back on how static/non-moving spatial data science tools have evolved over the last decades, we advocate for the following essential building blocks:

\textbf{A. Common movement data analysis concepts:} Whether through official standards or by adopting a pragmatic de facto standard, a common terminology would vastly improve the exchange of movement data methodology. However, so far, there are no strong contenders for a widely adopted standard. OGC's Moving Feature Standards lacks essential features that would facilitate its wide-spread adoption. A well-accepted standard would provide a common reference framework for development efforts from all sides of the mobility data science community. Furthermore, it has the potential to facilitate the exchange of movement analysis methods by establishing a shared language with a well defined terminology.

\textbf{B. An open general-purpose mobility (data) engine:} There are standard libraries for manipulating non-moving geometries (such as the JTS Topology Suite (Java Topology Suite) and its C++ port GEOS \cite{geos}) and numerous spatial data analysis tools (such as GeoTools, GDAL, R, PostGIS, QGIS, Google Earth, and NASA WorldWind) build on these libraries. Similarly, mobility data science would profit from a standard library for manipulating moving geometries. To avoid unnecessary duplication of spatial data handling functionality, such a mobility engine should build upon existing libraries for non-moving geometries. The mobility engine development should focus on the specific spatiotemporal aspects of movement data, while delegating the geometric manipulation to the underlying libraries. Language bindings for higher-level programming languages such as Python and R should be considered essential for a wider industry adoption. This general purpose mobility engine would enable faster development of more specialized movement data analysis tools while still keeping a common base framework.

\textbf{C. Focus on open science and reproducibility:} A push towards open scientific practices, including code sharing and replication of results by peers (before and after publication) would be one step towards incentivizing the development of better scientific software. Furthermore, it would speed up research and development by removing the need to re-implement the same things over and over again and instead build on each others work.

\bibliographystyle{apalike}


\section{Introduction}

\iffalse
In the last few years, a lot of works try to learn how people make decisions. However, those models usually simulate the process within latent space which is a black-box for humans. People can only passively accept the result without knowing the reasons. Among those decision-making algorithms, recommender system algorithms are broadly used in normal life. Moreover, researchers gradually realized the benefits brought by explanations in recommender systems. Through explaining why the system recommend such items, users are potentially telling why the algorithm recommends wrong items. The transparent system will increase the confidence or the trust of customers which can gain more profit in the long term. Not only academics are studying this issue, but the explainable recommendation systems have emerged in our daily lives. For example, the online shopping platform Amazon already started to tell consumers why recommended (Fig.\ref{fig:amazon}). \zfcomment{what behavior to humans?}

After the term \textit{explainable recommendation} was formally introduced in recent years\cite{zhang2014explicit}, a large number of
great works focus on this area. Many popular explainable recommendation approach is model-based method that one designed model can provide both of recommendation and explanation. The model-based explainable recommendation algorithms are closely related to recommendation model with machine learning such as matrix factorization, topic modelling, etc.

However, there are many much more complicated models cannot be explained. For example, in an e-commerce system, the recommendation lists are usually generated by a hybrid system containing many different methods. In this case, we rely on model-agnostic methods (also known as post hoc) which separate recommendation and explanation into two independent models. After an item is recommended by any recommendation model, the explanation model can generate a simple explanation based on statistic or data mining, for example, "because of another item you purchased". The main focus for these approaches is to find out the relationships between items.

Peake and Wang proposed an association rule mining approach\cite{peake2018explanation} for post hoc explainable recommendation. The authors consider a matrix factorization model as a black-box model and extract the association rules based on the data consisting of the user history and recommended items. The association rules are used to generate the explanation if the recommended item can also be found in an association rule, for example, an item pair $(A, B)$ can be explained as "Because you purchased A, the system recommends B".

An association rule is consist of an antecedent and a consequent which could be a kind of explanation. However, it is just co-occurrence instead of causality. Explanation by association rules could be a coincidence which cannot be a reliable explanation. One example is the "beer and diapers" story that young American males who buy diapers also have a predisposition to buy beer. The explanation with association rules may get that "because you purchased diapers, the system recommends beer". Obviously, it is not a trustworthy explanation. Instead, when explaining the recommendation with causality, this kind of mistakes due to high co-occurrence will disappear.

In this paper, we apply a method extracting causal relationships for post hoc explainable sequential recommendation. We try different sequential recommendation models as the black-box model and give explanations. We will theoretically prove that the relationship that extracted by our model is high probably a causality. Meanwhile, our model will get better performance against baselines in term of model fidelity which is defined as the percentage of explainable items in recommended items.
\fi
As widely used in decision-making, recommender systems have been recognized for its ability to provide high-quality services that reduces the gap between products and customers.
Nowadays, many state-of-the-art performances are achieved by neural network models.
Typically, deep learning models are used as black-box latent factor models accompanied with a high-dimensional latent space.
This allows them to achieve good expressiveness power and accuracy in various recommendation tasks.
However, it is also true that complex neural models
easily go beyond the comprehension of the majority of customers, since thousands or millions of parameters are involved.
Nevertheless, it is a natural demand of human-beings to understand why a model makes a specific decision, rather than blindly accepting the results without knowing the underlying reason. As a result, providing supportive information and interpretation along with the recommendation can be helpful for both the customers and the platform,
since it improves the transparency, trustworthiness, effectiveness, and user satisfaction of the recommendation systems, while facilitating system designers to refine the algorithms \cite{zhang2018explainable}. For example, a user may be inspired by the recommendation panel explained as ``you may also like'' in an e-commerce system, and thus decides to look around for more items that better satisfy his or her interests. In the meantime, the user preference would be captured more precisely, so that better services can be provided in further interactions.
On the other hand, system designers may also easily figure out the reasons of providing unsatisfied recommendation,
leveraging the result to take actions.

\begin{figure}[t]

    \centering
    \includegraphics[scale=0.5]{./fig/post-hoc.pdf}
    \vspace{-15pt}
    \caption{The model-intrinsic approaches directly generate the recommendations and explanations, while the model-agnostic approaches generate explanations after the black-box model has generated the recommendations.}
    \vspace{-20pt}
    \label{fig:posthoc}
\end{figure}

To address the explainability problem, researchers turned to explainable recommendation models, which are expected to not only generate effective recommendations but also intuitive explanations to humans.
Generally, the explainable models can be either model-intrinsic or model-agnostic (also known as post-hoc), as shown in Fig.\ref{fig:posthoc}.
The model-intrinsic approach takes advantage of the interpretable mechanism of the model, with the explanations directly provided as the intermediate stages of recommended decisions.
For instance, simple user-based collaborative filtering methods
pass an item from one user to another similar user, and thus it can directly output explanations like ``similar users also bought this item''.
However, the recommendation explanations are usually not obtained for free, and sometimes we have to trade-off with model accuracy~\cite{theodorou2017designing}. In many cases, outstanding recommendation performances are usually achieved by models that are less interpretable.
Thus, it is very challenging for current explainable recommendation methods to redesign a black-box model into an interpretable one whilst maintaining the recommendation performance \cite{zhang2018explainable}.

In contrast, post-hoc models make no assumption of the underlying recommendation model, and allow the decision mechanism to be a black-box, since it will provide explanations after a decision is made.
Although such explanations may not strictly follow the exact mechanism that generated the recommendations, they offer the flexibility to be applied to many different models.
Though it is still not fully understood what information is useful for generating explanations for a certain recommendation result, Peake ~\cite{peake2018explanation} argued that one can provide post-hoc item-level explanations.
In other words, interacted items (the causes) in a user's history can be used as explanations for the future item recommendations (the effect), since it answers ``what items you bought causes the recommendation of this item''.

However, existing work mostly use global association rule mining to discover the relationship between items, which relies on the item co-occurrence among all user transactions.
Therefore, the explanations are not personalized, i.e., different users would receive the same explanation as long as they are recommended with the same item and have overlapped histories.
This makes it incompatible with modern recommender systems, which aims to provide personalized services to users.
Moreover, the item-level explanation problem naturally involves causal analysis between a user's previous and future behaviors, which makes the problem even more challenging since it has to answer counterfactual questions such as ``what would happen if a different set of items were purchased''.

In this paper, we explore a causal analysis framework to provide post-hoc causal explanations for any sequential recommendation algorithm.
The goal is to design a model that generates post-hoc explanations for black-box recommendation models in order to reduce the accuracy-interpretability trade-off. Since the explanation model and recommendation model work separately, we obtain the benefit of explainability without hurting the prediction performance. Technically, we propose a Variational Auto-Encoder (VAE) based perturbation framework to create counterfactual examples for causal analysis, which extracts causal rules between a user's previous and future behaviors as explanations.

The key contributions of this paper are as follows:

\begin{itemize}
    \item We design and study a causal rule mining framework for sequential recommendation.
    \item We show that this framework can generate personalized post-hoc explanations based on item-level causal rules to explain the behaviors of a sequential recommendation model.
    \item We conduct experiments on real-world data to show that our explanation model outperforms state-of-the-art baselines.
\end{itemize}

In the following, we review related work in in Section 2, and introduce our model in Section 3. Experimental settings and results are provided in Section 4. We conclude this work in Section 5.

\iffalse
To date, researchers mainly focus on achieving the highest possible model performance and have been making great efforts to solve optimization problems in machine learning area. Most machine learning algorithms use black-box latent factor models to abstract their learning processes to a high-dimensional latent space, which goes beyond the direct comprehension of humans. Thus they can not provide an explanation of how a prediction has been made. With the broad use of machine learning algorithms in decision-making tasks, especially in some risk sensitive areas, such as health care, legal and fintech fields, the problems caused by the intransparency have emphasized the necessity of ensuring the explicability of models to develop people's trust in their decisions. In many cases, humans need to understand why a model make a specific decision rather than just can only decide whether to accept the results without knowing the reason. For example, in "human-in-the-loop" pattern where requires human interaction, the decisions given by models must be understandable by humans for achieving cooperation with machines. Therefore, it is very important yet challenging to establish approaches which can provide explanation for the behavior of the black-box models of decision-making algorithms without sacrificing their performance.

As a widely used decision-making system, state-of-the-art recommender systems, though having the ability to provide high-quality recommendations, usually can't explain their working process and results to humans. This lack of explainability for recommendation algorithms leads to critical problems in practice. For example, without letting users know the reason of providing specific results, the system may be less effective in persuading users to adopt the recommendations and further decrease the system's trustworthiness. What's more, providing supportive information and interpretation along with the recommendation is also potentially helpful for complicated decision making. For example, doctors may need comprehensive healthcare recommendations with interpretations to make best diagnoses.

Explainable recommendation algorithms focus on developing models which can generate not only effective recommendations but also intuitive explanations to help humans, including the users and system designers, to understand why a particular result is provided, and how to leverage the result to take actions, thus help to improve the transparency, trustworthiness, effectiveness and user satisfaction of the recommendation systems, as well as facilitates system designers to debug and refine the algorithms.

Typically, the explainable models can be either model-intrinsic or model-agnostic (also known as post-hoc) as shown in Fig.\ref{fig:posthoc}. The decision mechanism for model-intrinsic approach is transparent since they develop interpretable models, so as to provide explanations for the model decisions directly. Differently, the post-hoc models allow the decision mechanism to be a black-box and provide explanations for a decision after it has been made. Although such explanations may not follow the exact mechanism strictly, they have advantages that can be applied in many different models more flexibly.

In this paper, we design a model to generate post-hoc explanation for black-box recommendation models in order to reduce the accuracy-interpretability trade-off. In other words, since the explanation model and black-box model are working separately and the black-box model is fixed after trained, we can obtain the benefit of explainability without hurting the predictive performance.

The paper \cite{peake2018explanation} provides an approach to give the post-hoc interpretation for recommendation systems by mining the association rules among the purchased items. However, such explanation can't be satisfactory and even can be incomprehensible sometimes, since the association rules can only tell the co-occurrence of items, but do not represent the causation. What's more, too many spurious association rules among the data, such as the fake association brought by the sample selection bias, may usually cause the unreliable results of the systems. One classical example of association rule for the human purchase behavior is $\{$diapers$\rightarrow$ beer$\}$, which indicates that if a customer purchases diapers, it's likely that he also purchases beer. The explanation given by association rule is that "Because you purchased diapers, the system recommends you beer", which is obvious a confusing explanation. Thus it's naturally to think about using a kind of more stable and explainable relationship among the data to remove the fake associations and unreliable explanations caused by high co-occurrence rat
, which is the causation.

To the best of our knowledge, it is the first work to provide causal explanations for sequential recommendation systems. The causal relationships can be extracted from itemsets of user purchase history through our proposed two-step model. Firstly, for doing causal analysis, a perturbation model is introduced to generate data for achieving random experiments on black-box models. After that, causal relationships can be captured within a causal rule mining module. Our model is tested on real dataset to provide causal explanations for a state-of-the-art sequential recommendation system. The extracted causal rules are evaluated by showing their reliability, model fidelity and superior explanation performance comparing with association rules, which exhibits that causal rules can provide more reliable and intuitive explanations for black-box sequential recommender systems.

In the following, we review some related work in Section 2, introduce our model in Section 3, and conduct experiments and make discussions in Section 4. Finally, we conclude this work in Section 5.
\fi
\section{Related Work}
\subsection{Sequential Recommendation}
Sequential recommendation takes into account the historical order of items interacted by a user, and aims to capture useful sequential patterns for making consecutive predictions of the user's future behaviors.
Rendle et al. \cite{rendle2010factorizing} proposed Factorized Personalized Markov Chains (FPMC) to combine Markov chain and matrix factorization for next basket recommendation.
The Hierarchical Representation Model (HRM) \cite{wang2015learning} further extended this idea by leveraging representation learning as latent factors in a hierarchical model.
However, these methods can only model the local sequential patterns of very limited number of adjacent records.
To model multi-step sequential behaviors, He et al. \cite{he2016fusing} adopted Markov chain to provide recommendations with sparse sequences.
Later on, the rapid development of representation learning and neural networks introduced many new techniques that further push the research of sequential recommendation to a new level.
For example, Hidasi et. al. \cite{hidasi2015session} used an RNN-based model to learn the user history representation,
Yu et. al. \cite{yu2016dynamic} provided a dynamic recurrent model,
Li et. al. \cite{li2017neural} proposed an attention-based GRU model, Chen et. al. \cite{chen2018sequential} developed user- and item-level memory networks, and Huang et. al. \cite{huang2018improving} further integrated knowledge graphs into memory networks.
However, most of the models exhibit complicated neural network architectures, and it is usually difficult to interpret their prediction results.
As a result, we would like to generate explanations for these black box sequential recommendation models.

\subsection{Explainable Recommendation}
Explainable recommendation focuses on developing models that can generate not only high-quality recommendations but also intuitive explanations, which help to improve the transparency of the recommendation systems \cite{zhang2018explainable}.
Generally, the explainable models can either be model-intrinsic or model-agnostic as introduced in introduction.
As for model-intrinsic approaches, lots of popular explainable recommendation methods, such as factorization models \cite{zhang2014explicit,chen2018attention,wang2018explainable},
deep learning models \cite{seo2017interpretable,chen2018sequential,li2019capsule,costa2018automatic},
and knowledge graph models \cite{ai2018learning,xian2019reinforcement,huang2018improving,fu2020fairness,ma2019jointly,wang2019explainable}
have been proposed. A more complete review of the related models can be seen in \cite{zhang2018explainable}.
However, they mix the recommendation mechanism with interpretable components, which often results in a system too complex to make successful explanations.
Moreover, the increased model complexity may reduce the interpretability.
A natural way to avoid this dilemma is to rely on model-agnostic post-hoc approaches so that the recommendation system is free from the noises of the down-stream explanation generator.
Examples include \cite{mcinerney2018explore} that proposed a bandit approach, \cite{wang2018reinforcement} that proposed a reinforcement learning framework to generate sentence explanations, and \cite{peake2018explanation} that developed an association rule mining approach.
In their work, the transactions of all users, which consider user history as input and recommendation item as output, are used to extract association rules as the explanation for black-box models.
However, correlation is not reliable for its direction agnostic feature.
As will mention in the next part, our goal here is to find causal relationships in the user behaviour, which we can provide more stable explanations.

\subsection{Causal Inference in Recommendation}
Originated as statistical problems, causal inference \cite{pearl2000causality,imbens2015causal}
aims at understanding and explaining the causal effect of one variable on another.
While the observational data is considered as the factual world, causal effect inferences should be aware of the counterfactual world, thus often regarded as the question of "what-if".
The challenge is that it is often expensive or even impossible to obtain counterfactual data.
For example,
it is immoral to re-do the experiment on a patient to find out what will happen if we have not given the medicine.
Though the majority of causal inference study resides in the direction of statistics and philosophy, it recently attract attention from AI community for its great power of explainablity and bias elimination ability.
Efforts have managed to bring causal inference to several machine learning areas, including recommendation \cite{bonner2018causal},
learning to rank \cite{joachims2017unbiased},
natural language processing \cite{wood2018challenges},
and reinforcement learning  \cite{buesing2018woulda},
etc.
With respect to recommendation tasks, large amount of work is about how to achieve de-bias matrix factorization with causal inference.
The probabilistic approach ExpoMF proposed in \cite{liang2016modeling} directly incorporated user exposure to items into collaborative filtering, where the exposure is modeled as a latent variable.
Liang et. al.\cite{liang2016causal} followed to develop a causal inference approach to recommender systems which believed that the exposure and click data came from different models, thus using the click data alone to infer the user preferences would be biased by the exposure data.
They used causal inference to correct for this bias for improving generalization of recommendation systems to new data.
Bonner et. al.\cite{bonner2018causal} proposed a new domain adaptation algorithm which was learned from logged data including outcomes from a biased recommendation policy, and predicted recommendation results according to random exposure.
Differently, this paper focuses on learning causal rules to provide more intuitive explanation for the black-box recommendation models.
Additionally, we consider \cite{alvarez2017causal} as a highly related work though it is originally proposed for natural language processing tasks.
As we will discuss in the later sections, we utilize some of the key ideas of its model construction, and show why it works in sequential recommendation scenarios.

\vspace{-3pt}
\section{Proposed Approach}
In this section, we first define the explanation problem and then introduce our model as a combination of two parts: a VAE-based perturbation model that generates the counterfactual samples for causal analysis, and a causal rule mining model that can extract causal dependencies between the cause-effect items.

\subsection{Problem Setting}

We denote the set of users as $\mathcal{U}=\{u_1,u_2,\cdots, u_{|\mathcal{U}|} \}$ and set of items as $\mathcal{I}=\{i_1,i_2.\cdots,i_{|\mathcal{I}|}\}$.
Each user $u$ is associated with a purchase history represented as a sequence of items $\mathcal{H}^u$.
The $j$-th interacted item in the history is denoted as $H_j^u \in \mathcal{I}$.
A black-box sequential recommendation model $\mathcal{F}:\mathcal{H}\rightarrow\mathcal{I}$ is a function that takes a sequence of items (as will discuss later, it can be the permuted user history) as input and outputs the recommended item.
In practice, the underlying mechanism usually consists of two steps: a ranking function that scores all candidate items based on the user history, and then selects the item with best score as the final output.
Note that it only uses user-item interaction without any content or context information, and the scores predicted by the ranking function may differ according to the tasks (e.g. $\{1,\dots,5\}$ for rating prediction, while $[0,1]$ for CTR prediction).
Our goal is to find an item-level post-hoc model that captures the causal relation between the history items and the recommended item for each user.

\begin{mydef}
(Causal Relation) For two variables $X$ and $Y$, if $X$ triggers $Y$,
then we say that there is a causal relation $X \Rightarrow Y$, where $X$ is the \textbf{cause} and $Y$ is the \textbf{effect}.
\end{mydef}

When a given recommendation model $\mathcal{F}$ maps a user history $\mathcal{H}^u$ to a recommended item  $Y^u\in\mathcal{I}$, all items in $\mathcal{H}^u$ are considered as potential causes of $Y^u$.
Thus we can formulate the set of causal relation candidates as $\mathcal{S}^u = \{(H,Y^u)|H\in\mathcal{H}^u\}$.

\begin{mydef}
(Causal Explanation for Sequential Recommendation Model) Given a causal relation candidate set $\mathcal{S}^u$ for user $u$, if there exists a true causal relation $(H,Y^u) \in\mathcal{S}^u$, then the causal explanation for recommending $Y^u$ is described as ``Because you purchased $H$, the model recommends you $Y^u$'', denoted as $H\Rightarrow Y^u$.
\end{mydef}

Then the remaining problem is how to determine whether a candidate pair is a true causal relation.
We can mitigate the problem by allowing a likelihood estimation for a candidate pair to be causal relation.
In other words, we would like to find a ranking function that predicts the likelihood for each candidate pair.
In this way, causal explanations can be generated by selecting the most promising ones from these candidates.
One advantage of this formulation is that it allows the possibility that there is no causal relation between a user's history and the recommended item,
e.g. when algorithm recommends the most popular items regardless of the user history.
In the following sections, we will illustrate in detail how our model solves these problems.

\begin{figure*}[htbp]
    \vspace{-10pt}
    \centering
    \includegraphics[width=\linewidth]{./fig/Model.pdf}
    \vspace{-20pt}
    \caption{Model framework. $x$ is the concatenation of the item embeddings of the user history. $\tilde{x}$ is the perturbed embedding.

    }
    \vspace{-10pt}
    \label{fig:model}
\end{figure*}

\begin{algorithm}[t]
    \caption{Causal Post-hoc Explanation Model}
    \begin{flushleft}
    \textbf{Input:} users $\mathcal{U}$, items $\mathcal{I}$, user history $\mathcal{H}^{u}$, perturbation times $m$, black-box model $\mathcal{F}$, embedding model $\mathcal{E}$, causal mining model $\mathcal{M}$\\

    \textbf{Output:} causal explanations $H\Rightarrow Y^{u}$ where $H\in\mathcal{H}^{u}$\\
    \end{flushleft}
    \begin{algorithmic}[1]
        \STATE Use embedding model $\mathcal{E}$ to get item embeddings $\mathcal{E}(\mathcal{I})$
        \STATE Use $\mathcal{E}(\mathcal{I})$ and true user history

        to train perturbation model $\mathcal{P}$
        \FOR{each user $u$}
            \FOR{$i$ from $1$ to $m$}
                \STATE $\tilde{\mathcal{H}}^u_i \leftarrow \mathcal{P}(\mathcal{H}^{u})$; $\tilde{Y}^u_i\leftarrow\mathcal{F}(\tilde{\mathcal{H}}^u_i)$

            \ENDFOR

            \STATE Construct perturbed input-output pairs $\{(\tilde{\mathcal{H}}^u_i,\tilde{Y}^u_i)\}_{i=1}^m$
            \STATE
            $\theta_{\tilde{H}^u_{ij},\tilde{Y}^u_{i}} \leftarrow \mathcal{M}\Big(\{(\tilde{\mathcal{H}}^u_i,\tilde{Y}^u_i)\}_{i=1}^m \cup (\mathcal{H}^u,Y^u)\Big)$

            \STATE Rank $\theta_{\tilde{H}^u_{ij},Y^u}$ and select top-$k$ pairs $\{(H_j,Y^u)\}_{j=1}^k$
            \IF{$\exists H_{\min\{j\}}\in \mathcal{H}^{u}$}
                \STATE Generate causal explanation $H_{\min\{j\}}\Rightarrow Y^u$
            \ELSE
                \STATE No explanation for the recommended item $Y^u$
            \ENDIF
        \ENDFOR
        \RETURN all causal explanations $H\Rightarrow Y^{u}$
    \end{algorithmic}
    \label{alg:model}
\end{algorithm}

\subsection{Causal Model for Post-Hoc Explanation}

We thus introduce our causal explanation framework for recommendation.
Inspired by \cite{alvarez2017causal}, we divide our framework into two models: a perturbation model and a causal rule mining model.
The overview of the model framework is shown in Fig.\ref{fig:model}.
Before introducing our framework in detail, we define an important concept:
\begin{mydef}
(Causal Dependency): For a given pair of causal relation candidate $(H,Y^u)$, the causal dependency of the pair is the likelihood of the pair being a true causal relation.
\end{mydef}

\subsubsection{\textbf{Perturbation Model}}
\label{sec:perturbation-model}

To capture the causal dependency between items in history and the recommended items, we want to know what would happen if the user history had been different.
To avoid unknown influences caused by the length of input sequence (i.e., user history), we keep the input length unchanged, and only replace items in the sequence to create different sequences.
Ideally, for each item $H_j^u$ in a user's history $\mathcal{H}^u$, it will be replaced by all possible items in $\mathcal{I}$ to fully explore the influence that $H_j^u$ makes in the history.
However, the number of possible combinations will become impractical for the learning system, since recommender systems usually deal with hundreds of thousands or even millions of items.
Therefore, we pursue a perturbation-based method that samples the sequences, which randomly replaces items in the original user history $\mathcal{H}^u$.

There are various ways to obtain the perturbed version of user history, as long as they are similar to the original history.
The simplest solution is randomly selecting an item in $\mathcal{H}^u$ and replacing it with a randomly selected item from $\mathcal{I}$.
However, user histories are far from random selections.
Instead, a natural assumption is that the user behavior follows a certain distribution.
Here, we adopt VAE to learn such a distribution.
As shown in Figure \ref{fig:model}, we design a VAE-based perturbation method, which creates item sequences that are similar to but slightly different from a user's true behavior sequence, by sampling from a distribution in latent embedding space centered around the user's true sequence.

In detail, the VAE component consists of a probabilistic encoder $(\mu,\sigma) = \mathrm{ENC}(\mathcal{X})$ and a decoder $\tilde{\mathcal{X}} = \mathrm{DEC}(z)$.
The encoder encodes a sequence of item embeddings $\mathcal{X}$ into latent embedding space, and extracts the variational information for the sequence, i.e., mean and variance of the latent embeddings under independent Gaussian distribution.
The decoder generates a sequence of item embeddings $\tilde{\mathcal{X}}$ given a latent embedding $z$ sampled from the Gaussian distribution.
Here, both $\mathcal{X}$ and $\tilde{\mathcal{X}}$ are ordered concatenations of pre-trained item embeddings based on pair-wise matrix factorization (BPRMF) \cite{rendle2012bpr}.
We follow the standard training regime of VAE by maximizing the variational lower bound of the data likelihood \cite{kingma2013auto}.
Specifically, the reconstruction error involved in this lower bound is calculated by a softmax across all items for each position of the input sequence.
We observe that VAE can reconstruct the original data set accurately, while offering the power of perturbation.

After pretraining  $\mathrm{ENC}(\cdot)$ and $\mathrm{DEC}(\cdot)$, the variational nature of this model allow us to obtain perturbation $\tilde{\mathcal{H}}$ for any history $\mathcal{H}$.
More specifically, it first extracts the mean and variance of the encoded item sequences in the latent space, and the perturbation model samples $m$ latent embeddings $z$ based on the above variational information.
These sampled embeddings $z$ are then passed to the decoder $\mathrm{DEC}(\cdot)$ to obtain the perturbed versions $\tilde{\mathcal{X}}$.
For now, each item embedding in $\tilde{\mathcal{X}}$ may not represent an actual item since it is a sampled vector from the latent space, as a result, we find its nearest neighbor in the candidate item set $\mathcal{I}$ through dot product similarity as the actual item. In this way, $\tilde{\mathcal{X}}$ is transformed into the final perturbed history $\tilde{\mathcal{H}}$.
One should keep in mind that the variance should be kept small during sampling, so that the resulting sequences are similar to the original sequence.

Finally, the generated perturbed data $\tilde{\mathcal{H}}$ together with the original $\mathcal{H}$ will be injected into the black-box recommendation model $\mathcal{F}$ to obtain the recommendation results $\tilde{Y}$ and $Y$, correspondingly.
After completing this process, we will have $m$ perturbed input-output pairs: $\{(\tilde{\mathcal{H}}^u_i,\tilde{Y}^u_i)\}_{i=1}^m$, as well as the original pair $(\mathcal{H}^u,Y^u)$.

\subsubsection{\textbf{Causal Rule Learning Model}}
\label{sec:causal_rule_learning}
Denote $\mathcal{D}^u$ as the combined records of perturbed input-output pairs $\{(\tilde{\mathcal{H}}^u_i,\tilde{Y}^u_i)\}_{i=1}^m$ and the original pair $(\mathcal{H}^u,Y^u)$ for user $u$.
We aim to develop a causal model that first extracts causal dependencies between input and output items appeared in $\mathcal{D}^u$, and then selects the causal rule based on these inferred causal dependencies.

Let $\tilde{\mathcal{H}}^u_i=[\tilde{H}^u_{i1}, \tilde{H}^u_{i2}, \cdots, \tilde{H}^u_{in}]$ be the input sequence of the $i$-th record of $\mathcal{D}^u$.
Let $\tilde{Y}^u_i$ represent the corresponding output.
Note that this includes the original pair $(\mathcal{H}^u,Y^u)$.
The model should be able to infer the causal dependency $\theta_{\tilde{H}^u_{ij},\tilde{Y}^u_{i}}$ between input item $\tilde{H}^u_{ij}$ and output item $\tilde{Y}^u_{i}$.
We consider that the occurrence of a single output can be modeled as a logistic regression model on causal dependencies from all the input items in the sequence:
\begin{equation}
    P(\tilde{Y}^u_{i}|\tilde{\mathcal{H}}^u_{i}) = \sigma\Big(\sum_{j=1}^n \theta_{\tilde{H}^u_{ij},\tilde{Y}^u_{i}}\cdot\gamma^{n-j}\Big)\label{con:eq1}
\vspace{-5pt}
\end{equation}
where $\sigma$ is the sigmoid function defined as $\sigma(x) = (1+\exp(-x))^{-1}$ in order to scale the score to $[0,1]$.
Additionally, in recommendation task, the order of a user's previously interacted items may affect their causal dependency with the user's next interaction.
A closer behavior tends to have a stronger effect on user's future behaviors, and behaviors are discounted if they happened earlier~\cite{hidasi2015session}.
Therefore, we involve a weight decay parameter $\gamma$ to represent the time effect. Here $\gamma$ is a positive value less than one.

For an input-output pair in $\mathcal{D}^u$, the probability of its occurrence generated by Eq.\eqref{con:eq1} should be close to one. As a result, we learn the causal dependencies $\theta$ by maximizing the probability over $\mathcal{D}^u$.
When optimizing $\theta$, they are always initialized as zero to allow for no causation between two items.
By learning this regression model, we are able to gradually increase $\theta$ until they converge to the point where the data likelihood of $\mathcal{D}$ is maximized.

After gathering all the causal dependencies, we select the items that have high $\theta$ scores to build causal explanations.
This involves a three-step procedure.
\begin{itemize}[leftmargin=*]
\item[1.] We select those causal dependencies $\theta_{\tilde{H}^u_{ij},\tilde{Y}^u_{i}}$ whose output is the original output $Y^u$ (i.e., $\tilde{Y}^u_{i}=Y^u$). Note that these $(\tilde{H}^u_{ij},Y^u)$ pairs may come from either the original sequence or perturbed sequences, because when a perturbed sequence is fed into the black-box recommendation model, the output may happen to be the same as the original sequence $Y^u$.
\item[2.] We sort the above selected causal dependencies in descending order and take the top-$k$ $(\tilde{H}^u_{ij},Y^u)$ pairs.
\item[3.] If there exist one or more pairs in these top-$k$ pairs, whose cause item $\tilde{H}^u_{ij}$ appears in the user's original input sequence $\mathcal{H}^u$, then we pick such pair of the highest rank, and construct $\tilde{H}^u_{ij}\Rightarrow Y^u$ as the causal explanation for the given user.
Otherwise, i.e., no cause item appears in the user history, then we output no causal explanation for the user.
\end{itemize}

Note that the extracted causal explanation is personalized since the algorithm is applied on $\mathcal{D}^u$, which only contains records centered around the user's original record $(\mathcal{H}^u,Y^u)$, while collaborative learning among users is indirectly modeled by the VAE-based perturbation model. The overall algorithm is provided in Alg.\ref{alg:model}.

\section{Experiments}
In this part, we conduct experiments to show what causal relationships our model can capture and how they can serve as an intuitive explanation for the black-box recommendation model.

\subsection{Dataset Description}
We evaluate our proposed causal explanation framework against baselines on two datasets. The first dataset is MovieLens100k\footnote{https://grouplens.org/datasets/movielens/}.
This dataset consists of information about users, movies and ratings. In this dataset, each user has rated at least 20 movies, and each movie can belong to several genres. The second dataset is the office product dataset from Amazon\footnote{https://nijianmo.github.io/amazon/}, which contains the user-item interactions from May 1996 to July 2014. The original dataset is 5-core. To achieve sequential recommendation with input length of 5, we select the users with at least 15 purchases and the items with at least 10 interactions.

Since our framework is used to explain sequential recommendation models, we split the dataset chronologically.
To learn the pre-trained item embeddings based on BPRMF \cite{rendle2012bpr} (section \ref{sec:perturbation-model}), we take the last 6 interactions from each user to construct the testing set, and use all previous interactions from each user as the training set.
To avoid data leakage, when testing the black-box recommendation models and our VAE-based perturbation model, we only use the last 6 interactions from each user (i.e., the testing set of the pre-training stage).
Following common practice,
we adopt the leave-one-out protocol, i.e., among the 6 interactions, we use the last one for testing, and the previous 5 as input to recommendation models.
A brief summary of the data is shown in Table \ref{tab:data}.

\begin{table}[t]
    \centering

    \setlength{\tabcolsep}{1pt}
    \begin{tabular}{ccccccc}
        \toprule
        Dataset & \# users & \# items & \# interactions & \# train & \# test & sparsity \\\midrule
        Movielens & 943 & 1682 & 100,000 & 95,285 & 14,715 & 6.3\%\\
        Amazon & 573 & 478 & 13,062 & 9,624 & 3,438 & 4.7\%\\\bottomrule
    \end{tabular}
    \caption{Summary of the Datasets}
    \vspace{-10pt}
    \label{tab:data}
    \vspace{-10pt}
\end{table}

\begin{table*}[t]
    \centering
    \setlength{\tabcolsep}{4pt}
    \begin{tabular}{ccccccccccccc}
        \toprule
        Dataset & \multicolumn{6}{c}{Movielens 100k} & \multicolumn{6}{c}{Amazon}\\\cmidrule(lr){1-1}\cmidrule(lr){2-7}\cmidrule(lr){8-13}
        Method & \multicolumn{3}{c}{Causal} & \multicolumn{3}{c}{Association} & \multicolumn{3}{c}{Causal} & \multicolumn{3}{c}{Association}\\\cmidrule(lr){2-4}\cmidrule(lr){5-7}\cmidrule(lr){8-10}\cmidrule(lr){11-13}
        Parameter & $k=1$ & $k=2$ & $k=3$ & Support & Confidence & Lift & $k=1$ & $k=2$ & $k=3$ & Support & Confidence & Lift\\\midrule
        FPMC & 96.50\% & 99.25\% & 99.57\% & 14.32\% & 13.47\% & 13.47\%& 95.11\% & 98.95\% & 99.82\% & 12.57\% & 12.04\% & 12.04\%\\
        GRU4Rec & 98.51\% & 99.36\% & 99.68\% & 7.423\% & 7.635\% & 7.635\% & 95.94\% & 99.48\% & 99.65\% & 8.202\% & 8.551\% & 8.551\%\\
        Caser & 97.03\% & 99.15\% & 99.57\% & 8.801\% & 8.165\% & 8.165\% & 95.99\% & 99.30\% & 99.65\% & 9.250\% & 8.901\% & 8.901\%\\
        \bottomrule
    \end{tabular}
    \caption{Results of Model Fidelity. Our causal explanation framework is tested under the number of candidate causal explanations $k=1,2,3$. The association explanation framework is tested under support, confidence, and lift thresholds, respectively.}
    \vspace{-15pt}
    \label{tab:fidelity}
\end{table*}

\subsection{Experimental Settings}
\label{sec:experimental_setting}

We adopt the following methods to train black-box sequential recommendation models and to extract traditional association rules as comparative explanations. We include both shallow and deep models for experiment.

\textbf{FPMC}~\cite{rendle2010factorizing}: The Factorized Personalized Markov Chain model, which combines matrix factorization and Markov chains to capture user's personalized sequential behavior patterns for prediction\footnote{https://github.com/khesui/FPMC}.

\textbf{GRU4Rec}~\cite{hidasi2015session}: A session-based recommendation model, which uses recurrent neural networks -- in particular, Gated Recurrent Units (GRU) -- to capture sequential patterns for prediction\footnote{https://github.com/hungthanhpham94/GRU4REC-pytorch}.

\textbf{Caser}~\cite{tang2018personalized}: The ConvolutionAl Sequence Embedding Recommendation (Caser) model, which adopts convolutional filters over recent items to learn the sequential patterns for prediction\footnote{https://github.com/graytowne/caser\_pytorch}.

\textbf{Association Rule}~\cite{peake2018explanation}: A post-hoc explanation model, which learns the item-item association rules as item-level explanations\footnote{https://pypi.org/project/apyori/}.

For black-box recommendation models FPMC, GRU4Rec and Caser, we adopt their best parameter selection in their corresponding public implementation. For the association rule-based explanation model, we follow the recommendations in \cite{peake2018explanation} to set the optimal parameters: support = 0.1, confidence = 0.1, lift = 0.1, length = 2 for \textit{MovieLens100k}, and support = 0.01, confidence = 0.01, lift = 0.01, length = 2 for \textit{Amazon} dataset due to its smaller scale.

For our causal rule learning framework, we set the item embedding size as 16, both the VAE encoder and decoder are Multi-Layer Perceptrons (MLP) with two hidden layers, and each layer consists of 1024 neurons. The default number of perturbed input-output pairs is $m=500$ on both datasets. The default time decay factor is $\gamma=0.7$. We will discuss the influence of perturbation times $m$ and time decay factor $\gamma$ in the experiments.

In the following, we will apply both association rule learning and causal rule learning frameworks on the black-box recommendation models to evaluate and compare the association explanations and causal explanations. In particular, we evaluate our framework from two perspectives. First, we verify that the causal rules learned by our framework represent highly probable causal relationships (explanation quality). Second, we show that our model has the ability to offer explanations for most recommendations (explanation fidelity). Additionally, we shed light on how our model differs from other models on statistical metrics.

\iffals
\subsection{Implementation Details}
In this subsection, we introduce the implementation details to enhance reproducibility, including both our model and the baselines.

FPMC \cite{rendle2010factorizing} is implemented by public code\footnote{https://github.com/khesui/FPMC}. During the training process, all transactions in the training set are organized as sequences, where each sequence is the purchased items ordered by time for one user.
The original code put some users into training set while the rest users in the testing set. However, in our case, to keep consistent to all baselines, we divided the latest six interactions as testing set where the latest one is the target of prediction. For testing, the input is one user and a sequence with length five which represents the first five interacted items of the corresponding user. Since here we did not focus on the performance of recommendation, we use the same parameter setting for two datasets which is the same as the original setting in public code.

GRU4Rec\cite{hidasi2015session} is another sequential recommendation model. It is implemented by a pytorch version of public code from GitHub \footnote{https://github.com/hungthanhpham94/GRU4REC-pytorch}. This code is based on  paper\cite{hidasi2015session} and is extended on \cite{hidasi2018recurrent}. Here we group session by users. We consider the purchase history of one user is a session. The user id in the dataset is considered as the session id for GRU4Rec. Similar to FPMC, the dataset used in public code is split by session, however, considering that some sequential recommendation models may learn the representation of users, we split each user history into training and testing set to keep consistent to all baselines.

Caser\cite{tang2018personalized} is one of the state-of-the-art sequential recommendation models. This model is implemented by the code\footnote{https://github.com/graytowne/caser\_pytorch} published by authors. Caser also tests on \textit{Movielens 100k} dataset, therefore, we basically use the same parameter setting except the number of target. Compared with original one, we only predict the next one behavior instead of three. The parameter settings are same for two datasets.

To distinguish between association rule and causal rule, we firstly introduce association rule mining briefly. Association rules are extracted from items according to their strength in the collection. There are typically three scores to measure association rules: the support, confidence and lift. For two items $m$ and $n$, the support is a frequency measure to show how often we can expect to observe the pattern $\{ m \rightarrow n\}$; the confidence is an estimate of the conditional probability $P(m|n)$ to represent the predictive ability of the rule; and the lift measures the dependency between two items $m$ and $n$ . Specifically, a rule $\{ m \rightarrow n\}$ holds with support $supp$ if $supp\%$ of the records contain both $m$ and $n$; with confidence $conf$ if $conf\%$ of records that contain $m$ also contain $n$, which is defined as $supp(m \cup n) / supp(m)$. And the rule's lift is defined as $supp(m \cup n)$ divided by $supp(m) \times supp(n)$ if the items were independent.

The association rules are generally mined by setting specific thresholds for those three scores. Here we followed the pipeline offered by  \cite{peake2018explanation} and the association rules are extracted by apyori\footnote{https://pypi.org/project/apyori/} python package. For the association rules, we only involve global association rules and select top 30 for each interestingness as explanation. Peake's work\cite{peake2018explanation} also uses local association rule as explanation which requires to extract association rule based on neighborhood users with high cosine similarity. However, some sequential recommendation models may not have user representation which means that we cannot get the cosine similarity to get the neighborhood users. Based on that consideration, we only apply the global association rule as explanation. Meanwhile, we set the threshold same as \cite{peake2018explanation}: min\_supp = 0.1, min\_conf = 0.1, min\_lift = 0.1, max\_len = 2 for \textit{Movielens 100k}, and min\_supp = 0.01, min\_conf = 0.01, min\_lift = 0.01, max\_len = 2 for \textit{Amazon} due to its smaller scale.

For our framework, first we need to pretrain those black-box models and fix them after training. Before training the perturbation model, we also need to get the latent representation of each items. Here we simply use pair-wised \textit{Matrix Factorization} to learn the embedding of items. After fixing the embedding of items, we use the sequence with length five in testing set to train a perturbation model and fix it. The latent size for perturbation model is 16, and the encoder and decoder are two Multi-Layer Perceptron (MLP) with two hidden layer and each layer is consisted of 1024 neurons. For each user, we need several sequences which are similar to the true history sequence for achieving our random experiment. In \textit{Movielens 100k} dataset, we generate 500 unique perturbed sequence for each user. Meanwhile, we keep the original sequence to make sure all items in true history will appear at least once. We plug unique sequence (501 in total) into black box models to get corresponding recommendation as perturbed output which has the exactly same number (501 in total). In \textit{Amazon} dataset, the number of items is less than \textit{Movielens 100k}. Based on that observation, we only sample 250 times for each user which means that there are 251 unique sequences including the original one for each user. We will discuss about the influence of time of perturbation in later section. Same as \textit{Movielens 100k}, we plug the sequences into recommendation model to get corresponding predictions. After getting input-output perturbation pairs, we use our causal rule mining model to extract the causal relations set and select the causal explanation from the set. The time effect is set as 0.7 and we will discuss about the influence of it in later section.

In following part, we evaluate our framework from two perspectives. First, we verify that the causal rules which are mined by our framework are highly probably causal relationships. Second, we show that our model has the ability to give explanations for most recommendations. Additionally, we shed light on the characteristics of our model differs from other model on statistical metrics.

\f

\subsection{Causality Verification}
We first verify the quality of the extracted causal rules. Here, we adopt the following widely used definition of causation, which is introduced by Pearl \cite{pearl2000causality}:
\begin{equation}\label{eq:verification}
     \Pr(\text{effect}|do(\text{cause})) \textgreater \Pr(\text{effect}|do(\neg\text{cause}))
\end{equation}
where $do(c)$ represents an external intervention, which compels the truth of $c$, and $do(\neg c)$ compels the truth of not $c$. Actually, the conditional probability $\Pr(e|c)$ represents a probability resulting from a passive observation of $c$, which rarely coincides with $\Pr(e|do(c))$. We evaluate how many percentage of our extracted causal explanations really satisfy Eq.\eqref{eq:verification}.

Suppose the perturbation model (section \ref{sec:perturbation-model}) creates $m$ perturbed input-output pairs for each user $u$: $\{(\tilde{\mathcal{H}}^u_i,\tilde{Y}^u_i)\}_{i=1}^m$, plus the original pair $(\mathcal{H}^u,Y^u)$. Here $\tilde{\mathcal{H}}^u$ is created by our perturbation model (i.e. not observed in the original data), and thus observing $\tilde{\mathcal{H}}^u$ implies we have $do(\tilde{\mathcal{H}}^u)$ in advance. Let $H\Rightarrow Y^u$ be the causal explanation extracted by the casual rule learning model (section \ref{sec:causal_rule_learning}). Then we estimate the probability based on these $m+1$ total pairs as,
\begin{equation}
\small
\begin{aligned}
    &\Pr(\text{effect}|do(\text{cause}))=\Pr(Y^u|H)=\frac{\#\text{Pairs}(H\in\tilde{\mathcal{H}}^u\land Y=Y^u)}{\#\text{Pairs}(H\in\tilde{\mathcal{H}}^u)}\\
    &\Pr(\text{effect}|do(\neg\text{cause}))=\Pr(Y^u|\neg H)=\frac{\#\text{Pairs}(H\notin\tilde{\mathcal{H}}^u\land Y=Y^u)}{\#\text{Pairs}(H\notin\tilde{\mathcal{H}}^u)}
\end{aligned}
\end{equation}
and $H\Rightarrow Y^u$ is considered as a reliable causal rule if Eq.\eqref{eq:verification} is satisfied. The intuition here is that a causal rule should guarantee that the probability of seeing $Y^u$ when $H$ is purchased should be higher than the probability of seeing $Y^u$ when $H$ is not purchased.

We calculate the percentage of our extracted causal explanations that satisfy Eq.\eqref{eq:verification}. We tune the parameter $k$ from 1 to 3, where $k$ is the number of candidate causal explanations we consider for each user in the causal rule learning model (section \ref{sec:causal_rule_learning} step 2 and 3).
The results of causality verification on two datasets and three recommendation models are shown in Table \ref{tab:verification}.

Based on the results we can see that in most cases 90\% or more of the extracted causal explanations are reliable. We also see that with the increasing of $k$, the percentage tends to decease for all three recommendation models on both datasets. Intuitively, this means that unreliable causal relation candidates may be introduced if too many candidate pairs are considered in the causal rule learning model, which is reasonable. However, if too few candidates are considered (i.e., $k$ is too small), many users may not receive explanations at all (i.e., model fidelity will decrease). This will be further analyzed in the next subsection.

\begin{table}[t]
    \centering
    \setlength{\tabcolsep}{3pt}
    \begin{tabular}{ccccccc}
        \toprule
        Dataset & \multicolumn{3}{c}{Movielens} & \multicolumn{3}{c}{Amazon}\\\cmidrule(lr){2-4}\cmidrule(lr){5-7}
        $k$ & 1 & 2 & 3 & 1 & 2 & 3\\\midrule
        FPMC & 94.83\% & 94.23\% & 93.92\% & 91.56\% & 90.48\% & 90.03\%\\
        GRU4Rec & 97.84\% & 97.76\% & 97.75\% & 94.91\% & 94.20\% & 93.87\%\\
        Caser & 97.15\% & 96.89\% & 96.69\% & 98.91\% & 98.25\% & 98.14\%\\\bottomrule
    \end{tabular}
    \caption{The percentage of reliable causal explanations that satisfy the inequality in Eq.\eqref{eq:verification}.}
    \vspace{-27pt}
    \label{tab:verification}
\end{table}

\subsection{Model Fidelity}
An important evaluation measure for explanation models is model fidelity, i.e., how many percentage of the recommendation results can be explained by the model \cite{zhang2018explainable}.
The results of model fidelity are shown in Table \ref{tab:fidelity}. In this experiment, we still tune the number of candidate causal explanations $k$ from 1 to 3. For the association rule explanation model (section \ref{sec:experimental_setting}), we test three versions of the association rule learning algorithm as introduced in \cite{peake2018explanation}, i.e., the association rules are filtered by support, confidence, and lift thresholds, respectively.

We can see that on both datasets, our causal explanation framework can generate explanations for almost all of the recommended items, while the association explanation approach can only provide explanations for significantly fewer recommendations. The underlying reason is that association explanations have to be extracted based on the original input-output pairs, which limits the number of pairs that we can use for rule extraction. However, based on the perturbation model, our causal explanation framework is able to create many counterfactual examples to assist causal rule learning, which makes it possible to go beyond the limited original data to extract causal explanations.

Another interesting observation is that GRU4Rec and Caser have significantly ($p<0.01$) lower fidelity than FPMC when explained by the association model. This is reasonable because FPMC is a Markov-based model that directly learns the correlation of adjacent items in a sequence, as a result, it is easier to extract association rules between inputs and outputs for the model. However, it also means that the fidelity performance of the association approach highly depends on the recommendation model being explained. Meanwhile, we see that our causal approach achieves comparably good fidelity on all three recommendation models, because the perturbation model is able to create sufficiently many counterfactual examples to break the correlation of frequently co-occurring items in the input sequence. This indicates the robustness of our causal explanation framework in terms of model fidelity.

\iffalse
\begin{table}[htbp]
    \centering
    \footnotesize
    \begin{tabular}{ccccccc}
        \toprule
        \multicolumn{7}{c}{Movielens 100k}\\\midrule
        method & \multicolumn{3}{c}{Causal} & \multicolumn{3}{c}{Global Association}\\\cmidrule(lr){2-4}\cmidrule(lr){5-7}
        k & 1 & 2 & 3 & Support & Confidence & lift\\\midrule
        Caser & 85.365\% & 92.895\% & 95.546\% & 8.801\% & 8.165\% & 8.165\%\\
        GRU4Rec & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%}\\
        FPMC & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%}\\\bottomrule
        \multicolumn{7}{c}{Amazon}\\\midrule
        method & \multicolumn{3}{c}{Causal} & \multicolumn{3}{c}{Global Association}\\\cmidrule(lr){2-4}\cmidrule(lr){5-7}
        k & 1 & 2 & 3 & Support & Confidence & lift\\\midrule
        Caser & \fillblank{87.72\%}\% & \fillblank{87.72\%}\% & \fillblank{87.72\%}\% & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%}\\
        GRU4Rec & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%}\\
        FPMC & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%} & \fillblank{87.72\%}\\\bottomrule
    \end{tabular}
    \caption{Results of Model Fidelity}
    \label{tab:fidelity}
\end{table}
\fi

\iffals
\begin{table}[t]
    \centering
    \small
    \begin{tabular}{ccccccc}
        \toprule
        \multirow{2}{*}{Rules} & \multicolumn{3}{c}{Causal Rules} & \multicolumn{3}{c}{Global Association}\\\cmidrule(lr){2-4}\cmidrule(lr){5-7}
         & 1 & 2 & 3 & support & confidence & lift\\\midrule
        FPMC & 37.47\% & 37.39\% & 37.38\% & 49.62\% & 50.39\% & 50.39\%\\
        GRU4Rec & 40.90\% & 40.98\% & 40.97\% & 55.71\% & 61.11\% & 61.11\%\\
        Caser & 37.92\% & 37.75\% & 37.91\% & 51.80\% & 51.94\% & 51.94\%\\\bottomrule
    \end{tabular}
    \caption{The percentage of explanations that have same feature as recommendations in \textit{Movielens 100k}}
    \vspace{-25pt}
    \label{tab:feature}
\end{table}
\f

\iffals
\subsection{Explanation Observation}
How to evaluate explainability is an ill-defined problem. Our model tries to explain the behavior of the black-box models. In this case, the explanation may not be accepted by all of humans. However, we would like to show some characteristics of our model by some statistical observation comparing with baselines. Here we use \textit{Movielens 100k} as an example to show our observation.

We count the percentage of common features for explanation and recommendation. Here, for \textit{Movielens} dataset, we use the genres of movies as the feature of movies. The results of percentage of explanations that have the same feature as the recommendation are shown in Table.\ref{tab:feature}. According to the table, such percentage is higher for the explanations extracted from global association rule. However, combining with the model fidelity shown in Table.\ref{tab:fidelity}, we can see that global association rule can only explain a few portions of the recommendation, therefore, even it has a more percentage of common feature explanations, our causal model can give more explanations sharing common features with recommendations in terms of the total number.

Though sharing common features seems to make an explanation more reasonable and intelligible, it can not be guaranteed that such explanations are better or more accurate. Our causal model is capable of offering more diverse and inspiring explanations to users compared with the global association rule model, which will be shown through examples in subsection~\ref{sec:explanation}.

\f

\subsection{Influence of Parameters}
In this section, we discuss the influence of two important parameters. The first one is time decay parameter $\gamma$ -- in our framework, when explaining the sequential recommendation models, earlier interactions in the sequence will have discounted effects to the recommended item. A proper time decay parameter helps the framework to reduce noise signals when learning patterns from the sequence. The second parameter is the number of perturbed input-output pairs $m$ -- in our framework, we use perturbations to create counterfactual examples for causal learning, but there may exist trade-off between efficiency and performance. We will analyze the influence of these two parameters.

\begin{figure}[t]
\captionsetup[sub]{font=small,labelfont=normalfont,textfont=normalfont}
\centering
\hspace{-15pt}
\begin{subfigure}{0.24\textwidth}
\label{fig:GRU4Rec_ml100k_time_decay}
        \includegraphics[scale=0.3]{fig/GRU4Rec_fidelity_time_decay_ml100k.png}
        \vspace{-20pt}
        \caption{GRU4Rec on Movielens}
\end{subfigure}
\begin{subfigure}{0.24\textwidth}
\label{fig:GRU4Rec_amazon_time_decay}
        \includegraphics[scale=0.3]{fig/GRU4Rec_fidelity_time_decay_amazon.png}
        \vspace{-20pt}
        \caption{GRU4Rec on Amazon}
\end{subfigure}
\medskip
\hspace{-15pt}
\begin{subfigure}{0.24\textwidth}
\label{fig:k3_ml100k_time_decay}
        \includegraphics[scale=0.3]{fig/k3_fidelity_time_decay_ml100k.png}
        \vspace{-20pt}
        \caption{$k=3$ on Movielens}
\end{subfigure}
\begin{subfigure}{0.24\textwidth}
\label{fig:k3_amazon_time_decay}
        \includegraphics[scale=0.3]{fig/k3_fidelity_time_decay_amazon.png}
        \vspace{-20pt}
        \caption{$k=3$ on Amazon}
\end{subfigure}
\vspace{-10pt}
\caption{Model fidelity on different time decay parameters $\gamma$. $x$-axis is the time decay parameter $\gamma\in\{0.1, 0.3, 0.7, 1\}$ and $y$-axis is the model fidelity.
The upper two sub-figures plot the model fidelity for GRU4Rec on different number of candidate causal explanations $k=1,2,3$, and the lower two sub-figures plot the model fidelity for different recommendation models when fixing $k=3$.}
\label{fig:time}
\vspace{-10pt}
\end{figure}

\textbf{Time Decay Effect}:
Figure \ref{fig:time} shows the influence of $\gamma$ on different recommendation models and datasets. From the result we can see that the time decay effect $\gamma$ indeed affects the model performance on fidelity. In particular, when $\gamma$ is small, the previous interactions in a sequence are more likely to be ignored, which thus reduces the performance on model fidelity. When $\gamma$ is large (e.g., $\gamma=1$), old interactions will have equal importance with latest interactions, which also hurts the performance. We can see from the results that the best performance is achieved at about $\gamma=0.7$ on both datasets.

\begin{figure}[t]
\captionsetup[sub]{font=small,labelfont=normalfont,textfont=normalfont}
\centering
\hspace{-15pt}
\begin{subfigure}{0.24\textwidth}
\label{fig:GRU4Rec_ml100k_time_decay}
        \includegraphics[scale=0.3]{fig/k1_fidelity_time_pert_ml100k.png}
        \vspace{-20pt}
        \caption{Model Fidelity on Movielens}
\end{subfigure}
\begin{subfigure}{0.24\textwidth}
\label{fig:GRU4Rec_amazon_time_decay}
        \includegraphics[scale=0.3]{fig/k1_fidelity_time_pert_amazon.png}
        \vspace{-20pt}
        \caption{Model Fidelity on Amazon}
\end{subfigure}
\medskip
\hspace{-15pt}
\begin{subfigure}{0.24\textwidth}
\label{fig:k3_ml100k_time_decay}
        \includegraphics[scale=0.3]{fig/k1_verification_time_pert_ml100k.png}
        \vspace{-20pt}
        \caption{Percentage on Movielens}
\end{subfigure}
\begin{subfigure}{0.24\textwidth}
\label{fig:k3_amazon_time_decay}
        \includegraphics[scale=0.3]{fig/k1_verification_time_pert_amazon.png}
        \vspace{-20pt}
        \caption{Percentage on Amazon}
\end{subfigure}
\vspace{-10pt}
\caption{Model fidelity and the percentage of verified causal explanations (Eq.\eqref{eq:verification}) on different number of perturbed pairs $m$. $x$-axis is the number of perturbed pairs $m$. $y$-axis is model fidelity. The upper two sub-figures plot the model fidelity for different recommendation model on $k=1$, and the lower two sub-figures plot the percentage of verified causal explanations for different recommendation model on $k=1$.}
\label{fig:pert_time}
\vspace{-10pt}
\end{figure}

\textbf{Number of Perturbations}:
Figure \ref{fig:pert_time} shows the influence for the number of perturbed input-output pairs $m$.
A basic observation from Figure \ref{fig:pert_time} is that when $m$ increases, both model fidelity and the percentage of verified rules will decrease first and then increase. The underlying reason is as follows.

When $m$ is small, the variance of the perturbed input-output pairs will be small, and thus almost any change in the input sequence will be determined as a cause. For example, suppose the original input-output pair is $A, B, C \rightarrow Y$. In the extreme case where $m=1$, we will have only one perturbed pair, e.g., $A, \tilde{B}, C \rightarrow \tilde{Y}$. According to the causal rule learning model (section \ref{sec:causal_rule_learning}), if $\tilde{Y}\neq Y$, then $B\Rightarrow Y$ will be the causal explanation since the change of $B$ results in a different output, while if $\tilde{Y}=Y$, then either $A\Rightarrow Y$ or $C\Rightarrow Y$ will be the causal explanation since their $\theta$ scores will be higher than $B$ or $\tilde{B}$. In either case, the model fidelity and percentage of verified causal rules will be $100\%$. However, in this case, the results do not present statistical meanings since they are estimated on a very small amount of examples.

When $m$ increases but not large enough, then random noise examples created by the perturbation model will reduce the model fidelity. Still consider the above example, if many pairs with the same output $Y$ are created, then the model may find other items beyond $A,B,C$ as the cause, which will result in no explanation for the original sequence. However, if we continue to increase $m$ to sufficiently large numbers, such noise will be statistically offset, and thus the model fidelity and percentages will increase again. In the most ideal case, we would create all of the $|\mathcal{H}|^{|\mathcal{I}|}$ sequences for causal rule learning, where $|\mathcal{H}|$ is the number of item slots in the input sequence, and $|\mathcal{I}|$ is the total number of items in the dataset. However, $|\mathcal{H}|^{|\mathcal{I}|}$ would be a huge number that makes it computational infeasible for causal rule learning. In practice, we only need to specify $m$ sufficiently large. Based on Chebyshev's Inequality, we find that $m=500$ already gives >95\% confidence that the estimated probability error is <0.1.

\subsection{Case Study}
\label{sec:explanation}
In this section, we provide a simple qualitative case study to compare causal explanations and association explanations. Compared with the association explanation model, our model is capable of generating personalized explanations, which means that even if the recommendation model recommends the same item for two different users and the users have overlapped histories, our model still has the potential to generate different explanations for different users. However, the association model will provide the same explanation since the association rules are extracted based on global records.
An example by the Caser recommendation model on \textit{MovieLens100k} dataset is shown in Figure \ref{fig:case_study},
where two users with one commonly watched movie (\textit{The Sound of Music}) get exactly same recommendation (\textit{Pulp Fiction}). The association model provides the overlapped movie as an explanation for the two different users, while our model can generate personalized explanation for different users even when they got the same recommendation.

\iffalse
\begin{table}[h]
    \centering
    \begin{tabular}{l}
    \hline
    \textbf{Reliable Causal Rule Examples}\\\hline
    Dracula: Dead and Loving It $\Rightarrow$ Jeffrey\\
    Free Willy 3: The Rescue $\Rightarrow$ Mad City\\
    Crumb $\Rightarrow$ Evita\\\hline
    \textbf{Fake Causal Rule Examples}\\\hline
    Men in Black $\Rightarrow$ Twilight \\
    Beyond Bedlam $\Rightarrow$ Donnie Brasco\\
    Blood for Dracula $\Rightarrow$ Face/Off\\\hline
    \end{tabular}
    \caption{causal Rules Examples}
    \vspace{-10pt}
    \label{tab:results}
\end{table}
\fi

\vspace{-5pt}
\section{Conclusions}
Recommender systems are widely used in our daily life.
Effective recommendation mechanisms usually work through black-box models, resulting in the lack of transparency.
In this paper, we extract causal rules from user history to provide personalized, item-level, post-hoc explanations for the black-box sequential recommendation models.
The causal explanations are extracted through a perturbation model and a causal rule learning model.
We conduct experiments on real-world datasets, and apply our explanation framework to several state-of-the-art sequential recommendation models.
Experimental results verified the quality and fidelity of the causal explanations extracted by our framework.

In this work, we only considered item-level causal relationships, while in the future, it would be interesting to explore causal relations on feature-level external data such as textual user reviews, which can help to generate finer-grained causal explanations.

\section*{Acknowledgement}
This work was supported in part by NSF IIS-1910154. Any opinions, findings, conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect those of the sponsors.

\begin{figure}[t!]
    \centering
    \includegraphics[width=\linewidth]{./fig/CaseStudy.pdf}
    \vspace{-25pt}
    \caption{

    A case study on MovieLens by the Caser model. The first movie for $u_1$ is unknown in the dataset.
    }
    \vspace{-10pt}
    \label{fig:case_study}
\end{figure}

\vspace{-5pt}
\balance
\bibliographystyle{ACM-Reference-Format}


\section*{Availability}
This article may be downloaded for personal use only. Any other use requires prior permission of the author and AIP Publishing. This article appeared in \cite{Stehli2020} and may be found at \href{https://aip.scitation.org/doi/10.1063/5.0023533}{https://aip.scitation.org/doi/10.1063/5.0023533}. The data that support the findings of this study are available from the corresponding author upon reasonable request.

\section*{References}