Are Lorentz transformations linear transformations?

Minkowski space is a real affine space of dimension $4$ whose space of translations is equipped with a metric of Lorentzian type.

A (real) affine space is a triple $(\mathbb A, V, \vec{})$, where $\mathbb A$ is a set whose elements are said points, $V$ is a (real) vector space and $\vec{}$ is a map $\vec{} : \mathbb A \times \mathbb A \to V$ with the following properties,

$$\forall q \in \mathbb A\:, \forall v \in V\:, \exists \mbox{ and is unique } p \in \mathbb A\quad \mbox{such that}\quad \vec{qp} = v\:,\tag{1}$$

$$\vec{pq} + \vec{qr} = \vec{pr}\quad \forall p, q,r \in \mathbb A\:.\tag{2}$$

by definition, the dimension of the affine space is that of $V$, whose elements are said translations.

From now on, if $p,q\in \mathbb A$ a $v \in V$, $$p= q+ v$$

means

$$\vec{qp}=v\:.$$

Form (1) this notation is well posed. $q+v$ is the action of the translation $v$ on the point $q$. This action is transitive and free, its existence physically corresponds to homogeneity of both space and time in special relativity.

Assuming that $V$ is finite dimensional, if one fixes $o \in \mathbb A$ and a basis $e_1,\ldots, e_n \in V$, a Cartesian coordinate system on the affine space $\mathbb A$ with origin $o$ and axes $e_1,\ldots, e_n$ is the bijective map $$\mathbb R^n \ni (x^1,\ldots, x^n) \mapsto o + \sum_{j=1}^n x^je_j \in \mathbb A$$ Once again, using (1) one sees that, in fact, the map above is bijectve and thus identifies $\mathbb A$ with $\mathbb R^n$.

Changing $o$ to $o'$ and the basis $e_1,\ldots, e_n$ to the basis $e'_1,\ldots, e_n'$, one obtains a different Cartesian coordinate system $x'^1,\ldots, x'^n$. It is simply proved that the rule to pass form the latter coordinate system to the former has the form $$x'^a = c^a+ \sum_{j=1}^n {A^a}_j x^j \tag{3}$$ for $n$ constant coefficients $c^j$ and a nonsingular $n\times n$ matrix of coefficients ${A^a}_j$.

The said matrix verifies

$$e_k = \sum_{i=1}^n {A^i}_k e'_i\tag{3'}$$

whereas the coefficients $c^k$ are the components of the vector $\vec{oo'}$.

(As a matter of fact the affine structure gives rise to a natural differentiable real analytic structure on $\mathbb A$ of dimension $n$.)

A real affine space equipped with a (pseudo)scalar product in $V$ is called (pseudo)Euclidean space.

Minkowski spacetime $\mathbb M^4$ is a (real) four dimensional affine space equipped with a pseudo scalar product $g: V\times V \to \mathbb R$ of Lorentzian type.

"Of Lorentzian type" means that there exist bases, $e_0,e_1,e_2,e_3$, in $V$ such that (I adopt here the convention $-+++$)

$$g(e_0,e_0)=-1 \:,\quad g(e_i,e_i) = 1 \mbox{ if $i=1,2,3$}\:,\quad g(e_i,e_j) =0 \mbox{ if $i\neq j$.}\tag{4}$$

These bases are called Minkowskian bases. Lorentz group $O(1,3)$ is nothing but the group of matrices $\Lambda$ connecting pairs of Minkowskian bases. It is therefore defined by

$$O(1,3) := \left\{ \Lambda \in M(4,\mathbb R) \:|\: \Lambda \eta \Lambda^t = \eta \right\}$$

where $\eta = diag(-1,1,1,1)$ is the matrix representing the metric $g$ in (4) in every Minkowskian basis.

A Minkowskian coordinate system on $\mathbb M^4$ is a Cartesian coordinate system whose axes are a Minkowskian basis.

Lorentz transformations are transformations of coordinates between pairs of Minkowskian coordinate systems with the same origin (so that $c^k=0$ in (3)). Thus they have the form

$$x'^a = \sum_{j=1}^n {\Lambda^a}_j x^j $$

for some $\Lambda \in O(1,3)$. If we admit different origins we obtain the so-called Poincaré transformations

$$x'^a = c^a+ \sum_{j=1}^n {\Lambda^a}_j x^j \:.$$

When viewing Lorentz transformations as transformation of coordinates, their formal linearity does not play a relevant physical role, since it only reflects the arbitrary initial choice of the same origin for both reference frames. However, these transformations are also transformations of bases (3') in the space of translations (the tangent space), in this case linearity is natural because it reflects the natural linear space structure of the translations.

Strictly in the sense of coordinate transforms in special relativity (i.e. not general relativity), the Lorentz transforms are actually homogeneous, not linear. Linearity is as you rightly note, is a formal property of the transformation only in a certain coordinate system, the Cartesian system. There is no need to resort to identifying the space-time with a vector space.

So what does homogeneity mean? It means that the transformation does not spoil translational invariance. Translation of a set of points puts a set of parallel lines through every point and then moves every of the points along it's line by the same distance. Lorentz transformation has this nice property that arbitrary two parallel lines anywhere in space (i.e. not only the ones in the origin) remain parallel even after the transformation.

This is a geometrical statement independent of coordinates and physically a requirement of same velocity to have an invariant meaning. Remember that we are in a space-time, so a straight line is actually a particle moving at a constant velocity through space and a family of parallel lines is a family of objects at the same velocity. That is, two objects at the same velocity are observed to have equal velocity by any observer.

We now presume we can choose a coordinate system in which every family of parallel lines can be characterized by a unique coordinate-displacement slope - modulo multiples of that slope. We require that a Lorentz transformation preserves this structure - which eventually leads to the linearity of the transform in this very special coordinate system.

This set of coordinates in which we can identify parallel lines by coordinate slopes, are the Cartesian coordinates and a certain "natural time". The constructed families of parallel lines would then be a kind of a projective space (the definition of the introduction of the wiki article would not apply, but the structure is really the same).

One of the properties of a projective space is homogeneity, the property of being invariant to multiplication by a number. Thus the name of homogeneity of Lorentz transforms - it induces an isomorphism of projective ("homogeneous") spaces.

Special relativity takes place in Minkowski space, which is the vector space $\mathbb{R}^4$ equipped with the inner product given by the matrix

$$ G = \left(\begin{matrix} -1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1\end{matrix}\right) \text{ with inner product} \langle x,y \rangle := x^T G y \; \forall \; x,y \in \mathbb{R}^4 $$

We are not free to choose any basis of this vector space if we want to preserve the flat Minkowski metric $G$ in its given from, i.e. not all basis (or coordinate) transformations in $\mathrm{GL}(4,\mathbb{R})$ are suitable for leaving physics invariant. The set of basis transformations preserving the inner product is called isometries, and that is precisely what the Lorentz group $\mathrm{SO}(1,3)$ is - the group of transformations leaving the Minkowski inner product invariant, i.e. all matrices $M$ fulfilling

$$ M^T G M = M$$

If you look at special relativity in a slightly more general context, i.e. on arbitrary four-dimensional manifolds where every tangent space carries the Minkowski metric, then you must understand a Lorentz transformation as a coordinate change on the manifold whose Jacobian is an element of $\mathrm{SO}(1,3)$ at every point (i.e. on all tangent spaces), since coordinate changes act upon tangent spaces by their Jacobians.

Either way, the Lorentz transformations $\mathrm{SO}(1,3)$ are linear transformations on the tangent spaces, induced by a coordinate transformation.