\section{Tight-binding Hamiltonians}
For the calculation of the relaxation time and resistivity, we rely on the use of a tight-binding model. We adopt a general formulation of the tight-binding approach with Hamiltonian
\begin{align}
H = -\sum_{\langle i,j\rangle} t_{\parallel} (\r_i-\r_j) \; (a_{1,i}^{\dagger}a_{1,j}+h.c.) - \sum_{\langle i,j\rangle} t_{\parallel} (\r_i-\r_j) \; (a_{2,i}^{\dagger}a_{2,j}+h.c.) - \sum_{(i,j)} t_{\perp}(\r_i-\r_j) \; (a_{1,i}^{\dagger}a_{2,j}+h.c.)\;. 
\label{tbh}
\end{align}
The sum over the brackets $\langle...\rangle$ runs over pairs of atoms in the same layer (1 or 2), whereas the sum over the curved brackets $(...)$ runs over pairs with atoms belonging to different layers. $t_{\parallel} (\r)$ and $t_{\perp} (\r)$ are hopping matrix elements which have an exponential decay with the distance $|\r|$ between carbon atoms. A common parametrization is based on the Slater-Koster formula for the transfer integral\cite{Moon13} 
\begin{align}
-t(\d)=V_{pp\pi}(d)\left[1-\left(\frac{\d\cdot{\bm e}_z}{d}\right)^2\right]+V_{pp\sigma}(d)\left(\frac{\d\cdot{\bm e}_z}{d}\right)^2
\end{align}
with
\begin{align}
V_{pp\pi}(d)=V_{pp\pi}^0\exp\left(-\frac{d-a_0}{r_0}\right)\;,
V_{pp\sigma}(d)=V_{pp\sigma}^0\exp\left(-\frac{d-d_0}{r_0}\right)\;,
\end{align}
where $\d $ is the vector connecting the two sites, ${\bm e}_z$ is the unit vector in the $z$-direction, $a_0 $ is the C-C distance and $d_0$ is the distance between layers. A typical choice of parameters is given by $V_{pp\pi}^0=-2.7$ eV, $V_{pp\sigma}^0=0.48$ eV and $r_0=0.319 a_0$ \cite{Moon13}. In particular, we have taken these values to carry out the analysis reported in the main text. For an alternative comparison between the continuous and the tight-binding model, see Ref. \cite{Stauber18b}.

\begin{figure}[h]
\includegraphics[width=0.20\columnwidth]{fig5a}
\includegraphics[width=0.20\columnwidth]{fig5b}
\caption{Energy contour maps of the second highest valence band in the Moir\'e Brillouin zone of a twisted graphene bilayer with twist angle $\theta_{28}\approx 1.16^\circ$, showing the Fermi lines for filling levels shifted $-0.2$ meV (left) and $-1.5$ meV (right) below the level of the saddle points placed along the $\Gamma K$ lines.  
\label{fline}}
\end{figure}


\section{Transport decay rate of quasi-particles at the Fermi line} 
At low temperature, the transport decay rate is dominated by electron quasiparticles close to the Fermi line. In Fig. \ref{fline}, we show the Fermi line for two different chemical potentials $\Delta \mu$ taken with respect to the level of the van Hove singularity (vHS) arising from the saddle points at the $\Gamma K$ line. Also indicated are the discrete points on the Fermi line for which explicit calculations are here illustrated. 




In Fig. \ref{dr02}, we show the transport decay rate as function of temperature, computed according to the expression (2) in the main text, for the different points on the Fermi line indicated in Fig. \ref{fline}. As can be appreciated, the behavior depends crucially on the value of $\Delta \mu $, turning from linear to quadratic at low temperatures as the Fermi level deviates significantly from the vHS. 



\begin{figure}[h]
\includegraphics[width=0.40\columnwidth]{fig6a}
\hspace{1cm}
\includegraphics[width=0.40\columnwidth]{fig6b}
\\
 \hspace*{1.0cm} (a) \hspace{7.7cm} (b)
\caption{Plot of the temperature dependence of $1/\tau_{\rm tr}$ (weighted with the inverse of the square of the Fermi velocity to get dimensions of energy) when the Fermi level is $0.2$ meV (left) and $1.5$ meV (right) below the vHS, for six points along de Fermi line following the sequence shown in Fig. \ref{fline}, from the farthest position (1) to the closest location to the $K$ point (6). The on-site Hubbard interaction is taken as $U/(2\pi ) = 3$ meV $a_M^2$, $a_M$ being the lattice constant of the superlattice.    
\label{dr02}}
\end{figure}



 
\section{Quasi-particle properties at the Fermi line} 
Also for the self-energy, we can analyse the low-energy behaviour for different quasi-particles on the Fermi line. This is seen in Fig. \ref{imself}(a), which represents the imaginary part of the self-energy $\Sigma (\k, \omega )$ as function of $\omega $, computed according to the expression (5) in the main text, for the points on the Fermi line indicated in Fig. \ref{fline}. The linear behavior at low frequencies is consistent with the low-temperature dependence of the transport decay rate shown in Fig. \ref{dr02} for $\Delta \mu = -0.2$ meV.
\begin{figure}[h]
\includegraphics[width=0.43\columnwidth]{fig7a}
\hspace{1cm}
\includegraphics[width=0.40\columnwidth]{fig7b}
\\
 \hspace*{1.0cm} (a) \hspace{7.7cm} (b)
\caption{Plot of the frequency dependence of the imaginary (left hand side) and real (right hand side) part of $\Sigma (\k, \omega )$ for six points along de Fermi line following the sequence shown in Fig. \ref{fline}, from the farthest position (1) to the closest location to the $K$ point (6). The Fermi level is $0.2$ meV below the vHS and the on-site Hubbard interaction is taken as $U/(2\pi ) = 3$ meV $a_M^2$, $a_M$ being the lattice constant of the superlattice. 
\label{imself}}
\end{figure}


From the imaginary part of $\Sigma (\k, \omega )$, we can compute the real part of the self-energy by applying the Kramers-Kronig relation in Eq. (6) of the main text. The results corresponding to the different curves in Fig. \ref{imself}(a) are represented in Fig. \ref{imself}(b), which shows a clear logarithmic correction consistent with the linear dependence at low frequencies of the respective imaginary counterparts.

\section{Umklapp processes}
Umklapp processes define scattering events for which the final momentum lies outside the first Brillouin zone. The final momentum and its corresponding energy can be mapped back onto the first Brillouin zone, but this is not allowed for its wave function. Nevertheless, mapping also the eigenvectors onto the first Brillouin zone facilitates the numerical calculations and this approximation leads to a susceptibility that is periodic on the first Brillouin zone. 

The above approximation has been employed in the calculations of the main text and shall here be discussed for the continuum model, i.e., we compare it to the exact result. Another approximation would simply neglect all scattering processes that lie outside the first Brillouin zone. In Fig. \ref{Susceptibility}, we see that both approximations coincide in the case of the susceptibility for small wave numbers. All protocols are thus consistent with our main assumption, i.e., the marginal Fermi liquid behavior is caused by small momentum transfer along the quasi-one dimensional segments of the Fermi line. The same is true for the corresponding relaxation times.

\begin{figure}[t]
\includegraphics[width=0.99\columnwidth]{fig8}
\caption{Full susceptibility of the second highest VB (black line) for two different chemical potentials relative to the vHS $\Delta\mu=0.2$meV (left) and $\Delta\mu=1.5$meV (right). Also shown the susceptibilities that include approximate treatments of umklapp processes.
\label{Susceptibility}}
\end{figure}

\section{Relaxation time for effective dielectric media}
In the main text, we have assumed a strongly screened Hubbard interaction $U$, valid for gates close to the twisted bilayer sample. Here, we will outline the formalism including screening effects within the $G_0W$-approximation. We will first discuss the case of a dielectric function due to long-ranged Coulomb interaction and then also estimate the effect of localised plasmonic modes predicted in TBG.\cite{Stauber16}

\subsection{Relaxation time for long-ranged interaction.}
For gates further away, we expect also effects from the long-ranged Coulomb potential to become important. In this case, we calculate the scattering rate by incorporating the intrinsic screening effects within the $G_0W$-approximation of the self-energy, starting from the Coulomb potential, screened by a bottom and top gate at distance $D$:\cite{Cea19}
\begin{align}
 v_q=\frac{e^2}{2\epsilon_0\epsilon q}\frac{1-e^{-qD}}{1+e^{-qD}}\;.
\end{align} 
We will set $\epsilon=5$, the approximate value for hBN. 

The relaxation time within the $G_0W$-approximation at finite temperature for a quasiparticle (hole) state with $\Delta=E_{\bf p}\mu$ is given by\cite{Giuliani82}   
\begin{align}
\frac{1}{\tau(\Delta)}&=\int_{-\infty}^{\infty}\frac{d\omega}{2\pi}f(\omega)\frac{1}{A}\sum_{\bm{q}} v_q|\langle {\bf p}|{\bf p}+{\bm{q}}\rangle|^2\Im\left(\frac{1}{\epsilon({\bm{q}},\omega)}\right)\delta\left(\hbar\omega-(E_{\bf p}-E_{{\bf p}+{\bm{q}}})\right)\notag
\end{align}
which involves the dielectric function within the RPA
\begin{align}
\epsilon({\bm{q}},\omega)=1-v_q\chi({\bm{q}},\omega)
\end{align}
with the polarisability $(g_s=g_v=2)$
\begin{align}
\chi({\bm{q}},\omega)=\frac{g_sg_v}{A}\sum_\k|\langle \k|\k+{\bm{q}}\rangle|^2\frac{n_F(E_{\k})-n_F(E_{\k+{\bm{q}}})}{\hbar\omega-(E_{\k+{\bm{q}}}-E_{\k})+i0}\;,
\end{align}
and the temperature-dependent weight factor
\begin{align}
\label{f}
f(\omega)=\frac{\coth(\beta\hbar\omega/2)-\tanh(\beta(\hbar\omega-\Delta)/2)}{1+e^{-\beta\Delta}}\;.
\end{align}
We further defined the eigenstates $|{\bf p}\rangle$, the Fermi function $n_F(E)=(e^{\beta E}+1)^{-1}$, the inverse temperature $\beta=1/(k_BT)$. 

The $G_0W$-approximation requires the knowledge of the real and imaginary parts of the susceptibility and, in order to speed up the calculations, we have worked with the less time-consuming continuum model.\cite{Lopes07} In Fig. \ref{ScatteringRateLongRange}, we show the resulting scattering rate for  different chemical potentials relative to the vHS as function of the temperature. Interestingly, for larger gate distances $D\sim15$nm the results are close to the Planckian scattering rate $\hbar/\tau=0.086meV\cdot T[K]$ indicated as dashed line which is in good agreement to the experimental findings of Ref. \onlinecite{Cao19}. 



\begin{figure}[t]
\includegraphics[width=0.99\columnwidth]{fig9}
\caption{The scattering rate $\hbar/\tau$ of TBG with $i=29$ as function of the temperature for two chemical potentials around the vHS and screened long-ranged interaction with surrounding dielectric material $\epsilon=5$. $D$ denotes the distance of TBG to the top and bottom gate. The dashed line indicates the Planckian scattering rate $\hbar/\tau=0.086meV\cdot T[K]$.
\label{ScatteringRateLongRange}}
\end{figure}




What is seen independent of the gate distance $D$ is that for a chemical potential close to the van Hove singularity there is a linear low-temperature behaviour ($\Delta\mu\llless0.2$) in contrary to the quadratic low-temperature behaviour for $\Delta\mu\gsim1.5$meV. Also seen for all curves is the crossover to a different quasi-linear temperature regime for $T_{cr}\sim5$K.

\subsection{Relaxation time from collective modes}
For temperatures larger than the band-gap, the system is expected to reach the classical regime, characterised by a linear behaviour of the resistivity and dominated by the thermal charge fluctuations. In Fig. \ref{LossFunction}, we show the loss function $S(\omega)=$-Im$\epsilon^{-1}({\bm{q}},\omega)$ for $|{\bm{q}}|a=0.02$ in the $KK'$-direction for $\epsilon_0=4.8$ for various twist angles. The peak resembles a true plasmonic resonance as discussed in Ref. \cite{Stauber16} which shifts to smaller energies with decreasing twist angles. 

For $i=29$, also a Lorentzian fit is shown with 
\begin{align}
\tilde S(\omega)=\frac{2C\gamma}{(\omega-\omega_0)^2+\gamma^2}\to2\pi C\delta(\omega-\omega_0)\;,
\end{align}
where $\hbar C=3$meV, $\hbar\gamma=4$meV, and $\hbar\omega_0\sim8$meV. This yields the following scattering rate for $U=5meVa_M^2$ and $\mu$ close to the van Hove singularity, i.e., $U\rho(\mu)=4$:
\begin{align}
\frac{1}{\tau}=0.16\frac{k_BT}{\hbar}
\end{align}

With the plasmon energy $\hbar\omega_0/t\approx0.003$, the crossover temperature corresponds to $100$K which is clearly too high to explain the experiments of \cite{Cao19}. But this limit is imposed by the accuracy of our numerical solution and we expect a linear behaviour for the resonance as indicated by the red line in the inset. In fact, the plasmonic resonance is related to the band-width of the lowest valence/conduction bands which is around 1meV resp. 10K. 

The quasi-particle relaxation time will be mainly determined by the specific form of the loss function which was discussed in Ref. \cite{Stauber16} for small angle twisted bilayer graphene. There, a quasi-flat mode was found for small twist angles, independent of moderate doping-levels related to the localised states around the $AA$-region. As a first approach, we can thus approximate the loss function by the following analytical function:
\begin{align}
\Im\left(\frac{1}{\epsilon({\bm{q}},\omega)}\right)=2\pi C\delta(\omega-\omega_0)
\end{align}
with some suitable constant $C$ which permits for an analytical solution of the relaxation time neglecting the wave function overlap $|\langle {\bf p}|{\bf p}+{\bm{q}}\rangle|^2$. We get
\begin{align}
\frac{1}{\tau(\Delta)}=UCf(\omega_0)\rho(E_{\bf p}-\hbar\omega_0)
\end{align}
where $\rho(\omega)$ denotes the density of states. For large temperature, $\omega_0\ll k_BT$

We expect that the decay rate is dominated by the electron-plasmon coupling active for quasi-particle energies $\Delta\approx\omega_0$. Then, for sufficiently large temperatures $k_BT\gg\hbar\omega_0$, we obtain
\begin{align}
\frac{1}{\tau(\omega_0)}&\approx UC\rho(\mu)\frac{k_BT}{\hbar\omega_0}
\end{align}
Our approach thus yields the observed linear $T$-resistivity above some energy scale $\hbar\omega$. Furthermore, the prefactor is governed by the density of states which is decreasing as one approaches  the regime of half-filling from below. 

\begin{figure}[t]
\includegraphics[width=0.99\columnwidth]{fig10}
\caption{Loss function $S(\omega)=$-Im$\epsilon^{-1}({\bm{q}},\omega)$ for $|{\bm{q}}|a=0.02$ in the $KK'$-direction for $\epsilon_0=4.8$ for twist angles with $i=20-29$ (left). For $i=29$, also a Lorentzian fit is shown with $\tilde S(\omega)=2C\gamma[(\omega-\omega_0)^2+\gamma^2]^{-1}$ where $\hbar C=3$meV, $\hbar\gamma=4$meV, and $\hbar\omega_0\sim8$meV.
\label{LossFunction}}
\end{figure}





\section{Entropy of the electron liquid}



The electronic contribution to the entropy $S$ can be obtained by applying the formula\cite{Abrikosov63}
\begin{align}
\frac{S}{A} = \frac{i}{\pi}\frac{1}{T} \int \frac{d^2 k}{(2\pi)^2} \int_{-\infty }^{\infty } d\omega \: \omega \frac{\partial n_F(\omega)}{\partial \omega } \log \left( \frac{G_R (\k, \omega )}{G_A (\k, \omega )} \right)\;,
\end{align}
where $A$ is the area of the system and $G_R , G_A$ are the retarded and advanced electron Green's functions, respectively. When looking for the low-temperature dependence of the entropy, one can perform the momentum integral along the Fermi line by decomposing $\k $ into longitudinal $k_{\parallel }$ and transverse $k_{\perp }$ components. This leads to 
\begin{align}
\frac{S}{A}  \approx   \frac{1}{\pi^2 T} \oint  \frac{dk_{\parallel }}{2\pi } \int \frac{d\varepsilon_{\k}}{v_{\k}} \int_{-\infty }^{\infty } d\omega \: \omega \frac{\partial n_F(\omega)}{\partial \omega }    \arctan \left( \frac{{\rm Im} \: \Sigma (\k, \omega )}{\omega - {\rm Re} \: \Sigma (\k, \omega ) - \varepsilon_{\k}}  \right) 
\end{align} 
where the integral in $k_{\parallel }$ is carried out along the Fermi line. 
The integral over the energy variable $\varepsilon_{\k}$ can be computed by adopting a principal value prescription. Then we get
\begin{align}
\frac{S}{A}  \approx \frac{1}{2 \pi^2} \frac{1}{T} \oint  \frac{dk_{\parallel }}{v_{\k}}  \int_{-\infty }^{\infty } d\omega \: \omega \frac{\partial n_F(\omega)}{\partial \omega }  \left( \omega - {\rm Re} \: \Sigma (\k, \omega ) \right)
\label{entr}
\end{align}

The temperature dependence of the entropy can be estimated by absorbing $T$ into a dimensionless variable $\omega/T$ in the integrand of Eq. (\ref{entr}). In particular, when the real part of the electron self-energy has an anomalous logarithmic correction, we see that this is translated to the entropy, which gets a dominant scaling behavior $S \sim T \: |\log(T)|$.

\end{widetext}
