text:
\section{Introduction}

Discovering disease subgroups that share distinct disease mechanisms is fundamental to disease prevention, monitoring and treatment. For example, in autoimmune diseases, {specific autoimmune responses are associated with distinct disease phenotypes and trajectories} \citep{rosen2016autoantigens}. {Defining the molecular markers of these subgroups has value, as these markers are of diagnostic and prognostic significance, and guide management and therapy. For example, an immune response to RNA polymerase III in scleroderma is associated with cancer; this immune response arises in response to a mutation in RNA polymerase III in that patient's cancer.  While many prominent specificities recognized by the immune response have been defined, many remain to be discovered. Although modern measurement technologies are revolutionizing the ability to define specificities, each technique results in hundreds of possibilities, which are onerous and expensive to validate. A simple technique identifies patterns of antibody reactivity based on} the abundance of different weighted antigens immunoprecipitated by patient sera. {Defining similar antibody reactivity patterns prior to applying one of the new discovery technologies would greatly simplify validation and therefore reduce the cost and improve the speed of antigen identification.}

To identify the autoantibodies present in a patient's serum, scientists {mix serum collected from each patient with radiolabeled lysates made from cultured cells.  These lysates contain a representation of all the proteins expressed in that cell type.  Antibodies in each patient's serum recognize and bind tightly to the specific protein(s) in the lysate against which they are directed (termed "immunoprecipitation"). After further processing, electrophoresis is used to sort the immunoprecipitated} mixture of molecules using a crosslinked polymer or gel that separates the proteins by weight. Because different weighted molecules migrate with different speeds, the sorted molecules form distinct autoradiographed bands along the gel. By design, one gel can sort multiple samples on parallel lanes. Such experiments, referred to as gel electrophoresis autoradiography (GEA), serve to identify subsets that share one or more interesting observed bands. {It is noteworthy that the lysate proteins are present in their native conformation.  In our experience, many autoantibodies have epitopes that are conformationally dependent, giving GEA a powerful advantage over many of the new peptide-based (linear epitopes) sequencing technologies.} The method in this paper is designed to estimate a multivariate binary autoantibody signature for each sample that represents the presence or absence of autoantibodies over a grid of molecular weights, referred to as \textit{landmarks}.

To infer patient subsets, we  can cluster patients based upon the presence or absence of each band as well as other features of the radioactive intensities such as the peak scale and amplitude. There are two critical barriers to the successful implementation of this approach that we address. First, there are \textit{batch}, or \textit{gel effects} in the raw GEA data. By design, molecules of identical weight would migrate the same distance along the gel. This distance however varies by gel due to differential experimental conditions. Second, {gels are frequently slightly warped  as they electrophorese due to heating effects generated during the electrophoresis procedure and due to artifacts introduced during physical processing  of the gels.} As the size and complexity of GEA experiment database grow, the need for systematic, reproducible and scalable error correction has also grown.

In this paper, we introduce and illustrate a novel statistical approach to pre-process the high-frequency GEA data which we show improve our ability to compare and cluster band patterns across samples. The pre-processing involves peak detection and batch effect corrections.  In particular, we propose a local scoring algorithm for peak detection that is computationally efficient and performs well for minor peaks (Section \ref{sec:peakdetection}). The detected peaks then enter the image alignment method that corrects batch effects in two steps: reference alignment and spatial dewarping. First, reference alignment calibrates multiple gels towards a common molecular standard. We perform piecewise linear stretching/compression by placing knots at the \textit{marker} or \textit{reference} bands present on all the images (Section \ref{sec:reference.alignment}). The reference-aligned gel images produce a set of peak locations that are then fitted by a novel hierarchical Bayesian model. The proposed model assumes that the smooth spatial gel deformations have deviated the observed peaks from their true landmarks. We use Markov chain Monte Carlo to estimate both the smooth warping functions and, for each detected peak, the posterior probabilities over a grid of landmarks where it is aligned. The Bayesian framework has the advantage of incorporating inherent uncertainty in assigning a peak to a molecular weight landmark.

 { The \textit{aligned} high-frequency intensity data (Section \ref{sec::exact_matching}) may be the input of many methods including hierarchical clustering, latent class models and factor analyses. In this paper, we focus on illustrating the value of alignment for the standard hierarchical clustering applied to data with known and unknown clusters (Section \ref{sec:application}).}   At each iteration of the MCMC sampling, we obtain the multivariate binary signatures that represent autoantibody presence or absence over a grid of landmarks and align the gel images.  Upon hierarchically clustering the aligned intensities at each iteration, we obtain a collection of dendrograms. In particular, we use the standard correlation-distance based agglomerative hierarchical clustering to create nested subgroups. For $N$ samples, hierarchical clustering produces a dendrogram that represents a nested set of clusters. Depending on where the dendrogram is cut, between 1 and $N$ clusters result. We demonstrate through real data that pre-processing more clearly separates the estimated clusters and improves the accuracy of cluster detection compared to naive analyses done without alignment.

The rest of the paper is organized as follows. Section \ref{sec:preprocessing} introduces the importance of pre-processing GEA data followed by algorithmic details for peak detection in Section \ref{sec:peakdetection} and batch effect correction in Section \ref{sec:batch.effects}.  In Section \ref{sec:computing}, we describe model posterior inference by MCMC and the statistical property of the shrinkage priors. We demonstrate how the proposed methods function through an application to signature estimation and subgroup identification of scleroderma patients in Section \ref{sec:application}. The paper concludes with a discussion on model advantages and opportunities for extensions.

\section{Data Pre-Processing}
\label{sec:preprocessing}

\subsection{GEA Data and Pre-processing Overview}
\label{sec:gelexperiments}
Gel electophoresis for autoantibodies (GEA) is designed to separate {autoantibody} mixtures according to molecular weight and to radioactively map them as bands along the gel. Figure \ref{fig:raw_image} shows an example of a raw GEA image. We tested four sets of samples from scleroderma patients {with a malignancy; of note, these sera were pre-selected as being negative for the three most commonly found scleroderma autoantibodies (anti-topoisomerase 1, anti-centromere and anti-RNA polymerase III antibodies, which in aggregate are found in about $60\%$ of scleroderma patients).}

Each sample set consisted of $19$ patient sera plus one reference. In the middle panel of Figure \ref{fig:raw_data}, seven red vertical lines indicate the reference molecules of known weight (200, 116, 97, 66, 45, 31, 21.5) kDa. It also shows the band patterns that read out autoantibodies present in each of $19$ patient samples (lanes 2-20). The top of Figure \ref{fig:raw_data} shows the intensities from all the lanes; Seven clear spikes above the vertical lines again correspond to the reference molecules. The bottom of Figure \ref{fig:raw_data} shows the piecewise linear interpolation of the location-to-weight function using the seven reference weights as knots. The weight of an arbitrary peak (``*") can then be read from this interpolation. For example, the protein actin  (about 42kDa) produced the peaks immediately to the right of the 45kDa reference (the fifth vertical line from the left). The misalignment of the actin peaks is caused by non-rigid image deformation (Section \ref{sec:bayesian.warping}).

Identical reference molecules fail to align (empty circles, bottom panel of Figure \ref{fig:raw_data}) across multiple gels because of variation in experimental conditions such as the strength of the electric field. We correct such misalignment by matching the reference peak locations across gels and then piecewise-linearly stretch or compress each gel using the reference peaks as knots. The technique is referred to as \textit{piecewise linear dewarping} and was first used in human motion alignment anchored at body joints \citep[e.g.,][]{uchida2001piecewise}.

The autoradiographic process is also vulnerable to smooth non-rigid gel deformation. This is most evident from the bands of actin, a ubiquitous protein of molecular weight 42 kDa, present in all lanes at around $0.43$ (middle panel of Figure \ref{fig:raw_data}). The bands form a smooth curve from the top to the bottom. The curvature represents the gel deformation since actin has identical weight and should appear at identical locations across the $19$ lanes. Without correction, this deformation interferes with accurate sample comparisons even on the same gel.  In Section \ref{sec:bayesian.warping}, we propose a Bayesian hierarchical image dewarping model with shrinkage priors to correct the deformation and align the actin peaks.

To establish notations, let $(\mathbf{t}^{0}, \mathbf{M}^{0})=\left\{\left(t^{0}_{b}, M^{0}_{gib}\right)\right\}$ represent the standardized, high-frequency GEA data, for bin $b=1, \ldots, B$ on lane $i=1,\ldots, N_g$ from gel $g=1, \ldots, G$.  Appendix S1 describes the standardization of raw data. Here $\mathbf{t}^{0}$ is a equi-spaced grid over the unit interval $[0,1]$, where $t^{0}_{gb}=b/B\in [0,1]$, $b=1,\ldots, B$. $M_{gib}^{0}$ is the radioactive intensity scanned at $t^{0}_{b}$ for lane $i=1,\ldots, N_g$, gel $g=1,\ldots, G$. Let $N = \sum_{g}N_g$ be the total sample size.

For the rest of this section, we take the high-frequency data $(\mathbf{t}^{0}, \mathbf{M}^{0})$ and map it to multivariate binary data $\mathbf{Y}$ on a coarser common grid across gels. In Section \ref{sec:peakdetection}, we propose a general method to transform an arbitrary high frequency, nearly continuous intensity data into raw peak locations. We first apply the peak detection algorithm to $(\mathbf{t}^{0}, \mathbf{M}^{0})$ and obtain the peak locations $\mathcal{P}^{0}$. In Section \ref{sec:reference.alignment} we use the reference peaks, a subset in $\mathcal{P}^{0}$ from the first lane on each gel, to process $(\mathbf{t}^{0}, \mathbf{M}^{0})$ into reference-aligned data $(\mathbf{t}^{\sf R}, \mathbf{M}^{\sf R})$. In Section \ref{sec:bayesian.warping}, we transform the peaks $\mathcal{P}$ detected from $(\mathbf{t}^{\sf R}, \mathbf{M}^{\sf R})$ to a joint posterior distribution of a $N$ by $L$ binary matrix $\mathbf{Y}$ that represents presence or absence of a peak over a grid of $L$ landmarks for all the $N$ samples ($L=100$ in our application). In Section \ref{sec::exact_matching}, we will process the reference-aligned high-frequency data $(\mathbf{t}^{\sf R}, \mathbf{M}^{\sf R})$ into $(\mathbf{t}, \mathbf{M})$ where the peaks appear at the landmarks indicated by the ones in $\mathbf{Y}$.

\subsection{Peak Detection}
\label{sec:peakdetection}

This section presents a general algorithm for detecting peaks from intensity data. We illustrate the algorithm by detecting peaks $\mathcal{P}^{0}$  from data $(\mathbf{t}^{0}, \mathbf{M}^{0})$. The peaks may appear with varying background intensities. Because the occurrence of a local maximum is thought to be more important than the background level in autoantibody signature estimation, we design the algorithm to be insensitive to the absolute intensity level.

We adopted the following peak detection algorithm:

\begin{enumerate}

\item[i.] \textit{Local Difference Scoring}. For each bin $b=1,\ldots, B$, lane $i=1,\ldots, N_g$ of gel $g=1,\ldots, G$, calculate the local difference score by comparing the intensity at bin $b$ to its left and right neighbors exactly $h$ bins away and to the local minimum for locations in between $(t-h,t+h)$ (truncated at $1$ or $B$ if $b$ is near the endpoints). That is, we calculate
\begin{align}
\textsf{score}_{gi}(b)  & =
{\sf sign}\left\{M^{0}_{gib}-M^{0}_{gi,\ell(b)}\right\}+
{\sf sign}\left\{M^{0}_{gib}-M^{0}_{gi,r(b)}\right\}+\nonumber\\
& ~~~~ {\sf sign}\left\{M^{0}_{gib}-\min_{\ell(b) \leq b' \leq r(b)}{M^{0}_{gib'}}-C_0\right\},\label{eq:localscoring}
\end{align}
where ${\sf sign(a)}=1,0,-1$ indicates positive, zero, or negative values; $\ell(b)=\max\{b-h,1\}$ and $r(b)=\min\{b+h,B\}$ denote the left and right neighbors $h(=10)$ bins away, and $C_0$ denotes the minimum peak elevation. The tunning parameter $h$ controls the locality of the peaks and $C_0$ controls the minimum peak magnitude.

\item[ii.]  \textit{Peak Calling}. We look for the bins among peak candidates defined by $\{b\mid {\sf score}_{gi}(b)=3\}$ that maximize their respective local intensities (see  Appendix S2 for details and alternative peak calling methods). Let $\mathcal{P}^{0}_{gi}$ represent the collection of the peak locations for lane $i$ and gel $g$.

\end{enumerate}

\noindent \textbf{Remark 1:}  The score defined in (\ref{eq:localscoring}) depends only on the \textit{sign}s of the differences in local intensities. They can be computed in parallel across all the samples. A two-dimensional analogue has been used in astrophysics to find low grey-scale intensity galaxies from telescope images \citep{xu2016detection}.

\subsection{Batch Effect Correction}
\label{sec:batch.effects}

\subsubsection{Reference Alignment via Piecewise Linear Dewarping.}
\label{sec:reference.alignment}

Molecules with identical weight do not appear exactly at the same location in each lane of a single gel due to gel deformation or across gels due to variations in experimental conditions. We first align the reference peaks $\mathcal{P}^{0}_{g1}$, $g=1, \ldots, G$ via piecewise linear dewarping to address the gel-to-gel variation \citep{uchida2001piecewise}. In our application, we used seven reference molecules of known weight $(200, 116, 97, 66, 45, 31, 21.5)$ kDa.

We first match the reference peaks $\mathcal{P}^{0}_{g1}$ on a query gel $g$ to the reference peaks $\mathcal{P}_{g_01}$ on the template gel $g_0$, and then use the matched reference peaks and the endpoints as knots to linearly stretch or compress the gels. Quadratic or higher-order dewarping is also possible, but we found linear dewarping performs sufficiently well for our data.  Appendix S3 gives the details of the algorithm. We denote the high frequency, reference aligned data by $(\mathbf{t}^{\sf R}, \mathbf{M}^{\sf R})=\{(t^{\sf R}_{gib}, M^{\sf R}_{gib})\}$. Applying the peak detection algorithm in Section \ref{sec:peakdetection} to this data, we collect all the detected peaks in $\mathcal{P}=\{\mathcal{P}_g, g=1,\ldots, G\}$ where $\mathcal{P}_g$ represents the peaks from gel $g$.

\subsubsection{Bayesian Image Dewarping to Correct Gel Deformation.}
\label{sec:bayesian.warping}

Another source of error during autoradiographic visualization is the non-rigid, spatial gel deformation. The middle panel of Figure \ref{fig:raw_data} shows one such example. It also reveals three analytical challenges to be addressed before obtaining meaningful results from an automatic disease subsetting algorithm. First, some proteins, e.g., actin, are detected on multiple gels and must be aligned. The blue asterisks that denote the detected peaks near $0.43$, form a smooth but non-linear curve from the top to the bottom of the gel. Second, fewer bands appear on the right half of the image, because {these smaller proteins tend to contain fewer methionine residues for radiolabeling}. Higher estimation uncertainty of the dewarping function is therefore expected for the right half. Third, the observed locations of the peaks are likely random around their true locations as the result of the multiple sources of error.

To address these issues, we designed a hierarchical Bayesian dewarping algorithm for two-dimensional images. The algorithm simulates presence/absence data from the conditional distribution of protein occurrence over a grid of equi-spaced landmarks given the detected peaks $\mathcal{P}$ from the prior pre-processing. The stochastic model is defined on a coarser grid of landmark proteins, $\bm{\nu}=\{0=\nu_0 < \nu_1 < \ldots < \nu_L < \nu_{L+1}=1\}$ where $\nu_\ell = \ell/(L+1)$, $\ell=0, 1, \ldots, L+1$. In this paper, we will align peaks only to the internal knots $\{\nu_{\ell}, \ell\neq 0,L+1\}$; $\nu_0$ and $\nu_{L+1}$ will be used in the boundary constraint (\ref{eq:boundary_constraints}) to ensure endpoint alignment for all the sample lanes. We introduce a novel shrinkage prior to promote alignment of peaks to a common landmark. We also introduce shrinkage priors that regularize the overall smoothness of the spatial dewarping functions.

Let $(T_{gij}, u_{gij})$ denote the ({\sf location, lane number}) for peak $j=1, \ldots, J_{gi}$ on lane $i=1,\ldots, N_g$, gel $g=1, \ldots, G$. We fix $u_{gij}$ to take values in $\{1,2,\ldots, N_g\}$ and collect them in $\bm{u} = \{u_{gij}\}$ where $u_{gij}=u_{gi}$ if they belong to the same lane $i$. Let $P_g=\sum_i {J_{gi}}$ denote the total number of peaks on gel $g$ and $P = \sum P_g$. Let $\mathbf{T} = \{\bm{T}_g\}$, where $\bm{T}_g = (\ldots, T_{gi1},T_{gi2}, \ldots, T_{gi,{J_{gi}}}, \ldots)'$ collects the peak locations for gel $g=1, \ldots, G$. Both $\bm{u}$ and $\mathbf{T}$ are $P$-dimensional column vectors. For computational stability, without changing notation, we standardize $\mathbf{T}$, $\bm{u}$ and $\bm{\nu}$ by substract their means and dividing by their standard deviations. We now use $\mathcal{P}=\{\mathbf{T}, \bm{u}\}$ to denote the data for the Bayesian dewarping model.

\noindent \underline{Model Likelihood.} \textit{Peak-to-landmark indicators $\mathbf{Z}$.} Let $Z_{gij}$ take values in $\{1,\ldots,L\}$. For example, $Z_{gij}=3$ indicates that the $j$-th peak in lane $i$ on gel $g$ is aligned to landmark 3. Let $\mathbf{Z}=\{\bm{Z}_{g}, g=1, \ldots, G\}$ where $\bm{Z}_g=\{Z_{gij},j=1, \ldots, J_{gi}, i=1, \ldots, N_g\}$. Note that any $\mathbf{Z}$ can be converted to $N$ multivariate binary observations $\mathbf{Y} = \{(Y_{gi\ell}, \ell=1, \ldots, L)\}$ for the presence or absence of a landmark, where \({Y}_{gi\ell}=\mathbf{1}\left\{\ell \in \{{Z}_{gij}, j=1, \ldots, J_{gi}\}\right\}\), referred to as \textit{signature}.

\vspace{1em}
\noindent \textit{Gaussian mixture model for aligning observed peaks $\mathbf{T}$.} We model $\mathbf{T}$ as observations from a Gaussian mixture model with $L$ components, each representing one landmark. Given $\mathbf{Z} = \{Z_{gij}\}$ and the spatial dewarping function $\mathcal{S}_g$ to be discussed later, we assume
\begin{align}
p\left\{(\underbrace{T_{gij}=t}_{\substack{\text{peak}\\ \text{location} }}, \underbrace{u_{gi}}_{\substack{\text{lane}\\ \text{number} }})\mid \underbrace{Z_{gij}=\ell}_{\substack{\text{matched~to} \\
\text{landmark~}\ell}}, \underbrace{T_{gi,j-1}}_{\substack{\text{nearest~left} \\
\text{peak~location}}},
\underbrace{\mathcal{S}_g}_{\substack{\text{warping} \\
\text{function}}},
\underbrace{\sigma_{\epsilon}}_{\substack{\text{noise}\\ \text{level}}} \right\} = \begin{cases}
\phi\left(t; \mathcal{S}_g(\nu_\ell, u_{gi}), \sigma_{\epsilon}\right),&~t\in\mathcal{I}_{gij}(\nu_\ell, A_0);\\
0,&~\text{otherwise},
\end{cases} \label{eq:component_lkd}
\end{align}
$\ell=1,\ldots,L$, for peak $j=1, \ldots, J_{gi}$, lane $i=1, \ldots, N_g$, gel $g=1, \ldots, G$, where $\phi(\cdot; a,b)$ is the Gaussian density function with mean $a$ and standard deviation $b$, and $\mathcal{S}_g$ is an unknown smooth bivariate function that characterizes the deformation $(\nu, u) \mapsto \left(\mathcal{S}_g(\nu, u),u\right)$.

\vspace{1em}
\noindent \textbf{{{Remark 2:}}} The peak location $T_{gij}$ is assumed to follow a Gaussian distribution with mean equal to $\nu_\ell$ plus a horizontal displacement $\mathcal{S}_g(\nu_\ell, u_{gi})$ and noise level equal to $\sigma_\epsilon$. We assume $\sigma_{\epsilon}$ is independent of landmark and lane. The density function (\ref{eq:component_lkd}) is positive only in the set $\mathcal{I}_{gij}(\nu_\ell, A_0)\overset{\Delta}{=}\left\{t: |t-\nu_\ell|< A_0 \text{~and~} t > T_{gi,j-1}\right\}$. The first inequality prohibits $T_{gij} $ being matched to distant landmarks and limits the search space for $Z_{gij}$ in our algorithm; the second inequality places order constraints on the observed peak locations $T_{gij}>T_{gi,j-1}, j=2, \ldots, J_{gi}-1$. We will restrict $Z_{gij}>Z_{gi,j-1}$ to avoid reverse dewarping.

\vspace{1em}
\noindent \textit{Bivariate smooth warping functions $\mathcal{S}_g$.} For gel $g$, we model the warping function $\mathcal{S}_g:\mathbb{R}^2\rightarrow\mathbb{R}$ using the tensor product basis expansion
\begin{eqnarray}
\mathcal{S}_g(\nu, u) & = &\sum_{s=1}^{T_\nu}\sum_{t=1}^{T_{u}}\beta_{gst}B_{g1s}(\nu)B_{g2t}(u),\label{eq:tensor_product}
\end{eqnarray}
where $B_{g1s}(\cdot)$ and $B_{g2t}(\cdot)$ are the $s$-th and $t$-th cubic B-spline basis with intercept, and $\bm{\kappa}_{\nu}$ and $\bm{\kappa}_{u}$ are the knots along the two coordinate directions, respectively \citep[][Chapter 5]{friedman2001elements} and $T_{\nu}$ and $T_{u}$ are the total number of bases. In subsequent analyses, we choose $\bm{\kappa}_\nu$ with $T_\nu-4$ internal knots at the $s/(T_\nu-3)$-th quantile of $\{T_{gij}\}$, $s=1, \ldots, T_\nu-4$ and similarly for $\bm{\kappa}_u$. Let the two sets of B-spline basis functions along  $\nu$- and $u-$direction be $\bm{B}_{g1}(\cdot)=\left(B_{g11}(\cdot), \ldots, B_{g1T_\nu}(\cdot)\right)'$ and $\bm{B}_{g2}(\cdot)=(B_{g21}(\cdot), \ldots, B_{g2T_u}(\cdot))'$, respectively.

However, valid spatial gel deformations are limited to gel stretching, compression or shift along the $\nu$ direction. We thus constrain the shape of $\mathcal{S}_g, g=1, \ldots, G$ by
\vspace{-0.5em}
\begin{eqnarray}
\textit{Monotonicity}: &&\nu_{0}\leq \mathcal{S}_g(\nu, u)< \mathcal{S}_g(\nu', u)\leq \nu_{L+1}, \forall \nu< \nu', \forall u; \label{eq:monotonicity}\\
\textit{Boundary Constraint}:&& \mathcal{S}_g(\nu_0, u)=\nu_0, \mathcal{S}_g(\nu_{L+1}, u)=\nu_{L+1}. \label{eq:boundary_constraints}
\end{eqnarray}

\vspace{-0.5em}
The first constraint prevents reverse gel dewarping and the second assumes no gel shifting.  It can be relaxed to allow horizontal shifts by adding/substracting $\Delta$ for both equalities.  We implement both constraints by requiring the B-spline coefficients $\bm{\beta}_g=\{\beta_{gst}\}$ to satisfy: $\nu_{0}= \beta_{g1 t}< \beta_{g2t}<\ldots<\beta_{g,T_{\nu}-1,t}< \beta_{gT_\nu t} = \nu_{L+1}$, $\forall t=1,\ldots,T_u$. Although only sufficient for $\mathcal{S}_g$'s monotonicity and boundary constraints, the foregoing $\bm{\beta}_g$ constraints allow flexible and realistic warpings. Figure \ref{fig:example_warping} shows a member warping function that corrects for local $``L"$-, $``S"$- and $``7"$-shaped deformations.

The likelihood function (\ref{eq:component_lkd}) models the misaligned data $\mathcal{P}=\{\mathbf{T},\bm{u}\}$ in terms of the unknown spatial transformation $\mathcal{S}_g$ and the alignment $\mathbf{Z}$. Multiple raw gel images are then aligned by the model estimates accompanied by model-based uncertainty quantification. Importantly, coherent image registrations must align the universal actin peaks and hence require the borrowing of information among multiple misaligned observations. We accomplish this by sharing a set of intensity parameters $\{\lambda_\ell\}$ among the gels.

\noindent \underline{Prior}. \textit{Prior for} $\mathbf{Z}$. We describe a shrinkage prior for $\mathbf{Z}$ motivated by the need 1) to align the actin peaks (middle panel, Figure \ref{fig:raw_data}), and 2) to share the information about the location of actin peaks across multiple gels.

We specify the prior distribution based on a discretized, non-homogeneous Poisson process with extreme intensities at a small number of landmarks.  Let the total number of the observed peaks follow a Poisson distribution: $J_{gi}  \overset{d}{\sim} {\sf Poisson}(\Lambda_g)$, for sample $i=1, \ldots, N_g$, gel $g=1, \ldots, G$. Given $J_{gi}$, let $Z^*_{gij} \mid \{\lambda_\ell^*\} \overset{iid}{\sim} {\sf Categorical}\left(\left\{ \lambda^*_\ell\right\}_{\ell=1}^L\right)$ describe which landmarks are present in sample $i$ of gel $g$. For sample lane $i$, we then define $\{Z_{gij}\}$ as the increasingly sorted $\{Z^*_{gij}\}$.  That is, we impose the order restriction $Z_{gij}\leq Z_{gij'}$ whenever peak $j$ appears to the left of peak $j'$ ($T_{gij}\leq T_{gij'}$). For hyperpriors, let $\lambda^*_\ell={\lambda_{\ell}}/{\sum_{\ell'=1}^L\lambda_{\ell'}}$ where $\lambda_\ell \mid \tau \overset{iid}{\sim} {\sf Normal}(0,\tau),  \ell=1,\ldots, L$, and the hyperparameter $\tau \overset{d}{\sim} {\sf Inv}$-${\sf Gamma}(10^{-4},10^{-4})$. Integrating over $\tau$, we obtain a marginal $t$-distribution for $\lambda_\ell$.

\vspace{0em}
\noindent \textbf{Remark 3:} It is easy to calculate the prior probability of landmark $\ell$ present in a sample $\mathbb{P}(Y_{gi\ell}=1 \mid \lambda^*_\ell) \approx 1-\exp(-\lambda^*_\ell)$, $\ell=1, \ldots, L$, for large $L$. As shown by (A4) in Appendix, the ratio of the conditional posterior probabilities of assigning the peak $T_{gij}$ to landmark $\ell$ versus $\ell'$ is factorized into
\(\frac{\phi(T_{gij}; \mathcal{S}_g(\nu_{\ell}, u_{gi}),\sigma)}{\phi(T_{gij}; \mathcal{S}_g(\nu_{\ell'}, u_{gi}),\sigma)}\cdot\frac{1-\exp(-\lambda^*_\ell)}{1-\exp(-\lambda^*_{\ell'})}.\)
Suppose landmark $\ell$ is associated with a higher intensity, i.e., $\lambda^*_\ell>\lambda^*_{\ell'}$, the second ratio favors landmark $\ell$ given the likelihood ratio in the first term. Because the $\{\lambda^*_\ell\}$ are independent of $g$ and $i$, they globally modulate the probability of a landmark being present in all the gels. In our application, all the landmarks are \textit{in a priori} assumed to be equally likely by specifying independent $t$-distributed priors for the $\lambda_\ell$s. The $t$-distributions are heavy-tailed and can occasionally generate a large value of $\lambda_{\ell_0}$. Given $\lambda_{\ell_0}$, the posterior sampling algorithm will visit and then retain any configuration of $\mathbf{Z}$ that results in a large value of $\sum_{g,i}Y_{gi\ell_0}$ if the configuration substantially increases the joint posterior.

\vspace{1em}
\noindent \noindent \textit{Prior for} $\bm{\beta}_g$.
We incorporate the prior knowledge that large and abrupt image deformations are rare.  We first specify priors for the horizontal basis coefficients $\beta_{gst}$, $s=2, \ldots, T_{\nu}-1$ at the $u$-direction basis $t=1$. We use a first-order random walk prior \citep{lang2004bayesian}
\begin{equation}
\beta_{gst} - \beta^{\sf\footnotesize id}_{s}\overset{d}{\sim} {N}\left(\cdot; \beta_{g,s-1,t}-\beta^{\sf\footnotesize id}_{s-1}, \sigma^{-2}_{g1}\right)\mathbf{1}\left\{\beta_{gst}\in (\beta_{g,s-1,t},\nu_{L+1})\right\}, t=1, s=2, \ldots, T_{\nu}-1,\label{eq:horizontal_randomwalk}
\end{equation}
where $\bm{\beta}^{\sf \footnotesize id}=(\beta^{\sf \footnotesize id}_1, \ldots, \beta^{\sf \footnotesize id}_{T_{\nu}})'$ is the vector of coefficients to represent an identity function $I: \nu\mapsto \nu$ in terms of the bases $\{B_{gs1}(\cdot)\}_{s=1}^{T_\nu}$; The truncation of $\beta_{gst}$ is needed for monotonicity (\ref{eq:monotonicity}). The hyperparameter $\sigma_{g1}^2$ controls the similarity between $\{\beta_{gs1}\}_{s=1}^{T_{\nu}}$ and $\bm{\beta}^{\sf id}$ and hence the similarity between $\mathcal{S}_g(\cdot, u)$ and the identity function $I$; $\sigma^2_{g1}=0$ represents no warping. We refer to $\sigma_{g1}^{-2}$ as the smoothing parameter along the $\nu$-direction.

Next, for any $s=2, \ldots, T_{\nu}-1$, we specify another random walk prior for the vertical basis coefficients
\begin{equation}
\beta_{gst} \overset{d}{\sim} {N}\left(\cdot; \beta_{gs,t-1},\sigma^{-2}_{gs}\right)\mathbf{1}{\{\beta_{gst}\in (\beta_{g,s-1,t},\nu_{L+1})\}}, t=2, \ldots, T_u.\label{eq:vertical_randomwalk}
\end{equation}
Similarly, the hyperparameter $\sigma^{-2}_{gs}$ controls the smoothness of $\mathcal{S}_g$ along the vertical or $u-$ direction; $\sigma^2_{gs}=0$ produces identical amounts of warping for all the lanes. Details about the hyperpriors for $\{\sigma_{gs}^2,s=1,2,\ldots, T_\nu-1\}$ are provided in  Appendix S4.

\label{sec:joint_distribution}
{\noindent {\it Joint Distribution.} The joint distribution of all the unknowns is
\begin{eqnarray}
&  & \prod_{g=1}^{G} \biggr\{ \underbrace{\prod_{i=1}^{N_g}\biggr[\prod_{j=1}^{J_{gi}} N\left(T_{gij}; \bm{B}_{g1}(\nu_{Z_{gij}})'\bm{\beta}_{g}\bm{B}_{g2}(u_{gi}), \sigma_{\epsilon}^{-2}\right)\mathbf{1}\{T_{gij}\in \mathcal{I}_{gij}(\nu_{Z_{gij}},A_0)\}}_{\sf likelihood~(\ref{eq:component_lkd})} \nonumber\\
&\times & \underbrace{J_{gi}!\prod_{j=1}^{J_{gi}}{\sf Categorical}(Z_{gij}; \bm{\lambda})\mathbf{1}\{Z_{gij}\leq Z_{gi,j+1},j=1,\ldots, J_{gi}-1\}\biggr]}_{\sf prior~of~\mathbf{Z}}
\nonumber \\
& \times & \underbrace{{N}_{T_{\nu}-1}\left(\{\beta_{gs1}\}_{s=1}^{T_{\nu}-1}; \bm{\beta}^{\sf id}_{[-T_{\nu}]}, \sigma_{g1}^{-2}\Delta_1'\Delta_1\right)\mathbf{1}\{\nu_0=\beta_{g11}<\ldots <\beta_{gs1} <\ldots < \beta_{g,T_\nu-1,1}<\nu_{L+1}\}\cdot p(\sigma_{g1}^2)}_{\sf prior~(\ref{eq:horizontal_randomwalk})~and~hyperprior~of~the~smoothing~parameter} \nonumber \\
& \times &\underbrace{ \prod_{s=2}^{T_{\nu}-1} \left[{N}_{T_u}\left(\bm{\beta}_{gs\bullet}; \mathbf{0}, \sigma_{gs}^{-2}\Delta_2'\Delta_2\right) \mathbf{1}\{\nu_0=\beta_{g1t} < \beta_{g,s-1,t} <\beta_{gst} < \nu_{L+1}, \forall t \geq 2\}
\cdot
p(\sigma^2_{gs}, \rho_{g})
\right]}_{\sf prior~(\ref{eq:vertical_randomwalk})~and~hyperpriors~of~the~smoothing~parameters}\biggr\}  {p(\bm{\lambda})},\label{eq:joint_distn}
\end{eqnarray}
where $p(\bm{\lambda})$, $p(\sigma_{g1}^2)$ and $p(\sigma^2_{gs}, \rho_{g})$ are the priors and hyperpriors and $N_d(\cdot; \bm{\mu}, \bm{\Lambda})$ denotes the $d$-dimensional multivariate normal density with mean $\bm{\mu}$ and precision matrix $\bm{\Lambda}$ (can be degenerate). The matrix $\Delta_1$ maps a column vector to its first-order differences (used in \ref{eq:horizontal_randomwalk}): $\Delta_{1kk'}=\delta(k+1,k')-\delta(k,k')$, $k=1, \ldots, T_{\nu}-2$, $k'=1, \ldots, T_\nu-1$, where $\delta(a,b)=1$ if $a=b$ and equals $0$ otherwise; Similarly we define $\Delta_2$ with $T_\nu$ replaced by $T_u+1$.
}

\section{Model Estimation and Implementation}
\label{sec:computing}

\subsection{Posterior Sampling}
\label{sec::posterior_sampling}
We use Markov chain Monte Carlo (MCMC) to simulate samples from the joint posterior distribution of all the unknowns \citep[e.g.,][]{gelfand1990sampling} and then draw posterior inferences about chosen functionals of the model parameters. Of special interest are the gel warping functions $\{\mathcal{S}_g(\cdot,\cdot; \bm{\beta})\}$ and the peak-to-landmark alignment $\mathbf{Z}$. Appendix S5 describes the sampling algorithm and discusses conditions for statistical identifiability of the warping functions. All the model estimation and visualization are performed by the \texttt{R} package \texttt{spotgear} (\url{https://github.com/zhenkewu/spotgear}).

Turning to dewarping a new GEA image $g^*$, we perform reference alignment (Section \ref{sec:reference.alignment}) and then obtain the peaks $\mathcal{P}_{g^*}$ (Section \ref{sec:peakdetection}). We approximate the joint posterior of $(\bm{\beta}_{g^*}, \bm{Z}_{g^*})$ by
\begin{eqnarray}
p(\bm{\beta}_{g^*}, \bm{Z}_{g^*} \mid \mathcal{P}, \mathcal{P}_{g^*}) & =
& \int p(\bm{\beta}_{g^*}, \bm{Z}_{g^*} \mid \bm{\lambda}, \mathcal{P}_{g^*})p(\bm{\lambda} \mid \mathcal{P}, \mathcal{P}_{g^*}) \mathrm{d} \bm{\lambda}   \approx \int p(\bm{\beta}_{g^*}, \bm{Z}_{g^*} \mid \bm{\lambda},\mathcal{P}_{g^*})p(\bm{\lambda} \mid \mathcal{P}) \mathrm{d} \bm{\lambda},\nonumber
\end{eqnarray}
where the first term of the integrand is an one-sample conditional posterior and the second term is the posterior of $\bm{\lambda}$ given the old peaks $\mathcal{P}$.  Given $\mathcal{P}_{g^*}$, the first term can be derived from the joint distribution (\ref{eq:joint_distn}) with $G=1$. The integral can then be approximated by $K^{-1} \sum_k p(\bm{\beta}_{g^*}, \bm{Z}_{g^*} \mid \bm{\lambda}^{(k)},\mathcal{P}_{g^*})$ where $\{\bm{\lambda}^{(k)},k=1,\ldots, K\}$ are the stored posterior samples.

\subsection{Exact Peak Alignment for $\mathbf{M}$}
\label{sec::exact_matching}

We now describe the high-frequency data $(\mathbf{t}, \mathbf{M})=\left(\mathbf{t}, \mathbf{M}(\mathbf{Z})\right)$ that have peaks \textit{exactly} aligned according to $\mathbf{Z}$. We simply perform piecewise linear dewarping for each sample lane so that the detected peak $(T_{gij}, u_{gi})$ is horizontally adjusted to match its landmark $\nu_{{Z}_{gij}}$, $ j=1, \ldots, J_{gi}$. To do so, we apply the algorithm in Appendix S3 with $\{\nu_0, T_{gi1}, \ldots, T_{giJ_{gi}}, \nu_{L+1}\}$ as a query and $\left\{\nu_0, \nu_{{Z}_{gi1}}, \ldots, \nu_{Z_{giJ_{gi}}}, \nu_{L+1}\right\}$ as a template. In the following analyses, we will create aligned data using either 1) $\mathbf{Z}=\mathbf{Z}^{(k)}, k=1,\ldots, K$, the stored MCMC samples when calculating the posterior distributions of the parameters that are functions of $\mathbf{Z}$, or 2) $\mathbf{Z}=\widehat{\mathbf{Z}} :=\{\hat{{Z}}_{gij}= \arg \max_{\ell=1, \ldots, L} p(Z_{gij}=\ell \mid \mathcal{P})\}$, the \textit{maximum a posteriori} (MAP) alignment to obtain $\mathbf{M}(\hat{\mathbf{Z}})$.

\noindent {\bf Remark 4}. Because ${\mathcal{S}}_g$ is monotonic in $\nu$ given $u$, let ${\mathcal{S}}_g^{-1}(\cdot; u)$ denote its inverse.  One might be tempted to dewarp the images so that $(t,u)$ is horizontally aligned to $\left({\mathcal{S}_g^{-1}}(t; u), u\right)$. However, because a peak $T_{gij}$ varies around its mean $\mathcal{S}_{g}(\nu_{Z_{gij}})$, unless $\sigma_{\epsilon}^2=0$, the inverse mapping cannot guarantee that ${\mathcal{S}_g^{-1}}(T_{gij};u_{gi})$ is equal to $\nu_{{Z}_{gij}}$.

\section{Applications to Scleroderma Patient Subsetting}
\label{sec:application}

Our methodology is motivated by the long-term clinical objective of finding an autoantibody signature that subsets autoimmune disease patients into groups with more homogeneous phenotypes and disease trajectories. The first step is to use the GEA data to cluster patients into subgroups with potential to have different outcomes. {We used sera from well-characterized patients with scleroderma and an associated cancer identified through the IRB-approved Johns Hopkins Scleroderma Center database \citep{Shah2017}.}                                                                                To test our algorithms, we first analyze two GEA replicates each of 20 samples. Compared to the results of hierarchical clustering without pre-processing, we show our pre-processing method creates more clearly separated clusters. We also show our pre-processing improves the accuracy of cluster detection evaluated against the true matching. As a second test, we apply the pre-processing method  to GEA measurements on 76 patients with unknown clustering. We observe that the use of the pre-proecessing algorithm identifies clusters that are clearly separated and scientifically meaningful.

\subsection{Outline of Analyses}

\label{sec:analysis_outline}

For subsequent analyses, this section describes the steps of pre-processing, clustering and three metrics that evaluate the obtained clusters.

\noindent \underline{\it Pre-processing}. We apply the peak detection algorithm in Section \ref{sec:peakdetection} followed by batch effect corrections as described in Section \ref{sec:batch.effects}. We exclude the reference lane on each gel when performing the two dimensional Bayesian dewarping.  We used $T_\nu=10$ and $T_u=6$ cubic B-spline basis functions in the horizontal and vertical directions, respectively. The dewarping functions are estimated by $\{\hat{\mathcal{S}}_g=\mathcal{S}_g(\cdot, \cdot; \hat{\bm{\beta}}_g)\}$ where $\hat{\bm{\beta}}_g$ is the posterior mean estimated by the empirical average of the MCMC samples. We also obtain the MAP estimate $\widehat{\mathbf{Z}}=\{\hat{Z}_{gij}\}$. The choice of the number of bases is crucial for the estimation of the warping functions \citep[e.g.,][]{lang2004bayesian}. For example, larger values of $(T_\nu, T_u)$ define a richer class of functions that can accommodate abrupt local image deformations. Visual inspection of the alignment of the actin peaks makes clear that more parsimonious models are preferred. Further improvements in knot selection is possible using knots on a nonequidistant grid so that more knots are placed where the spatial image warping is severe and the peaks are dense.

\noindent \underline{\it Clusterings}. Given the peak-to-landmark alignment $\mathbf{Z}$, we follow Section \ref{sec::exact_matching} to obtain peak-aligned images $\mathbf{M}=\mathbf{M}(\mathbf{Z})$ and then obtain clustering solutions. For example, let $\mathbf{M}=\mathbf{M}(\hat{\mathbf{Z}})$ where $\hat{\mathbf{Z}}$ is the MAP alignment. For each pair of sample $i$ and $i'$, we calculate the pairwise distances \(d(i,i') = 1-\widehat{\sf cor}(\mathbf{M}_{gi\cdot},\mathbf{M}_{gi'\cdot})\) where $\mathbf{M}_{gi\cdot} = (M_{gi1}, \ldots, M_{giB})'$ and $\widehat{\sf cor}(\cdot, \cdot)$ is the Pearson's correlation coefficient. Denote the $N$ by $N$ matrix of pairwise distances by $\hat{D}=\{d(i,i')\}$. We use $\hat{D}$ in the standard agglomerative hierarchical clustering with complete linkage to produce a dendrogram $\widehat{\mathcal{T}}={\mathcal{T}}(\hat{D})$. By varying the level of cutting the dendrogram $\widehat{\mathcal{T}}$, we obtain a nested set of clusterings $\hat{\mathcal{C}}(n)$, $n=2, \ldots, N$. We similarly denote the dendrogram produced without pre-processing by $\widehat{\mathcal{T}}^0 = {\mathcal{T}}(D^0)$ where $D^0$ is the correlation-based distance matrix computed from $\mathbf{M}^{0}$. We denote the nested clusters by  $\hat{\mathcal{C}}^{0}(n)$, $n=2, \ldots, N$. We will evaluate the obtained clusters by three criteria below.

\noindent \underline{\it Adjusted Rand Index}. We assess the agreement between two clusterings of the identical set of observations using the adjusted Rand index ({\sf aRI}; \citet{hubert1985comparing}). {\sf aRI} is defined by
\begin{eqnarray}
{\sf aRI}(\mathcal{C}, \mathcal{C}') & = & \frac{\sum_{r,c}{n_{rc}\choose 2}-\left[\sum_r {n_{r\cdot}\choose 2}\sum_c {n_{\cdot c}\choose 2}\right]/{N\choose 2}}{0.5\left[\sum_r {n_{r\cdot}\choose 2}+\sum_c {n_{\cdot c}\choose 2}\right]-\left[\sum_r {n_{r\cdot}\choose 2}\sum_c {n_{\cdot c}\choose 2}\right]/{N\choose 2}},
\end{eqnarray}
where $n_{rc}$ represents the number of observations                                                                                placed in the $r$th cluster of the first partition $\mathcal{C}$ and in the $c$th cluster of the second partition $\mathcal{C}'$, $\sum_{r,c}{n_{rc}\choose 2} (\leq 0.5\left[\sum_r {n_{r\cdot}\choose 2}+\sum_c {n_{\cdot c}\choose 2}\right])$ is the number of observation pairs placed in the same cluster in both partitions and $\sum_r {n_{r\cdot}\choose 2}$ and $\sum_c {n_{\cdot c}\choose 2}$ calculates the number of pairs placed in the same cluster for the first and the same cluster for second partition, respectively. {\sf aRI} is bounded between $-1$ and $1$ and corrects for chance agreement. It equals one for identical clusterings and is on average zero for two random partitions; larger values indicate better agreements between the two clustering methods.

\noindent \underline{\it Average silhouette}. We also evaluate the strength of each clustering method using the \textit{average silhouette} \citep{rousseeuw1987silhouettes}. For observation $i$, its silhouette $s(i)$ for a partition $\mathcal{C}$ compares the within- to the between-cluster average distances: $s(i)=[b(i)-a(i)]/\max\{a(i),b(i)\}$ where $a(i)$ is the average distance of $i$ to all other observations within the same cluster and $b(i)=\min_{C\in \mathcal{C}: i\notin C} \frac{\sum_{i'\in C} d(i,i')}{|C|}$ is the minimum average distance between $i$ and a cluster not containing $i$. $s(i)$ lies in $[-1,1]$ where a large value indicates observation $i$ is in a tight and isolated cluster. A larger \textit{average} silhouette $\bar{s}(\mathcal{C})=N^{-1}\sum_i s(i)$ indicates more clearly separated and tighter clustering $\mathcal{C}$.

\noindent \underline{\it Confidence levels of clusters}. In addition to the alignment uncertainty addressed by the posterior distribution $[\mathbf{Z}\mid \mathcal{P}]$, another source of uncertainty is the clustering of the aligned high-frequency intensity data $(\mathbf{t}, \mathbf{M}(\mathbf{Z}))$ given $\mathbf{Z}$. In this paper we chose not to specify the full probability distribution for the continuous intensities $(\mathbf{t}, \mathbf{M}(\mathbf{Z}))$. Following \citet{shimodaira2004approximately} and \citet{efron1996bootstrap}, we use bootstrap resampling to assess the confidence in the estimated dendrogram $\widehat{\mathcal{T}}$ (setting $\mathbf{Z}=\widehat{\mathbf{Z}}$, the \textit{MAP}). The bootstrap method perturbs the data by randomly sampling the columns of $(\mathbf{t}, \mathbf{M}(\widehat{\mathbf{Z}}))$ with replacement  and assesses the confidence levels for the presence of each subtree in $\widehat{\mathcal{T}}$. We calculate the frequency with which a subtree appears in an estimated dendrogram across all the bootstrap iterations where a large value (e.g., $> 0.95$) indicates strong evidence. We similarly bootstrap $(\mathbf{t}^0, \mathbf{M}^0)$ to assess the confidence in the dendrogram $\widehat{\mathcal{T}}^0$ estimated without alignment.

\subsection{Replication Experiments}
\label{sec::replicates}

Each of 20 biological samples were tested with two different lengths of exposure to autoradiographic devices: long (two-week) versus short (one-week) exposure. We ran 40 lanes on two gels that form 20 replicate pairs. Each gel image has 20 sample lanes: 19 serum sample lanes plus one reference lane comprised of molecules with known weights. The posterior dewarping results are shown in Appendix Figure S2.

We assess the agreement between the estimated clustering $\hat{\mathcal{C}}^{(k)}(n)$ and the true replication-based clusters $\mathcal{C}^*$ by {\sf aRI($\hat{\mathcal{C}}^{(k)}(n), \mathcal{C}^*$)}, for the number of clusters $n=2, \ldots, 20$ and the stored MCMC iteration $k=1, \ldots, K$. At iteration $k$, $\hat{\mathcal{C}}^{(k)}(n)$ is the clustering solution obtained by cutting the dendrogram that hierarchically clusters the peak-matched data $\mathbf{M}(\mathbf{Z}^{(k)})$ where $\mathbf{Z}^{(k)}$ is drawn from the posterior $[\mathbf{Z}\mid \mathcal{P}]$.

The pre-processing enhances the hierarchical clustering to produce clusters closer to the true replicate pairs. In Figure \ref{fig:replication_aRI}, the posterior mean of the adjusted Rand indices based on $K=5,000$ saved MCMC samples (solid line, $K^{-1}\sum_k$ {\sf aRI($\hat{\mathcal{C}}^{(k)}(n), \mathcal{C}^*$)}) are uniformly higher than the adjusted Rand indices based on data without pre-processing (dashed line, {\sf aRI}($\hat{\mathcal{C}^0}(n), \mathcal{C}^*$)). In the bottom panel, for every $n$, the posterior distribution for the difference between the two {aRI}s excludes zero increases with the numbers of clusters. In addition, the pairwise distances in $\hat{D}$ (obtained from $\mathbf{M}(\widehat{\mathbf{Z}})$) decreased by between $6.2$ and $66.4\%$ (mean $26.9\%$) relative to those in $D^0$. These decreases in the distances result in a dendrogram $\widehat{\mathcal{T}}$ that puts 13 replicate pairs at the terminal leaves as compared to 8 in  $\widehat{\mathcal{T}}^0$.

We also observed uniformly increased confidence levels of the presence of true replicate pairs upon pre-processing. Appendix Figure S3 examines the confidence levels associated with each subtree with (hierarchically clustering the MAP-aligned data $\left(\mathbf{t}, \mathbf{M}(\widehat{\mathbf{Z}})\right)$) and without pre-processing. For example, for pair 18, the estimated confidence level increases from 0.71 to 1 after pre-processing; The confidence levels for detecting the pairs 2, 8 and 11 see similar increases from $0.67, 0.79, 0.66$ to $0.97, 0.86, 1$, respectively.  The increase in confidence levels is partly explained by the tighter clusters obtained after data pre-processing: the average silhouette computed from the MAP clustering, $\bar{s}\left(\hat{\mathcal{C}}(n)\right)$, increased $14.2-117.6\%$ ($0.03-0.18$ in magnitude) for $n=2, \ldots, 20$ clusters.

\subsection{Scleroderma GEA Data without Replicates}
\label{sec::second_data}

We ran 4 GEA gels, each with 19 patient sera and one reference lane. The sera are from scleroderma patients with cancer who are all negative for common autoantibodies to RNA polymerase III, topoisomerase I and centromere proteins. We had no other prior knowledge about known or novel autoantibodies at the time the study was conducted.  The sera were loaded in random order on each gel; the reference sample comprised of known molecules was always in the first lane. In the following, we describe the estimated dewarping, alignment and the resulting clusters.

\noindent \underline{\it Dewarping}. We pre-process the four gel sets by estimating the dewarping functions $\{\mathcal{S}_g, g=1, 2, 3, 4\}$ and the peak-to-landmark alignment $\mathbf{Z}$. We first removed a few spots on the right of the gels caused by localized gel contamination and assumed absence of peaks at these spots. The posterior dewarping results are shown in Figure \ref{fig:warped_all}. Each detected peak $\{T_{gij}\}$ (blue dot) is connected to its matched MAP landmark $\hat{Z}_{gij}$ (red triangle). The vertical bundle of black curves, one per landmark, visualizes the global shape of the estimated warping functions $\hat{\mathcal{S}}_g$. Along each estimated vertical curve, the locations $\left\{\left(\hat{\mathcal{S}}_g(\nu_\ell, u), u\right), \forall u\right\}$   represent identical molecular weights.

\noindent \underline{\it Alignment to landmarks}. The marginal posterior probabilites of each landmark in a sample are shown at the bottom of Figure \ref{fig:warped_all}. For example, the posterior probability is $0.59$ for landmark 50 (about $43.4$ kDa, actin): the MAP estimate $\hat{\mathbf{Z}}$ shows that $73$ out of $76$ lanes. The marginal posterior probability is expected to further increase when more samples containing actin are analyzed via hierarchical Bayesian dewarping. Landmark 46 (about $46.6$ kDa) is another autoantibody hotspot where $54$ out of $76$ lanes have matched peaks. On the other hand, only $18$ and $1$ out of $76$ are matched to Landmarks 36 (about $59.8$ kDa) and 89 (about $23.4$ kDa), respectively. Their marginal posterior probabilities are hence low at $0.21$ and $0.01$.

An animation of  the continuous dewarping process is available at \url{https://github.com/zhenkewu/spotgear}. It matches the detected peaks $T_{gij}$ to their MAP landmarks $\hat{Z}_{gij}$ and morphs the posterior mean dewarping $\hat{\mathcal{S}}_g$ into the constant function $\mathcal{I}: (\nu, u) \mapsto (\nu, u)$. Also shown is the pre-processed high-frequency data $(\mathbf{t}, \mathbf{M}(\hat{\mathbf{Z}}))$ with exactly matched peaks as described in Section \ref{sec::exact_matching}.

\noindent \underline{\it Clusters}. Our pre-processing method removed global warping phenomena and revealed a few strong clusters. The clusters with $0.95$ confidence levels or higher are shown in red boxes in Figure \ref{fig:bootstrap_dendro_foursets} for the analyses done with pre-processing (top) and without pre-processing (bottom). A comparison of the two clustering solutions favors the pre-processing approach. For example, within the dendrogram at the top, the first cluster from the right (number 44) consists of seven sample lanes ((Set, Lane): (1,19), (4,3), (1,18), (3,8), (4,10), (2,4), (2,13)) that are enriched at roughly $32.7$ and $27.9$ kDa. This group is split into two clusters (numbers 47 and 14) for the analyses done without pre-processing. In a second example, the cluster 46 at the bottom and cluster 40 at the top are comprised of identical samples (enriched at about $103.4$ kDa). We observe the confidence level increases from 0.97 to 1 after pre-processing. Pre-processing again produced more clearly separated clusters and eliminated many large clusters that are otherwise formed at the bottom of Figure \ref{fig:bootstrap_dendro_foursets}; We observed $8.8-39.5\%$ increases in the average silhouette based on the MAP alignment.

\section{Discussion}
\label{sec:discussion}

In this article, we have developed a novel statistical  approach to pre-processing and analyzing two-dimensional image data obtained from gel electrophoresis autoradiography (GEA). Our objective is to eliminate artifactual data patterns that can confound our ability to use standard clustering algorithms such as hierarchical clustering to detect subsets of autoimmune disease patients. The hierarchical Bayesian image dewarping model provides a natural framework for assessing uncertainty in the estimated alignment and warping functions and allows us to make inferences about many functions of parameters.

In Section \ref{sec:application}, we analyzed two sets of data from scleroderma patients. For the data with replication, we showed that the adjusted Rand indices increased if we perform pre-processing prior to standard hierarchical clustering. Based on the MAP alignment, the average silhouette that measures the strength of clustering increased by $14$ to $118\%$. The pre-processing also increased the confidence levels for detecting true replicates.

For the data without replicates, we showed that our pre-processing method successfully aligned the actin peaks. It also increased the confidence levels for the clusters that appeared in both clusterings (one with  pre-processing and the other without pre-processing).

We conclude that there is added benefits of applying the pre-processing procedure prior to estimating disease subsets. We expect marginal though worthwhile gains to be achievable by using more carefully designed and tested tuning parameter selection procedure for local scoring (Section \ref{sec:peakdetection}).

In the analysis of data with out replicates (Section \ref{sec::second_data}), we grouped the samples by creating a \textit{single} dendrogram given a fixed $\mathbf{Z}=\hat{\mathbf{Z}}$, i.e., the MAP alignment. Uncertainty exists in both the alignment and the dendrogram obtained by hierarchical clustering. We have addressed the former by the posterior distribution $[\mathbf{Z} \mid \mathcal{P}]$ and the latter by bootstraping. Future work is needed to assume a likelihood involving the unknown dendrogram structure to obtain and represent its posterior uncertainty \cite[e.g.,][]{chakerian2012computational}.

Two extensions based on prior biological knowledge are the current subject of further research. First, in our hierarchical Bayesian dewarping model, we assumed that the intensity parameters $\{\lambda^*_\ell\}$ are shared among the samples. However, the prevalence of autoantibodies may differ by subpopulation. For example, cancer versus non-cancer patients may have distinct distributions for the abundance of certain autoantibodies. We can either add another hierarchy on top of $\{\lambda_\ell^*\}$ or develop regression models for $\{\lambda^*_\ell\}$ to incorporate disease phenotype information and covariates such as age and gender.

Second, proteins in the cells tend to work in complexes, so multiple autoantibodies are likely to be produced {against} a particular protein complex. This mechanism can be represented by a binary matrix $\mathbf{E}_{C\times L}$ where the $c$-th row $\left(E_{c1}, \ldots, E_{cL}\right)$ is a multivariate binary vector with $1$ for presence of landmark $\ell$ in complex $c$ and $0$ otherwise. The complexes are then assembled via $\bm{\eta}_{N\times L}=\mathbf{A}\mathbf{E}$ to produce the actual presence or absence of the landmarks for every patient, where $\mathbf{A}$ is a $N\times C$ binary matrix where each row represents the presence or absence of the $C$ complexes. Prior biological knowledge can be readily implemented via constraints on $\mathbf{A}$ or $\mathbf{E}$. For example, $A_{i1}=1$ for all the samples acknowledges the universal presence of autoantibodies produced against complex $1$, e.g., actin and likely others.  $\mathbf{A}$ and $\mathbf{E}$ can be inferred from the alignment indicators $\mathbf{Y}$ and continuous intensities. We may use regularization or shrinkage priors in a Bayesian framework to encourage a few maximally different complexes \citep[e.g.,][]{broderick2013}. One practical advantage of the Bayesian factorization approach lies in its convenient accommodation of repeated GEA on the same unknown sample by placing equality constraints on the rows of $\mathbf{A}$. Finally, our latent variable formulation $\bm{\eta}=\mathbf{A}\mathbf{E}$  makes it easy to incorporate multiple sources of patient lab and phenotype data that inform $\bm{\eta}$,  facilitate subgroup definition by $\mathbf{A}$ and perform individual predictions via the posterior predictive distributions of $\bm{\eta}$ \citep[e.g.,][]{BIOM:BIOM12577, wu2016partially, wu2017nplcm}.

\vspace{-3em}
\section*{Supplementary Materials}
Supplementary Material is available at
\url{http://biostatistics.oxfordjournals.org}.

\vspace{-2em}
\section*{Acknowledgments}
We thank the Editors, Associate Editor and two referees for their constructive suggestions that improved the presentation of our method and results. Research reported in this work was partially funded through a Patient-Centered Outcomes Research Institute (PCORI) Award (ME-1408-20318), NIH grant K23 AR061439 (to AS) and generous grants from the Jerome L. Greene Foundation and the Donald B. and Dorothy L. Stabler Foundation. The Johns Hopkins Rheumatic Disease Research Core Center, where the sera were processed and banked, and the antibody assays were performed, is supported by NIH grant P30 AR-070254.

\vspace{-2em}
\begingroup
    \setlength{\bibsep}{0.4em}
\bibliographystyle{biorefs}

redpajama_set_name: RedPajamaArXiv


text:
Home EntertainmentMusicLive PerformanceTBD Fest 2015: Saturday Recap
TBD Fest 2015: Saturday Recap
Posted on 23 September 2015 by Tan The Man
With one day removed from TBD Festival's opening day, Saturday began as a fresh start of music enjoyment and alcohol consumption that unfortunately coincided with the upward trend of a heatwave. I'm a wuss when it comes to high temperatures, which explains why I have never been to Coachella or Bonnaroo nor will I probably ever.
And considering the challenges I had with today being in the mid-90s, I think not rocking out in the desert or the humid southeast would be a smart decision. I honestly don't know how everyone handles the heat without going insane. It was like a thirty-minute wait to use the water refill station; fortunately, beer and good beats temporarily diverted most people's attention.
I was excited to start my day seeing Chicago-based electronic duo Autograf at the main stage. The duo's mellow tracks were perfect for mid-afternoon, which included crowd-favorite "Dream" as well as remixes of Lorde's "Team" and Stevie Wonder's "Superstition." One of the benefits of going to these festivals is that bands are eager to perform new music, and Autograf was more than eager to unleash the new track "Heartbeat" toward the end of its set.
Like Cathedrals the day before, I pushed to see New York City-based electronic duo ASTR after having missed them at SXSW earlier this year. I was not disappointed, as Adam ASTR provided the infectious beats and Zoe ASTR pumped up the crowd with her voice and dance moves. ASTR shifted naturally between old and new hits ("We Fall Down" and "Operate" to "How Did We Get So High?") with "Activate Me" being the biggest highlight.
Manchester-native Chad Valley (joined by Pamela Martinez) took it down a few notches with his mellower dance tunes. At one point, someone in the audience shouted "play for an hour" as a dig to the thirty-minute allotted set time given to Chad Valley and most of the other afternoon artists. Big hit "Now That I'm Real (How Does it Feel)" and "Fall 4 U" were crowd pleasers.
Rapper A-1 gave an entertaining set, who started off with an instrumental-less version of "Turlian" that sounded a lot more profound as a freestyle. A-1 sampled the theme song for Nickelodeon's old sketch comedy-variety show All That that was very unexpected despite a long story about Canadian television not offering Nickelodeon during his recent trip up north. There was even a sample of Lana Del Rey's "Video Games" and the entertaining "Doin the Most."
The best show of the day was Vérité who gave a very rousing performance at the 51 Fifty stage to a relatively large crowd. Vérité glibly apologized for using expletives, but the enthusiastic crowd quickly forgave the young singer after her performances of the passionate "Wasteland" and "Weekend." As the sun started to set, the live versions of "Colors" and "Sentiment" produced some warm fuzzy feelings that I hope were not heat-related.
Joywave also had the honor of performing to the sunset. The indie rock quintet from Rochester, New York probably had the most energy of any music act I had seen so far, which also included conversing and interacting with the crowd. Vocalist Daniel Armbruster definitely seemed to get a kick out of getting the crowd involved and bouncing during "Tongues" and "Something New."
Australian-natives Cut Copy put on a hits show, performing the new wave dance music to a nonstop dance fest with the outback-tinged "Take Me Over" (think Men at Work) and "We Are Explorers." Always cognizant of where they came from, the quartet also reached back to perform "Saturdays" from their freshman album Bright Like Neon Love. That first album was released almost ten years ago, which is somewhat hard to fathom. What's even harder to fathom is that Cut Copy has been making music for about as long as Brooklyn-based duo Ratatat who performed later to a more industrial light show (i.e., there were holograms and lasers that made up for the otherwise very dark stage) that helped shift the focus towards the music (like "Montanita" and "Loud Pipes").
A lot of people lined up early to get a good spot for 23-year-old wunderkind Porter Robinson who definitely satisfied the crowd with his unique sound. Robinson's set included his hits "Sad Machine" and "Lionhearted" with remixes (Nero's "The Thrill") and many loads of confetti strings.
Tycho (formerly of Sacramento proper, but now San Francisco residents) was sandwiched between Robinson and Pretty Lights, which proved to be a perfect calm in-between two very energetic artists. Tycho's ambient electronic music ("Awake" and "Spectre") helped reset the mood and provided perfect background tunes for a nearby group of fire performers. I still wonder how watching the fire breathers and fire twirlers would have felt with the more industrial electronica from Pretty Lights. I don't think we'll ever know.
NOTE: Originally published on Blogcritics.
Categories: Live Performance, MusicTags: A-1, ASTR, Autograf, Chad Valley, Cut Copy, Joywave, Porter Robinson, Pretty Lights, Sacramento, TBD Fest, TBD Festival, Tycho, Verite, West Sacramento
Previous Previous post: TBD Fest 2015: Friday Recap
Next Next post: TBD Fest 2015: Sunday Recap

redpajama_set_name: RedPajamaCommonCrawl


This quilt is one of my Quarter 1 Goals for the 2015 Finish Along. HERE is my original post.
Grand Illusion Mystery Quilt - the Dotted Spotted Version - is finished!
Best use of a clothesline - quilt pictures!
It was a windy photoshoot!
It matches her room perfectly!
Turquoise backing and binding for pop!
Backing: Robert Kaufman Spot On Yardage from Craftsy.
Quilting: Allover small meandering stitch in white thread on Mrs.Singer, my long-arm quilting machine.
Binding: Double fold, Spot On, attached by machine.
Craftsy quilt kits are a great way to build your stash even when you don't care for the pattern! It's my favorite way to get my hands on an entire fabric collection plus background fabrics for amazing prices!
My quilt was entry #72 in the original post at the beginning of March, and you can find my March goal post HERE.
Loving the borders, can't wait to see it!
So fun - I love dots!
What a gorgeous quilt, it really does match your daughter's room, it picks up the colours in her books! Lovely!
This is absolutely gorgeous - and I love that you stuck with all dots for it!
I'm happy to see another Grand Illusion finished! Mine is waiting to be quilted. Maybe this month. The Dots are great.
What a fun quilt! Your color choices are really good and look fantastic on your daughter's bed.
What a fun quilt! The colors and layout are really striking, and look great in your daughter's room! I totally agree--clotheslines are the best for quilt photos!
Wow! What a stunning quilt. It is so exciting to have a quilt finished after all of that work. I am a fan of the clothesline for showing off quilts - outdoors in the garden makes them look so pretty.
Wow - so bright and cheerful and I love all those dots!


We got some snow! I am happy! Now wasn't that just beautiful? One beautiful snowy day and it was gone in two days. The roads weren't bad for long, it didn't stay long enough to get dirty and it was good for the crops. It puts nitrogen in the soil, you know. I loved it!
And did you see the picture of the snow dinosaur on the channel 8 news. Right here in Rector, it was!
The weekend before, of course we were honored with a beautiful preview of spring. What a glorious weekend that was! Daffodils were blooming everywhere and the crocuses were beginning to open, too.
Speaking of performances, I hope you attended "And Then They Came for Me" presented Saturday by the Rector High School drama department. It was really a heart wrenching presentation. The videos of Eva Geiringer-Schloss and Helmut (Ed) Silverberg as they narrated their stories put so much emotion behind it and brought the whole thing to life in a stark and powerful manner.
The cast seemed to really meld into their characters who were near their own ages as they showed us the story of the Holocaust and a very different lifestyle. There are some statements we dare not forget. The simple writings of Anne Frank fall into that category. There is definitely some reason for the diary of a young girl to survive against all odds for more than 60 years. The words and emotions of these two people who knew Anne and survived what she did not speak to us all and ask that we remember and not allow this history to keep repeating. Thank you Gail Burns and cast for a wonderful performance.
Rosemarie Polk spent several days in Cape Girardeau recently with her mother. Ella Ruth King died on Feb. 22. She will be missed, but her family and friends are comforted to know she is now at peace. Rest in God's love, Mrs. King.
Keat Cohn was in town this weekend for a late birthday celebration at the home of his grandmother, Jessie Malin. Rob Louder fried up some wonderful catfish and Cathy Louder made the hush puppies, with Jessie's help, of course, and everyone enjoyed the feast. Of course, Keat did not come alone. He was accompanied by his wife Heidi, daughters Mason and Hayden and son Cooper. Also there to help him celebrate were David Cohn, Meagan Louder, Margo Malin, Conner Langille, KK Manchester, Paula Benton, Lane Wozniak, Paul and Cheryl Schwager and me.
Spring training has begun, and Jay Beacham left Sunday for Florida, where he will be training his heart out and looking forward to an awesome season. He played at Williams Baptist College on Feb. 28. Mom, Tracy Horton, was there to see it all with her great-nephew, Lawson Ford. Can you just imagine the cheering that was going on?
Jay even received the game ball. As pitcher he had 6 K's, two walks and earned one run. They won 6-2. Yeah!! Sounds like a great day, a great game and a great time. The next game will be March 18 at CRC, weather permitting. Tracy says anyone wanting to go and cheer Jay and the team on is welcome to ride with her. Just give her a call and let her take you out to the ball park.
Merley Moseley and Eunice Morris have been under the weather and we want to send them our most positive thoughts. We want to see you back to normal soon. Mae Ruth and Ralph Hardin are both recuperating, and we send them our best wishes, too.
Time to close for this week. I suppose I will sit on the porch and see what the weather will be like in the next 15 minutes. Don't forget to stay in touch. You can give me a call or email me. My phone number is 595-3173 or you can drop me a note at 514 West Third Street. My email address is porchswing@newwavecomm.net. We all want to know what has made you happy or sad. Please share. Bye now!


You have a hot sauce connoisseur on your gift list? Satisfy their taste buds, with these three unique flavors from Cin Chili. Sweet, hot and smokey with our Cin Chili Hot Sauce. Cindy's Cin-Fully Hot Hot Sauce, doubles as a Jalapeño marinade for chicken or pork. Just pour it on, for an on-the-edge flavor, that will tantalize your taste buds. Bloody Mary's will never be the same with our Hell's Passion hot sauce. Watch out though; the cool up front, ends up biting you in the back. This will guarantee to wake up your senses.


All questions, concerns or comments should be directed to shippantrot@gmail.com. We can't wait to hear from you. But, be sure to check out our Frequently Asked Questions page first.
Click here if your question is about registering for the race.
Oh! And in case you were wondering how to give, just go ahead and smash that DONATE button right below.
And don't forget to friend us on Facebook before the Russians take over the account!


The Pure Nautica Discovery ($49 – $65) was created for the adventurous male in mind since it is energetic and unapologetically bold without being overbearing. It's truly one of the best men's fragrances I've smelled so far. I'm glad that I was given the opportunity to give this to my dad. He has been loyal to one fragrance for years, which is Black Suede by Avon. However, he even had to agree with me that this is one of the best fragrances out there for men.
For a limited time receive an MP3 compatible speaker your Pure Nautica Discovery fragrance purchase. Limit 1 per customer, while supplies last.


text:
Metal Roof | Roofing materials
The Components of a Residential Metal Roof System
November 23rd, 2022 | 5 min. read
A residential metal roof is much more than the metal panels you see from the street. It's actually different components coming together to form a complete metal roof system.
Metal panels make up the bulk of a metal roof, but each metal roof component is important and has a specific function to keep your home leak-free. So, before diving into a residential metal roof investment, you need to know what components you'll pay for and what they do.
This ensures you'll know what needs to be included in your metal roof estimate, so you get exactly what you pay for.
For over 30 years, the team at Bill Ragan Roofing has helped homeowners learn everything they need to know about their metal roofing project. Because of this, I'll break down the different components of a residential metal roof system.
Just know that these components are specific to a standing seam metal roof (the recommended residential metal roof system).
Roof decking
Roof decking acts as the foundation of a roof system. These wooden boards (plywood, OSB, plank) are what the rest of your metal roof and its components are installed on.
When getting a metal roof replacement, you may not have to replace your roof decking. But if there's any rotten decking, your roofing contractor will have to replace it before they can begin installing your metal roof.
Unfortunately, your roofing contractor won't know the exact amount that needs to be replaced until after they finish tearing off your old roof.
Lock strip
Lock strips hold (lock) your metal roof panels to your roof edges and keep them from lifting from strong wind. They're installed at the eaves (bottom edge of your roof at the gutter line) and rakes (the sloping side of your roof).
Lock strips will be the same color and gauge as the metal roof panels. Just like steel metal panels, the lock strips will come with a Kynar 500 painted finish.
The painted finish protects your lock strips from the elements and keeps the color you choose from fading.
High-heat ice and water shield as underlayment
Roof underlayment is a felt (15 lb, 30 lb) or synthetic material used over roof decking to provide additional protection. For a metal roof, the underlayment will be high-heat ice and water shield.
Your metal roof's high-heat ice and water shield acts as the last line of defense for your roof decking. Just know that the ice and water shield must be installed over the lock strip at the eaves.
If the lock strip is installed over the underlayment, the condensation under the ribs runs down the fascia board behind the lock strip, and eventually, it rots out. I can't stress enough that the underlayment must be installed over the lock strip on your metal roof.
The metal panels make up the bulk of your metal roof system. Standing seam metal panels are attached to the roof deck by lock strips and clips at the ribs, so they're free of any holes from nails or other fasteners.
(Metal panels before installation)
The ribs (standing seams) on the metal panels are usually 12 to 24 inches apart, depending on how long the metal panel is. Because the ribs are locked together at the seam, the metal panels can glide back and forth freely as they go through thermal expansion.
The most common type of metal panels is a 24 gauge pre-finished steel with the Kynar 500 finish. However, you can choose a different gauge, and there are multiple colors options to choose from.
The clips are fasteners that lock (clip) your metal panels together. These clips are installed every 12 to 24 inches on the seam's vertical part.
These clips are screwed to the vertical part of the panel using two gimlet screws. Once the clip is installed, the next panel snaps on top of the clip fastening it to the other panel.
These clips allow the panels to glide up and down your underlayment (high heat ice and water shield) without wearing it out.
Z bar
A z bar is a piece of metal bent in the shape of a "Z" that is 1/8th of an inch taller than the ribs on the metal panels. After the metal panels are installed, this piece of metal flashing is installed to fill in the spots between the ribs and the ridge.
This protects the panels if water comes through the rib gaps at the ridge of your roof. Once the z bars are installed, the entire metal roof is ready to lock in at the ridge capping.
Ridge capping or transition metal
Ridge capping is installed at the peak where two slopes of a roof meet. Like every other metal component, it will be the same type and color as your metal panels.
But if the metal panels are coming from a low slope roof to a steeper slope, you'll need transition metal instead of ridge capping. The transition metal will come down the roof and out 6 inches to ensure there's a straight plane across, and then it'll lock onto the Z bar.
Once your ridge capping or transition metal is installed, your metal panels will be completely locked in, and you'll have a new standing seam metal roof.
What are the pros and cons of metal roofing?
Now you know the components of a residential metal roof and their function. Make sure each of these components is included in your metal roof estimate before signing anything.
But before going any further, there's still much more to learn about residential metal roofing to determine if it's the right choice. The best way to do this is by learning about the pros and cons of a residential metal roof system.
That's why I wrote another article breaking down the main advantages and disadvantages when getting a metal roof for your home.
The team at Bill Ragan Roofing has provided high-quality metal roof replacements to homeowners in Nashville and surrounding Middle Tennessee areas since 1990. Our experience and workmanship ensure you get a metal roof that lasts for decades.
Check out The Pros and Cons of Metal Roofing to learn if a metal roof is the right roofing material for your upcoming residential roof replacement.
Does Metal Roofing Need Underlayment?
October 17th, 2022|4 min read
What's the Best Metal Roof for Residential Roofing?
What is a Hidden Fastener Metal Roof?
How Much Does Metal Roofing Cost?
The 7 Main Parts of a Metal Roof (& What They Do)
Is Adding a Metal Accent to Your Roof Replacement Worth It?
August 5th, 2022|5 min read
Is a Residential Standing Seam Metal Roof Worth It?
How Much Does a Metal Roof Cost Per Square Foot?
5 Problems With Installing a Metal Roof Over Shingles
What's the Most Common Gauge of Metal Roofing?
Everything You Need to Know About Metal Roofing
June 7th, 2022|1 min read
The Pros and Cons of a Standing Seam Metal Roof
June 1st, 2022|6 min read
The 5 Best Metal Roofing Companies in Brentwood, Tennessee
The 4 Best Metal Roofing Contractors Gallatin, TN
May 9th, 2022|3 min read
What's the Most Common Gauge Used in Residential Metal Roofing?
4 Best Metal Roof Replacement Companies in Hendersonville, Tennessee
April 25th, 2022|3 min read
What is Oil Canning on a Metal Roof and How Do Your Prevent It?
The Pros and Cons of Standing Seam Metal Roofing
Can You Do a Metal Roof and Asphalt Shingle Roof Combination?
Can a Standing Seam Metal Roof be Copper? (4 Things to Know About It)
How Much Does a Metal Roof Replacement Cost?
February 7th, 2022|6 min read
What is an Exposed Fastener Metal Roof? (& 3 Things to Know About It)
December 15th, 2021|4 min read
The 5 Best Metal Roofing Contractors in Franklin, Tennessee
What are the Disadvantages of a Metal Roof?


text:
Howard's Folly launches Touriga Nacional 2015 - Howard's Folly
PORTUGAL (Prices in Euro €) UNITED KINGDOM (Prices in GB £) HONG KONG (Prices in HK $)
Our Winery
Visiting Estremoz
Restaurant & Visitor Centre
Switch the Country Site
Howard's Folly launches Touriga Nacional 2015
Howard's Folly Wines has announced the addition of Touriga Nacional 2015 to its portfolio.
Produced from grapes from the Portalegre region, one of the eight official sub-regions of the Alentejo located along the Spanish border, Howards's said the wine shows the "typical characteristics of the region which is known for its powerful spicy red wines".
Known for its blended wines, the single varietal Touriga Nacional is a change of direction for Howard's, though David Baverstock and Howard Bilton, the duo behind the business, said the new wine would be both "exceptional" and "popular".
"I am excited to bring this single varietal red wine to the market and am confident that it has all the qualities to sell extremely well," said winemaker David Baverstock.
"The taste of dark berry fruit and plum aromatics with cumin spice notes, good balancing acidity and firm tannins on the finish will make it a popular choice".
Howard Bilton said the addition "demonstrates why we have invested heavily in our own winery and vineyards: we want the freedom to create exceptional, small batch wines and to diversify into single varietals and increase production of reserve and super- premium wines. You can expect to see more of these from us in future."
Bilton also enthused about the new labels for the wine."I am honoured and delighted that such a famous artist as Sir Michael Craig-Martin created the totally new artworks that appear on our labels," he said.
"The collaboration started as an off-the-cuff request from me when we met in Hong Kong. We shared a bottle of the Howard's Folly Howard Hodgkins reserve and he obviously thought that our wines had sufficient merit to be associated with our brand. His work will be a hard act to follow."
The new urban winery in Estremoz started production in September 2018, when the partnership announced their long-term investment in the Alentejo with the signing of a ten-year lease on vineyards in Portalegre.
"We will be investing in additional wine making equipment, including lagars, amphorae and small tanks," said Baverstock. "I am really excited about the future of the Howard's Folly project: we have all the facilities and equipment we need to make truly exceptional wines."
Howard's Folly Sponsors London Art Fair
Howard's Folly Wine
Rua General Norton de Matos,
7100-107 Estremoz
susana@howardsfollywine.com
You must be of legal drinking age to purchase alcohol.


text:
Parker Markel
PARKER MARKEL, RIGHT HANDED PITCHER
Born: September 15, 1990
Drafted: 39th Round, 1181st overall pick, 2010 (Rays)
How Acquired: Waiver Claim (from Mariners)
College: Yavapai College (CC)
Agent: KVA Sports
Markel was drafted out of junior college as a strong-armed pitcher with a fastball that reached 93 mph, along with a change and slider that both showed promise. Scouts had concerns about his delivery, though, and believed he might be limited to relief. The Rays did in fact move him to the bullpen in 2014. He generally had weak walk and strikeout rates during his time with them and they did not add him to the 40-man roster, instead letting him leave as a free agent. He spent a year out of baseball and then a year in independent ball, before signing with Seattle. He's put up impressive numbers in AAA, but had a rough time during some brief opportunities in the majors. He showed a fastball that averages 95-96 mph and relied very heavily on a slider. The Pirates claimed him off waivers.
R: 2-0-0, 1.74 ERA, 1.07 WHIP, 10.1 IP, 2.6 BB/9, 11.3 K/9
Markel dominated in a brief debut in rookie ball.
A-: 3-4-0, 3.14 ERA, 1.13 WHIP, 57.1 IP, 3.6 BB/9, 6.9 K/9
Pitching as a starter in the New York-Penn League, Markel impressed scouts enough with his stuff that Baseball America rated him the third best prospect in the league and 16th best in the Rays' system. This was in spite of concerns about his mechanics and inconsistent command. He also didn't miss a lot of bats.
A: 11-5-0, 3.53 ERA, 1.26 WHIP, 120.0 IP, 2.6 BB/9, 7.2 K/9
The Rays kept Markel in the rotation in low A all year and he got good results, with better control. He improved dramatically about a third of the way through the season due to a better change. BA rated him 25th in the Rays' system after the season.
Markel had a rough time in the Florida State League. Despite the low-offense environment, opponents batted 298/363/464 against him.
The Rays sent Markel back to the FSL and moved him full-time to the bullpen. He didn't show much improvement, but in mid-July, they moved him up to AA. Despite the solid ERA, he didn't pitch well there, either, with control being a particular problem.
AAA: 0-0-0, 4.91 ERA, 1.36 WHIP, 7.1 IP, 6.1 BB/9, 8.6 K/9
Markel returned to AA and again had an ERA that was better than his other stats. He got a late-season promotion to AAA.
The Rays sent Markel back to AA to open the season. He didn't pitch especially well there, despite greatly improved walk and K rates, but the team moved him up to AAA after a month anyway. He got much better results there, although with unimpressive walk and K rates. Tampa did not add him to their 40-man roster and he became a free agent after the season.
Markel did not appear with any team.
IND: 1-1-4, 4.14 ERA, 1.43 WHIP, 41.1 IP, 4.6 BB/9, 13.5 K/9
Markel pitched in relief for Sioux City of the independent American Association. He wasn't especially effective, with continued control problems, but he struck out a lot of hitters. In late September, Seattle signed him to a minor league contract for 2019.
AA (Sea): 2-0-1, 0.00 ERA, 0.52 WHIP, 7.2 IP, 2.3 BB/9, 21.1 K/9
AAA (Sea): 2-0-8, 2.60 ERA, 1.23 WHIP, 27.2 IP, 6.8 BB/9, 14.3 K/9
AAA: 1-0-0, 0.00 ERA, 0.50 WHIP, 6.0 IP, 4.5 BB/9, 19.5 K/9
MLB (Sea): 0-0-0, 15.43 ERA, 3.00 WHIP, 4.2 IP, 7.7 BB/9, 5.8 K/9
MLB: 0-1-0, 5.71 ERA, 1.67 WHIP, 17.1 IP, 6.8 BB/9, 10.9 K/9
Seattle sent Markel to AA, but quickly moved him up to AAA. At both levels he was nearly unhittable, although he had problems with walks in AAA. He got brief callups to the Mariners in May and July, but those obviously didn't go well at all. In late July Seattle designated Markel for assignment and the Pirates claimed him off waivers. They optioned Markel to AAA after claiming him, but called him up soon after. He was sent down once more, then spent the season's final six weeks in the majors. He posted a good K rate but had trouble throwing strikes.
When the Pirates acquired Markel, Seattle had the third worst pitching in MLB and dfa'd him anyway, so the pickup was a good indication of the clown show their pitching had become. Markel has two options left, so the Pirates may see him as useful depth. They've become obsessed with relievers who have mid-90s fastballs and poor control, and have had no success with such pitchers recently.
2019: Major league minimum
Signing Bonus: $75,000
MLB Debut: 5/12/2019
MiLB FA Eligible: N/A
MLB FA Eligible: 2025
Rule 5 Eligible: N/A
Added to 40-Man: 11/20/2015
Options Remaining: 2 (USED: 2019)
June 11, 2009: Drafted by the Detroit Tigers in the 32nd round, 960th overall pick.
June 9, 2010: Drafted by the Tampa Bay Rays in the 39th round, 1181st overall pick; signed on July 31.
November 7, 2016: Became a free agent.
September 28, 2018: Signed as a minor league free agent by the Seattle Mariners.
May 12, 2019: Contract purchased by the Seattle Mariners.
July 23, 2019: Designated for assignment by the Seattle Mariners.
July 27, 2019: Claimed off waivers from the Seattle Mariners by the Pittsburgh Pirates.


text:
Donald Trump, "I don't need GOP leaders to Win in November"
By Paul Ebeling on June 19, 2016 Comments Off on Donald Trump, "I don't need GOP leaders to Win in November"
Donald Trump asks reluctant Republican leaders to join his Presidential campaign, but says he does not need them to win in November.
In a TV interview Sunday Donald Trump said it would be "nice" if Republicans all got on board in the race to the White House
"We are going to have a great convention and I think we're going to go on to a great victory. It would be nice if the Republicans stuck together," he said Sunday. "I think that I win either way. I can win one way or the other."
"I obviously won the primaries without them," he added. "You know, I'm an outsider and I won the primaries. I do believe that we can win either way, but it would be nice if we stuck together."
"I think that honestly, they should go about their business and do a wonderful job, work on budgets, get the budgets down, get the military the kind of money they need, lots of other things," Donald Trump said. "And they should not be talking so much.
"They should go out and do their job, let me do my job. I have tremendous support from both politicians and the people. Tremendous support," he continued. "Unfortunately the media just likes to cover really a small number of people that maybe have something to say. I think they should go about their work. Let me run for President. I think I'm going to do very well."
House Speaker Paul Ryan and Senate Majority Leader Mitch McConnell have both declared their support.
The fallen elite GOP leaders including former Presidents George W. Bush (43) and George H.W. Bush ($!) have yet to endorse Donald Trump, and former GOP Presidential nominee and loser Mitt Romney is strident in his criticism of Donald Trump.
Tensions flare regularly between Donald Trump and some establishment figures within the party since he became the presumptive nominee in late May when he defeated every candidate put up against him with record Republican primary voter turnout.
The W-P reported Friday that Mr. Trump's critics led by Low Energy Jeb hope to challenge him in Cleveland next month by making changes to rules governing the convention. Do not expect anything to come of that, they have about 12 aboard Vs over 2,000 delegates from primary states won by Donald Trump.
Donald Trump said former Florida Governor Jeb Bush, 1 of the 16 Republicans vanquished during the party's primary elections, may be among those attempting to undermine him. "Jeb is working on the movement, just so you understand," he said in Las Vegas.
In Las Vegas, Donald Trump said he raised $12 to 13-M for the RNC in the past 2 days alone, and that if Republicans "don't want to help out as much, I will fund my own campaign."
Donald Trump predicted success with groups including Women and Latinos.
Jan Brewer, the former Arizona governor who has endorsed Mr. Trump, was at Saturday's event in Phoenix, AZ.
"Look at all these Women for Trump," he said at the Veteran's Memorial Coliseum at the the Arizona State Fairgrounds on a day when the outside temperatures reached 110F, and a heat warning was in effect.
"We are going to do unbelievably well with Women, Mexicans, Hispanics, and with the Latinos," he said. "It is a movement. I've been on the cover of Time magazine so many times, I feel like a super-model."
Trumpeting for Trump
Campaign, Donald, Donald Trump, GOP, Hispanics, Jeb, Latinos, Mexicans, president, Republicans, Trump, turnout, voter, White House, win, women
Donald Trump, "I don't need GOP leaders to Win in November" added by Paul Ebeling on June 19, 2016
"Women Vote Trump" Raising $30-M for Presidential Campaign
Donald Trump Extends Lead Over Hillary Clinton
Donald Trump, " When I am President, it will Always be America First."
Donald Trump Now Out Ahead of Rival Hillary Clinton by 4%


text:
HomeLesson PlansBury My Heart at Wounded Knee: Day 5: Discussion of Thought Questions
When the decision came down in Standing Bear v. Crook, Indians for the first time were granted rights as people. Explain how that could have changed the way that they could advocate for themselves, if the order had been carried out.
Students should be encouraged to think in terms of rights that accrue to citizens that are otherwise unavailable to non-citizens, including the right to sue in court.
What false promise were the Cheyennes given when they went to inspect the Indian Territory?
(5 minutes)
They were told that they could return if they didn't like it. Students should understand that the promise to return could have been either a misunderstanding or a betrayal...


text:
Haunting New TV Series Law of Perdition Premieres on Amazon Prime
SGL Entertainment and Potent Media are pleased to announce the worldwide premiere of the new chilling TV Series LAW OF PERDITION, Now Available Exclusively on Amazon Prime TV. A cross between police crime drama and supernatural horror. If you will, picture "Chicago PD" and "Law & Order" Intertwined with "Ghost Whisperer" and the "Twilight Zone." LAW OF PERDITION "Season One" is 9 episodes that will keep you on the edge of your seat and dying for more.
The TV Series is Written and Directed by Christian Grillo and Carmela Hayslett and Produced by Christian Grillo, Carmela Hayslett, David Gechman and John Martineau
Starring: Andrew Hunsicker, Michael Everett Johnson, Evangeline Young, Carmela Hayslett, Chuck Maher, Daniel Plaza, Marti Keegan, B.O. Rand, Daquan Wright, Lilly Cassimore, Aaralyn Anderson, Joe Polito, Hannah Amacker and John Martineau
Watch the Official Trailer:
SYNOPSIS: Law of Perdition sets the creepy factor sky high as it follows a detective whose cold cases are returning from the dead as ghosts to help him solve their horrific murders.
LAW OF PERDITION the TV Series is Now Streaming Worldwide on Amazon Prime:
WATCH THE MOVIE ON AMAZON PRIME (USA)
WATCH THE MOVIE ON AMAZON PRIME (UNITED KINGDOM)
WATCH THE MOVIE ON AMAZON PRIME (GERMANY)
WATCH THE MOVIE ON AMAZON PRIME (JAPAN)
STILLS FROM THE TV SERIES:
For More Info on LAW OF PERDITION Go To:
www.sglmoviestore.com/tv-shows/law-of-perdition
For More Info on Potent Media Go To:
www.facebook.com/PotentMediaLLC
Watch APPARITION Today with VIMEO ON DEMAND:
By sglentertainment • Posted in Uncategorized	• Tagged Amazon, chicago pd, ghost whisperer, law & order, Potent Media, prime, SGL Entertainment, TV, twilight zone
Laura Dern's winning streak continues! The actress won the award for Outstanding Performance by a Female Actor in a Supporting Role in a Motion Picture at tonight's 2020 SAG...
Hair accessories are having a big moment tonight at the 2020 SAG Awards. From America Ferrera's pearl headband, Michelle Williams' classic hair brooch and Priah Ferguson's...
Everything You Didn't See on TV at the 2020 SAG Awards
Couldn't snag a ticket to the 2020 SAG Awards? You're in luck. E! News is inside tonight's star-studded spectacle, which takes place at the Shrine Auditorium in Los Angeles...
The Marvelous Mrs. Maisel Cast Honors Late Co-Star Brian Tarantina After 2020 SAG Awards Win
Just marvelous. The ensemble cast of The Marvelous Mrs. Maisel won one of the biggest awards of the night at the 2020 SAG Awards on Sunday, beating out the ensembles of Barry, Fleabag,...
Fleabag's Phoebe Waller-Bridge Wins Her First SAG Award for Best Female Actor in a Comedy Series
Did anyone hit Hollywood with a bigger bang this past year than Phoebe Waller-Bridge? The British star completed the prestige trifecta Sunday night, adding a SAG Award for Outstanding...
David Harbour and Lily Allen Make Red Carpet Debut as a Couple at the 2020 SAG Awards
And just like that, it's officially official. After months of dating, Stranger Things actor David Harbour and his lady friend Lily Allen have made their red carpet debut at the 2020...
Tom Hanks Once Again Steals the Show at 2020 SAG Awards
Tom Hanks is charming as ever at the 2020 SAG Awards. Per usual, the star is wowing the crowds with his funny, relatable and utterly charismatic reactions to the various award moments....
Tony Shalhoub Dedicates His SAG Awards Win to the Late Brian Tarantina
It's a great night for the marvelous Mr. Shalhoub. Tonight, at the 2020 Screen Actors Awards, the biggest actors in Hollywood are honoring one another for their tremendous work in TV...
Tonight's the night for the 2020 SAG Awards and your favorite celeb couples are making their way down the red carpet giving us all the feels. Making a rare public appearance on the...
Leave it to Chrissy Teigen to always speak her mind. Our favorite outspoken supermodel took to Twitter on Sunday, Jan. 19 as much of Hollywood was hitting the red carpet at the 2020...
Jennifer Aniston Proves Less Is More at the 2020 SAG Awards in Elegant Satin Gown
Best Dressed at 2020 SAG Awards: Jennifer Garner, Lupita Nyong'o and More


Elite Precast Concrete were approached by Midlands-based contractor Currall, Lewis and Martin (Construction) Ltd (CLM) to provide their high strength interlocking Legato blocks as part of the Salt Barn refurbishment project being carried by Kier Construction on behalf of the Highways Authority.
The Legato blocks were to be used to line the walls of the salt barn creating additional safe storage for the road salt which is kept on hand for gritting during the winter months.
Key elements of the project included the difficult access into the building, the fact that the building floor was sloping and issue of there being thousands of tonnes of salt already stored in the building.
CLM took on the task of casting level plinths onto which the blocks could be placed and organised the timing of installation of the blocks by making sure the salt already stored on the barn was moved over to one side.
The off-loading and installing of the blocks was contracted to the specialist lifting and contracting business owned and operated by Andy Fones – AA Installations. Andy has many years of experience of operating at the highest level of safety standard whilst completing some of the most challenging lifting / placing / installing projects involving heavy precast concrete elements.
By close liaison between CLM, AA Installations and Elite Precast Concrete the blocks were safely installed in two parts – one side of the salt barn and then the other – thus insuring minimal disruption to the day-to-day running of the facility at Hilton Park.
Elite Precast Concrete carry large stocks of their high strength Legato blocks so AA Installations were able to order them straight 'off-the-shelf' and the fact that the blocks are designed to a minimum 100 year working life means that the client will have benefit of the newly refurbished salt barn as part of their annual challenge of keeping the motorways open during the winter months for many years to come.


The world's first electrogroove. In the roadbed built rail connected to the movable connector of the battery of the car. When he drives off the track, the connector is automatically disengaged.
Sweden opened the world's first electrified road, which in the process of moving can be charged commercial vehicles and passenger cars, reports The Guardian.
The plot connects a transport node in Rosersberg with the cargo terminal of the airport Stockholm-Arlanda.
On electrified road will be used special electrogrooving developed in the framework of the project eRoadArlanda. How will charging a serial commercial vehicles and passenger cars, are not reported.
In the roadbed built rail connected to the movable connector of the battery of the car. When he drives off the track, the connector is automatically disengaged.
Electricity is supplied to the third rail only in the moments when it is connected to the electric car.
According to the team members, the price of laying a road 50 times lower than the costs of tram lines, and to equip the contact rail sufficient for only 5% of roads, because of the distance between the tracks does not exceed 45 km, and a period of electric cars can travel without recharging.
eRoadArlanda — part of the initiative of the Swedish Ministry of transport aimed at achieving "carbon neutrality" by 2050.
As reported, in Ukraine, the growing demand for electric vehicles.
Previously, the autopilot Tesla experienced collisions on the child.


Pedaltown Bicycle Company considers that each task is different and demands personalized solutions. To accomplish the ideal results, the bicycle company takes time to hear out people's suggestions. Pedaltown Bicycle Company has everything you need for all your bicycle services needs.
The main goal of launching a bike shop is to earn money. Don't forget to focus on the important details. Once you are willing to work at it, you can pick up the fundamentals of operating a bike repair services pretty easily. Use these simple recommendations to start down the right path.
Irrespective of whether this is your first entrepreneurial journey or not, it's bound to be one of the most challenging things you will ever do in life. Learn as much as you could about your competition before beginning to develop your bike repair services. Plan early and have a great strategy in order to build a high performing bike shop. Do not forget to take advantage of all the available resources offered on the web.
Many customers consult comments and ratings on review websites before clicking on the website of a bike repair services. Approach several of your star customers and ask to provide feedback on your bike repair services. It's a great idea to advertise feedback that highlights your most quality aspects and selections. Always thank those customers that do take the time to leave a review, even when they suggest things that you could improve on.
Do you own Pedaltown Bicycle Company?
Due to the unique way that the eBusinessPages Business Directory is structured a Premium + Verified Listing is the most powerful way to get more clients to contact Pedaltown Bicycle Company. Not only will your listing be more noticeable and more trustworthy, it will also be displayed more prominently in our category (Bicycles and related parts and Parent Categories) and location (Memphis, TN and USA Wide) listings. Pedaltown Bicycle Company will be displayed in your competitor's listings while no ads will be shown in your own listing.
Starting at $2.78 for a Premium Verified Listing, there is every reason to give Pedaltown Bicycle Company the eBusinessPages Premium treatment today.