What would be a Matlab function with two arguments that returns a Matrix (as shown) ?

I would like to create a Matlab function (with two arguments, N = Number of parameter, PO = polynomial order) that returns a matrix as shown in the linked file. Desired output Matrix

[Rik]

Have a read here and here. It will greatly improve your chances of getting an answer.

[Image Analyst]

I have no idea what those matrices are. It looks vaguely similar to meshgrid output, but it's not. I imagine, from the names of things it might have something to do with polyfit() and polyval(), but I can't figure it out. There is obviously some information missing, either you didn't tell us, or your instructor didn't tell you. So you either get it, or start guessing and trying to figure it out.

[Salaijobhaka]

Thank you for your answer ! There is no instructor for this task. I can provide more information,

1. Number of columns = N = Number of parameter
2. Total sum of each row is always less than or equal to polynomial order = PO
3. I have created those matrix manually in such a way that 1) and 2) are valid. The number in the matrix indicates the power of each parameter when one creates a higher dimensional polynomial function.

For example if I would like to create a polynomial function of order 2 with two different parameter, my function would be like this f(X1, X2) = a1+ a2*X1 + a3*X2 + a4*X1^2 + a5*X2^2 + a6*X1*X1, The matrix that I have created will be used to create polynomial function exactly like this. If you look at a matrix with N=2 and PO = 2, the first row gives X1^0*X^0 = 1 second row gives X1^0*X2^1 = X2, third row gives X1^0*X2^2 = X2^2, fourth row gives X1^1*X2^0 = X1, fifth row gives X1^1*X2^1 = X1*X2, sixth row gives X1^2*X2^0 = X1^2, this way, thus created matrix helps to build a polynomial function of order PO for N number of parameters. My problem is with programming skill, I am struggling to create a function that takes N and PO as argument an gives the corresponding matrix as shown in the excel file.

Answer 1

[Salaijobhaka]

I just found a existing solution.

https://se.mathworks.com/matlabcentral/answers/2128-multinomial-matrix

Thanks

[Roger Stafford]

The example given in Bruno Luong's answer at:

https://www.mathworks.com/matlabcentral/answers/2128

which presumably corresponds to your N=4, PO=3 case, may have the same set of four integers between 0 and 3 whose sum is 3, but the ordering of these is very different from that you give in your manually computed "Desired output Matrix". There does appear to be a comparatively simple algorithm which would reproduce the ordering which you give here. However, I am reluctant to work this out if ordering is not signifcant.

Also note that in your matrix as you show it, there seems to be one line missing. Line 53 with 1 1 1 1 appears to have been left out so that the N=4, PO=4 case has only 69 lines whereas it ought to have 70.

[Salaijobhaka]

Hi Roger,

Thank you for your comment, In my case ordering doesn't matter. I create a polynomial function and fit to the experimental data and figure out the coefficients of the polynomial function, i.e, {C} = [M]^(-1){F}, where {F} is a vector of experimental data, and [M] is 70 X 70 matrix for N=4, PO =4. Number of column of [M] is equal to the length of multinomial matrix (the matrix for which I want to create a function) and number of row of [M] is equal to number of sample points /design points that we need to create the polynomial function. Once we have {C} vector I can easily create the fitted polynomial function using the same order that was used to create the multinomial matrix.

The algorithm created by Bruno, should work for me as well, but his algorithms are a bit complicated and has more features than needed for my purpose. I would greatly appreciate if we can create simple algorithm.

Answer 2

[Roger Stafford]

@Salaijobhaka: I finally managed to write the script I mentioned earlier in the comment above. It produces the same ordering as you described in your "Desired output Matrix" file five days ago. It wasn't quite as simple as I thought it would be. It makes important use of the 'cumsum' function for appropriate addressing. The 'round' function is supposed to help protect against round-off errors for very large numbers of rows in the resulting matrix, M.

   PO = 4; N = 7;  % <-- Choose PO and N
   M = zeros(round(prod(N+(1:PO))/prod(1:PO)),PO);
   c = 1:N+1;
   M(1:N+1,1) = (0:N)';
   for ip = 2:PO
     s2 = repmat(c(N+1),1,N+1);
     s1 = c(N+1)-[0,c(2:N+1)-1];
     c = cumsum(c);
     d2 = c(N+1)-[0,c(1:N)];
     d1 = [d2(2:N+1)+1,1];
     for in = 1:N+1
       M(d1(in):d2(in),3:ip) = M(s1(in):s2(in),2:ip-1);
       M(d1(in):d2(in),2) = M(s1(in):s2(in),1)-M(s1(in),1);
       M(d1(in):d2(in),1) = repmat(M(s1(in)),d2(in)-d1(in)+1,1);
     end
   end

[Roger Stafford]

@Salaijobhaka: Note: I see that unfortunately I have inadvertently interchanged N and PO, so that my N is your PO and my PO is your N, as you have defined them in your comment. I will leave the task of reversing these two variables in the script to you, since I am rather sleepy at the present moment and might make a mistake in that undertaking.

[Salaijobhaka]

Hi Roger,

This is an excellent solution. I have already interchanged the representation of N and PO. Thank you very much. I will be using the "M" every now and then.

Answer 3

[Salaijobhaka]

The following code works for any polynomial order PO but for only N = 4, if I want less than 4 I need decrease the no of for loop and if I want N more than 4 I need to increase for loop manually, would be nice if we could parametric the for loop (may be a recursive function).

PO = 4;
a = 0;
  for i = 0:PO
      for j = 0:PO-i
          for k = 0:PO-i-j
              for l = 0:PO-i-j-k
                  a = a + 1;
                  v(a, 1) = i;
                  v(a, 2) = j;
                  v(a, 3) = k;
                  v(a, 4) = l;
              end
          end
      end
  end