fitnlm の modelfun の解説を読むと分かる通り、modelfun は入力次元数を持つ列ベクトルを返すように設定しなくてはいけません。つまり modelfun(b,x) の引数 x の次元が n であったら、(nx1) の出力を返すように設定しないといけません。
modelfun = @(b,t)b(1)+b(2)*t+(b(3)+b(4)*exp(-(t-t_eq)/b(5))).*heaviside(t-t_eq);
に関しては、heaviside(t-t_eq) が返す値は1次元なので、このままで問題有りません。逆に、
modelfun = @(b,t)b(1)+b(2)*t+(b(3)+b(4)*exp(-(t-t_eq)/b(5)))*heaviside(t-t_eq);
で良いです。" .* " を入れる意味は無いと思われます。この関数は1次元入力で1次元出力であることにご注意ください。
一方、
modelfun = @(b,t)b(1)+b(2)*t+(b(3)+b(4)*exp(-(t-t_eq)/b(5))).*h;
では、t はスカラー、b は係数ベクトルです。hがサイズ (1x212) なのでこの出力も (1x212) になってしまいます。入力の次元が1 であるのに、出力が 212 次元になってしまって約束に反しています。
恐らくですが、この無名関数の書き方に t が入っているので、ベースワークスペースにある変数 t と混同してしまったのかなと思います。
modelfun = @(b,x)b(1)+b(2)*x+(b(3)+b(4)*exp(-(x-t_eq)/b(5))).*h;
こういうことですので、ご注意ください。このシンタックスで t (1x212) を代入すると、(212 x 212) を返すことになってしまって、大問題です。
具体的にやってみます
t_eq = 70;
h = zeros(1,212);
h(1,41:70) = 1/30*(1:30);
h(1,71:end) = 1;
beta0 = [0.1 0.0001 -0.2 -0.1 20];
% .*h なし
modelfun1 = @(b,x)b(1)+b(2)*x+(b(3)+b(4)*exp(-(x-t_eq)/b(5)));
modelfun1(beta0,1)
% .*h あり
modelfun2 = @(b,x)b(1)+b(2)*x+(b(3)+b(4)*exp(-(x-t_eq)/b(5))).*h;
modelfun2(beta0,1)
この状態では、modelfun2 は 1次元 --> 212次元 への写像関数になります。これを修正する必要があるかと。
何をされたいのか分からないですが、h を一次元のデータを戻すような関数に書き換えてあげるひつようがあると思います。