Derivation of Toeplitz Matrix
This reference should help.  If you carry out the matrix-vector multiplication, you will see they are the same. 

Also note that it is only possible to express the convolution as a matrix multiplication when the filter impulse response is of finite duration, i.e. a FIR filter. 

The equivalence $h_{n,m} = h_{n+q,m+q}$ comes from the "TI" part of LTI, that is, time-invariance. If you delay the input by $q$ samples, the output should be the same, just delayed in time by $q$.