Improved estimates of $n$ quantities via $n$ measurements $\newcommand{\de}{\delta}
\newcommand{\De}{\Delta}
\newcommand{\ep}{\epsilon}
\newcommand{\ga}{\gamma}
\newcommand{\Ga}{\Gamma}
\newcommand{\la}{\lambda}
\newcommand{\si}{\sigma}
\newcommand{\Si}{\Sigma}
\newcommand{\thh}{\theta}
\newcommand{\R}{\mathbb{R}}
\newcommand{\F}{\mathcal{F}}
\newcommand{\E}{\operatorname{\mathsf E}} 
\newcommand{\PP}{\operatorname{\mathsf P}}
\newcommand{\ii}[1]{\operatorname{\mathsf I}\{#1\}}$ 
One wants to estimate $n$ unknown real numbers $x_1,\dots,x_n$. Toward this end, one can measure $n$ linear combinations 
\begin{equation}
 a_{11}x_1+\dots+a_{1n}x_n,\ \dots,\ a_{n1}x_1+\dots+a_{nn}x_n
\end{equation}
of $x_1,\dots,x_n$ with coefficients $a_{ij}\in\{-1,0,1\}$. Each measurement involves a random error, which has mean zero and the same standard deviation $\si\in(0,\infty)$; all measurements are independent or, more generally, non-correlated. 
So, one knows the values of the coordinates $y_1,\dots,y_n$ of the random column vector 
\begin{equation}
 y:=Ax+\xi,
\end{equation}
where $A$ is the $n\times n$ matrix with entries $a_{ij}$, $x:=[x_1,\dots,x_n]^T$, and $\xi=[\xi_1,\dots,\xi_n]^T$ is the column vector of the errors of the $n$ measurements. 
Thus, one obtains the unbiased estimate 
\begin{equation}
 \hat x:=A^{-1}y=x+A^{-1}\xi
\end{equation}
of $x$, provided that the matrix $A$ is nonsingular; the unbiasedness means that $\E\hat x=x$. The covariance matrix of the estimation error $A^{-1}\xi$ is 
\begin{equation}
 \E A^{-1}\xi\xi^T(A^{-1})^T=\si^2(A^T A)^{-1}=:(b_{A;\,ij})_{i,j=1}^n, 
\end{equation}
so that the standard error of this estimation of $x_i$ is $\sqrt{b_{A;\,ii}}$. 
The question is this: 

What can be said about 
  \begin{equation}
 \si_n^2:=\min\big\{\max_1^n b_{A;\,ii}\,\colon A\in\{-1,0,1\}^{n\times n},\ A\text{ is nonsingular}\big\}, 
\end{equation}
  the minimum of $\max\limits_1^n\mathsf{Var}\,\hat x_i$ 
  over all choices of nonsingular $n\times n$ matrices $A$ with entries $a_{ij}\in\{-1,0,1\}$? 


For small enough $n$, the problem can be solved by direct calculation. In particular, for $n=2$ we have $\si_n^2=\si^2/2$ -- attained e.g. at $A=\begin{bmatrix}1&1\\1&-1\end{bmatrix}$. For $n=3$, $\si_n^2$ is $\si^2/2$ as well -- attained e.g. at 
$A=\begin{bmatrix}1&1&1\\1&1&-1\\1&-1&0\end{bmatrix}$ and at $A=\begin{bmatrix}1&1&1\\1&1&-1\\1&-1&1\end{bmatrix}$; the first of these two choices of $A$ may be considered a better one, because for it $\{b_{11},b_{22},b_{33}\}=\{\frac38,\frac38,\frac12\}\si^2$, whereas for the latter choice of $A$ we have $\{b_{11},b_{22},b_{33}\}=\{\frac12,\frac12,\frac12\}\si^2$. 

Of course, choosing $A=I_n$, we see that $\si_n^2\le\si^2$ for all $n$. It would be interesting to find the asymptotics of $\si_n^2$ as $n\to\infty$. 
Consider the (hopefully nearly orthogonal) imitation of the Rademacher orthogonal system given by formula 
\begin{equation}
 a_{n,s;ij}=\text{sgn}_s\,\sin\frac{\pi i j}{n+\pi/3}
\end{equation}
for $i,j=1,\dots,n$, where $\text{sgn}_s\,x:=\ii{|x|\ge s}\,\text{sign}\,x$ and $\ii\cdot$ denotes the indicator function. Then we have the following (connected) plot of the points $\Big(n,\dfrac{\tilde\si_n^2}{\si^2/n}\Big)$ for $n=2,\dots,100$, where 
\begin{equation}
 \tilde\si_n^2:=\max_1^n b_{A_{n,0.2};\,ii}
\end{equation}
and $A_{n,s}$ is the $n\times n$ matrix with entries $a_{n,s;ij}$. 

So, one may conjecture that $\si_n^2$ behaves as $\si^2/n$ for large $n$, possibly with an extra logarithmic factor. That would be an (almost) $n$ times as good as with the "naive" choice $A=I_n$. 

This question was sparked by the somewhat similar one at Optimal linear measurement operator, which in turn was inspired by the question at Can I really double my accuracy?, which in turn goes back to Mosteller.  
 A: This is certainly a well established topic. Probably originating at the NIST (which would then have been called the National Bureau of Standards.) Here is a 1958 paper on weighing designs. 
A Weighing Matrix of weight $w$ and order $n$ is an $n \times n$ $0,\pm1$ matrix $A$ with $A^tA=wI_n.$ The case $w=n$ is a Hadamard Matrix with $A^tA=nI_n.$ Other than $n=1,2$ this requires $n=4m.$ It is an open question if they exist for every such $n.$ According to the linked article $668, 716, 892, 1004, 1132, 1244, 1388, 1436, 1676, 1772, 1916, 1948,$ and  $1964$ are the only multiples of $4$ up to $2000$ for which such a design is unknown. 
So $\sigma^2_n=\frac{\sigma^2}{n}$ for $n$ such that there is a Hadamard matrix of order $n.$
I would have guessed that for such an $n$ the optimal thing (in whatever sense) to do for $n-1$ is take a Hadamard matrix and delete a row and column. You showed that in some sense that is not optimal for $4-1=3.$
I really wonder what the case is for $n=8-1=7.$ The underlying design is a Fano Plane which seems as if it must be optimal (yet it isn't). For example here is a circulant matrix: $$\left[ \begin {array}{ccccccc} -1&-1&1&-1&1&1&1\\  1
&-1&-1&1&-1&1&1\\  1&1&-1&-1&1&-1&1
\\  1&1&1&-1&-1&1&-1\\  -1&1&1&1&-1&
-1&1\\  1&-1&1&1&1&-1&-1\\  -1&1&-1&
1&1&1&-1\end {array} \right] 
$$
Here $AA^t=7I-(J-I)$ is all $-1$ except $7$'s on the diagonal. The inverse of this is $\frac14I+\frac18(J-I).$
If the top left entry is changed to $0$ then the maximum entry on the diagonal is still $\frac14$  but in four of the seven diagonal positions one has $\frac{13}{64}=.203125.$  So this is, in some sense, even better.
Nicer in another way is the circulant matrix with first row 01011-10 with $AA^t=4I$ 
Definitely better is making the entire main diagonal $0.$ So first row $0-11-1111.$ Then the entire diagonal of $(AA^t)^{-1}$ is $\frac{155}{788}\approx .1967$
A systematic look at circulant matrices seems worthwhile (and has probably been done.)
