Mean values of skew diagonals of a $(n+1,n)$ matrix

Question

How mean values of skew diagonals of a $(n+1,n)$ matrix can be computed efficiently?

Here is my naive implementation:

ClearAll[build] ;
build[matrix_] := Block[
    {col,row,signal},
    {col,row} = Dimensions[matrix] ;
    signal = ConstantArray[0,2*row] ;
    Do[
        signal[[i]] = Table[If[q+p==i+1,matrix[[q,p]],Nothing],{q,1,col},{p,1,row}] ;
        signal[[i]] = Mean[Flatten[signal[[i]]]] ;
        ,{i,1,2*row,1}
    ] ;
    signal
]

Looks like it's time complexity is $O(n^3)$, can it be reduced?

Example:

n = 4 ;
ncols = n + 1 ;
nrows = n ;
matrix = Array[m,{ncols,nrows}] ;
matrix
build[matrix]
(* {{m[1,1],m[1,2],m[1,3],m[1,4]},{m[2,1],m[2,2],m[2,3],m[2,4]},{m[3,1],m[3,2],m[3,3],m[3,4]},{m[4,1],m[4,2],m[4,3],m[4,4]},{m[5,1],m[5,2],m[5,3],m[5,4]}} *)
(* {m[1,1],1/2 (m[1,2]+m[2,1]),1/3 (m[1,3]+m[2,2]+m[3,1]),1/4 (m[1,4]+m[2,3]+m[3,2]+m[4,1]),1/4 (m[2,4]+m[3,3]+m[4,2]+m[5,1]),1/3 (m[3,4]+m[4,3]+m[5,2]),1/2 (m[4,4]+m[5,3]),m[5,4]} *)

n = 4 ;
ncols = n + 1 ;
nrows = n ;
data = Range[1,2*n] ;
data = Partition[data,n,1] ;
data
build[data]
(* {{1,2,3,4},{2,3,4,5},{3,4,5,6},{4,5,6,7},{5,6,7,8}} *)
(* {1,2,3,4,5,6,7,8} *)

Answer 1

We can use the Diagonal function, but first we must "rotate" the matrix. Start by constructing a matrix with $n$ rows and $n+1$ columns:

n = 3;
mat = Table[m[irow, jcol], {irow, n}, {jcol, n + 1}];
mat // MatrixForm

$$\left( \begin{array}{cccc} m(1,1) & m(1,2) & m(1,3) & m(1,4) \\ m(2,1) & m(2,2) & m(2,3) & m(2,4) \\ m(3,1) & m(3,2) & m(3,3) & m(3,4) \\ \end{array} \right)$$

"Rotate" the elements like this

rot = Reverse @ Transpose @ mat;
rot // MatrixForm

$$\left( \begin{array}{ccc} m(1,4) & m(2,4) & m(3,4) \\ m(1,3) & m(2,3) & m(3,3) \\ m(1,2) & m(2,2) & m(3,2) \\ m(1,1) & m(2,1) & m(3,1) \\ \end{array} \right)$$

The diagonals of rot are the skew diagonals of mat. So the mean of the skew diagonals of mat can be obtained by

diags = Table[Diagonal[rot, k], {k, 1-n,n-2];
Mean /@ diags  // Column

$$\begin{array}{l} \frac{1}{2} (m(1,2)+m(2,1)) \\ \frac{1}{3} (m(1,3)+m(2,2)+m(3,1)) \\ \frac{1}{3} (m(1,4)+m(2,3)+m(3,2)) \\ \frac{1}{2} (m(2,4)+m(3,3)) \\ \end{array}$$

For $n+1$ rows and $n$ columns in the original matrix, use

n = 3;
mat = Table[m[irow, jcol], {irow, n + 1}, {jcol, n}];
rot = Reverse @ Transpose @ mat;
diags = Table[Diagonal[rot, k], {k, 2 - n, n - 1}];
Mean /@ diags

Answer 2

The following is long to comment.

Numeric: power 2 behavior

t0 = Table[H = RandomReal[{-1, 1}, {n, n + 1}];
   {n, AbsoluteTiming[HR = Reverse /@ H;
      Table[
       Total@Diagonal[HR, i], {i, -Length@HR + 
         1, +Length@HR}];][[1]]}, {n, 1000, 10000, 1000}];
ff = a x^n /. FindFit[t0, a x^n, {a, n}, x]

array: around power 2 behavior

t0 = Table[H = Array[aa, {n, n + 1}];
   {n, AbsoluteTiming[HR = Reverse /@ H;
      Table[
       Total@Diagonal[HR, i], {i, -Length@HR + 
         1, +Length@HR}];][[1]]}, {n, 100, 1000, 100}];

Answer 3

Here is a procedural style solution.

ClearAll[build] ;
build[matrix_] := Block[
  {row,col,array, signal, shift, start, count, i, j},
  {row,col} = Dimensions[matrix] ;
  array = Flatten[Transpose[matrix]]    ;
  signal = Table[0,2*col]  ;
  shift = 1 ;
  Do[
  {
    start = i ;
    If[
      i > row,
      start = i + col*shift ;
      shift++ ;
    ] ;
    count = 0 ;
    Do[
    {
      signal[[i]] += array[[start+j*col]] ;
      count++ ;
    },
    {j,0,Min[i,col]-shift,1}
    ] ;
    signal[[i]]/=count ;
  },
  {i,1,2*col,1}
  ] ;
  signal
] ;

Compare scaling with Diagonal solution:

ClearAll[LouisB] ;
LouisB[matrix_] := Block[
  {row, col, rotate, count},
  {row, col} = Dimensions[matrix] ;
  rotate = Reverse[Transpose[matrix]] ;
  Table[Mean[Diagonal[rotate,count]],{count,1-col,row-1}]
] ;

ClearAll[test] ;
test[function_,seed_:"SkewSum"] := Block[
    {data, fit},
    SeedRandom[seed] ;
    data = Table[
        Block[
            {matrix, time},
            matrix = RandomReal[{-1, 1}, {n+1, n}] ;
            time = First[AbsoluteTiming[function[matrix]]] ;
            {n, time}
        ],
        {n, 500, 5000, 500}
    ] ;
    fit = NonlinearModelFit[data,a*x^n,{a,n},x] ;
    fit = fit["Function"] ;
    Show[
        ListPlot[data, PlotTheme -> "Detailed", ImageSize -> 400, PlotStyle -> Directive[{PointSize[Large], Red}], AspectRatio -> 1/2],
        Plot[fit[x],{x,500,10000}, PlotStyle -> Black],
        PlotLabel -> fit[x]
    ]
] ;

Row[{test[LouisB],test[build]}]