Main Content

mvregress

다변량 선형 회귀

설명

예제

beta = mvregress(X,Y)는 설계 행렬 X를 기반으로 Y에 있는 d차원 응답 변수를 다변량 정규 회귀로 추정한 계수를 반환합니다.

예제

beta = mvregress(X,Y,Name,Value)는 하나 이상의 이름-값 쌍의 인수로 지정된 추가 옵션을 사용하여 추정한 계수를 반환합니다. 예를 들어, 회귀에 대한 추정 알고리즘, 초기 추정값 또는 최대 반복 횟수를 지정할 수 있습니다.

예제

[beta,Sigma] = mvregress(___)는 위에 열거된 구문의 입력 인수 중 하나를 사용하여 Y에 대해 추정한 d×d 분산-공분산 행렬도 반환합니다.

예제

[beta,Sigma,E,CovB,logL] = mvregress(___)는 잔차로 구성된 행렬 E, 추정한 회귀 계수의 분산-공분산 행렬 CovB 및 마지막 반복 후의 로그 가능도 목적 함수 값 logL도 반환합니다.

예제

모두 축소

서로 다른 절편을 가지면서 동일한 기울기를 갖는다고 가정하고 다변량 회귀 모델을 패널 데이터에 피팅합니다.

표본 데이터를 불러옵니다.

load('flu')

dataset형 배열 flu는 Google® 쿼리 데이터를 기반으로 하는 CDC의 전국 독감 추정값과 9개 개별 지역의 추정값을 포함합니다.

응답 변수와 예측 변수 데이터를 추출합니다.

Y = double(flu(:,2:end-1));
[n,d] = size(Y);
x = flu.WtdILI;

Y의 응답 변수는 9개 지역의 독감 추정값입니다. 1년 동안의 매 주별 관측값이 존재하므로 n = 52입니다. 응답 변수의 차원은 지역에 대응되므로 d = 9입니다. x의 예측 변수는 주별 전국 독감 추정값입니다.

독감 데이터를 지역별로 그룹화하여 플로팅합니다.

figure;
regions = flu.Properties.VarNames(2:end-1);
plot(x,Y,'x')
legend(regions,'Location','NorthWest')

Figure contains an axes object. The axes object contains 9 objects of type line. One or more of the lines displays its values using only markers These objects represent NE, MidAtl, ENCentral, WNCentral, SAtl, ESCentral, WSCentral, Mtn, Pac.

다변량 회귀 모델 yij=αj+βxij+ϵij를 피팅합니다. 여기서 i=1,,n이고 j=1,,d이며 지역 간 동시 상관관계 COV(ϵij,ϵij)=σjj를 가집니다.

9개의 절편 항과 1개의 공통된 기울기를 가지므로 추정할 회귀 계수의 개수 K = 10입니다. 입력 인수 Xd×K 크기의 설계 행렬로 구성된 요소를 n개 가진 셀형 배열이어야 합니다.

X = cell(n,1);
for i = 1:n
	X{i} = [eye(d) repmat(x(i),d,1)];
end
[beta,Sigma] = mvregress(X,Y);

betaK차원 계수 벡터 (α1,α2,,α9,β)의 추정값을 다음과 같이 포함합니다.

Sigma는 지역 간 동시 상관관계에 대한 d×d 분산-공분산 행렬 (σij)d×d(i,j=1,,d)의 추정값을 포함합니다.

피팅된 회귀 모델을 플로팅합니다.

B = [beta(1:d)';repmat(beta(end),1,d)];
xx = linspace(.5,3.5)';
fits = [ones(size(xx)),xx]*B;

figure;
h = plot(x,Y,'x',xx,fits,'-');
for i = 1:d
	set(h(d+i),'color',get(h(i),'color'));
end
legend(regions,'Location','NorthWest');

Figure contains an axes object. The axes object contains 18 objects of type line. One or more of the lines displays its values using only markers These objects represent NE, MidAtl, ENCentral, WNCentral, SAtl, ESCentral, WSCentral, Mtn, Pac.

플롯을 통해 각 회귀선이 절편은 서로 다르지만 기울기는 동일함을 알 수 있습니다. 시각적으로 검토해 보면 일부 회귀선이 다른 회귀선보다 데이터를 더 잘 피팅하는 것으로 보입니다.

서로 다른 절편과 기울기를 갖는다고 가정하고 최소제곱을 사용하여 다변량 회귀 모델을 패널 데이터에 피팅합니다.

표본 데이터를 불러옵니다.

load('flu');

dataset형 배열 flu는 Google® 쿼리를 기반으로 하는 CDC의 전국 독감 추정값과 9개 개별 지역의 추정값을 포함합니다.

응답 변수와 예측 변수 데이터를 추출합니다.

Y = double(flu(:,2:end-1));
[n,d] = size(Y);
x = flu.WtdILI;

Y의 응답 변수는 9개 지역의 독감 추정값입니다. 1년 동안의 매 주별 관측값이 존재하므로 n = 52입니다. 응답 변수의 차원은 지역에 대응되므로 d = 9입니다. x의 예측 변수는 주별 전국 독감 추정값입니다.

다변량 회귀 모델 yij=αj+βjxij+ϵij를 피팅합니다. 여기서 i=1,,n이고 j=1,,d이며 지역 간 동시 상관관계 COV(ϵij,ϵij)=σjj를 가집니다.

9개의 절편 항과 9개의 기울기 항을 가지므로 추정할 회귀 계수의 개수 K = 18입니다. Xd×K 설계 행렬로 구성된 요소를 n개 가진 셀형 배열입니다.

X = cell(n,1);
for i = 1:n
    X{i} = [eye(d) x(i)*eye(d)];
end
[beta,Sigma] = mvregress(X,Y,'algorithm','cwls');

betaK차원 계수 벡터 (α1,α2,,α9,β1,β2,,β9)의 추정값을 다음과 같이 포함합니다.

피팅된 회귀 모델을 플로팅합니다.

B = [beta(1:d)';beta(d+1:end)'];
xx = linspace(.5,3.5)';
fits = [ones(size(xx)),xx]*B;

figure;
h = plot(x,Y,'x',xx,fits,'-');
for i = 1:d
	set(h(d+i),'color',get(h(i),'color'));
end

regions = flu.Properties.VarNames(2:end-1);
legend(regions,'Location','NorthWest');

Figure contains an axes object. The axes object contains 18 objects of type line. One or more of the lines displays its values using only markers These objects represent NE, MidAtl, ENCentral, WNCentral, SAtl, ESCentral, WSCentral, Mtn, Pac.

플롯을 통해 각 회귀선이 서로 다른 절편과 기울기를 가짐을 알 수 있습니다.

모든 응답 변수 차원에 대해 단일 n×P 설계 행렬을 사용하여 다변량 회귀 모델을 피팅합니다.

표본 데이터를 불러옵니다.

load('flu')

dataset형 배열 flu는 Google® 쿼리를 기반으로 하는 CDC의 전국 독감 추정값과 9개 개별 지역의 추정값을 포함합니다.

응답 변수와 예측 변수 데이터를 추출합니다.

Y = double(flu(:,2:end-1));
[n,d] = size(Y);
x = flu.WtdILI;

Y의 응답 변수는 9개 지역의 독감 추정값입니다. 1년 동안의 매 주별 관측값이 존재하므로 n = 52입니다. 응답 변수의 차원은 지역에 대응되므로 d = 9입니다. x의 예측 변수는 주별 전국 독감 추정값입니다.

n×P 설계 행렬 X를 만듭니다. 회귀에 상수항을 포함시키기 위해 1로 구성된 열을 추가합니다.

X = [ones(size(x)),x];

다음과 같은 다변량 회귀 모델을 피팅합니다.

yij=αj+βjxij+ϵij,

여기서 i=1,,n이고 j=1,,d이며, 다음과 같은 지역 간 동시 상관관계가 존재합니다.

COV(ϵij,ϵij)=σjj.

9개의 절편 항과 9개의 기울기 항을 가지므로 추정할 회귀 계수는 18개입니다.

[beta,Sigma,E,CovB,logL] = mvregress(X,Y);

betaP×d 계수 행렬의 추정값을 포함합니다. Sigma는 지역 간 동시 상관관계에 대한 d×d 분산-공분산 행렬의 추정값을 포함합니다. E는 잔차로 구성된 행렬입니다. CovB는 회귀 계수의 추정된 분산-공분산 행렬입니다. logL은 마지막 반복 후의 로그 가능도 목적 함수 값입니다.

피팅된 회귀 모델을 플로팅합니다.

B = beta;
xx = linspace(.5,3.5)';
fits = [ones(size(xx)),xx]*B;

figure
h = plot(x,Y,'x', xx,fits,'-');
for i = 1:d
    set(h(d+i),'color',get(h(i),'color'))
end

regions = flu.Properties.VarNames(2:end-1);
legend(regions,'Location','NorthWest')

Figure contains an axes object. The axes object contains 18 objects of type line. One or more of the lines displays its values using only markers These objects represent NE, MidAtl, ENCentral, WNCentral, SAtl, ESCentral, WSCentral, Mtn, Pac.

플롯을 통해 각 회귀선이 서로 다른 절편과 기울기를 가짐을 알 수 있습니다.

입력 인수

모두 축소

다변량 회귀에 대한 설계 행렬로, 행렬 또는 행렬로 구성된 셀형 배열로 지정됩니다. n은 데이터에 있는 관측값 개수이고, K는 추정할 회귀 계수의 개수이고, p는 예측 변수의 개수이고, d는 응답 변수 행렬 Y에 있는 차원 수입니다.

  • d = 1인 경우 X를 단일 n×K 설계 행렬로 지정하십시오.

  • d > 1이고 모든 d개의 차원이 동일한 설계 행렬을 갖는 경우 X를 (셀형 배열에 포함되지 않은) 단일 n×p 설계 행렬로 지정할 수 있습니다.

  • d > 1이고 모든 n개의 관측값이 동일한 설계 행렬을 갖는 경우 X를 단일 d×K 설계 행렬을 포함하는 셀형 배열로 지정할 수 있습니다.

  • d > 1이고 n개의 관측값이 동일한 설계 행렬을 갖지 않는 경우 X를 d×K 설계 행렬을 포함하는 길이가 n인 셀형 배열로 지정하십시오.

회귀 모델에 상수항을 포함시키려면 각 설계 행렬이 1로 구성된 열을 포함해야 합니다.

mvregressX에 있는 NaN 값을 누락값으로 처리하고 X에서 누락값이 있는 행을 무시합니다.

데이터형: single | double | cell

응답 변수로, n×d 행렬로 지정됩니다. n은 데이터에 있는 관측값 개수이고, d는 응답 변수에 있는 차원 수입니다. d = 1이면 mvregressY에 있는 값을 n개의 독립된 응답 변수 값처럼 처리합니다.

mvregressY에 있는 NaN 값을 누락값으로 간주하고 이름-값 쌍의 인수 algorithm을 사용하여 지정된 추정 알고리즘에 따라 처리합니다.

데이터형: single | double

이름-값 인수

선택적 인수 쌍을 Name1=Value1,...,NameN=ValueN으로 지정합니다. 여기서 Name은 인수 이름이고 Value는 대응값입니다. 이름-값 인수는 다른 인수 뒤에 와야 하지만, 인수 쌍의 순서는 상관없습니다.

R2021a 이전 릴리스에서는 쉼표를 사용하여 각 이름과 값을 구분하고 Name을 따옴표로 묶으십시오.

예: 'algorithm','cwls','covar0',C는 공분산 행렬 C를 사용하는 공분산 가중 최소제곱 추정을 지정합니다.

추정 알고리즘으로, 'algorithm'과 함께 다음 중 하나가 쉼표로 구분되어 지정됩니다.

'mvn'일반 다변량 정규 최대가능도 추정.
'ecm'ECM 알고리즘을 통한 최대가능도 추정.
'cwls'공분산 가중 최소제곱 추정.

디폴트 알고리즘은 누락값의 존재 여부에 따라 달라집니다.

  • 전체 데이터의 경우 디폴트 값은 'mvn'입니다.

  • NaN으로 표시되는 누락값인 응답이 있으면 표본 크기가 모든 모수를 추정하기에 충분할 경우 디폴트 값은 'ecm'입니다. 그렇지 않을 경우 디폴트 알고리즘은 'cwls'입니다.

참고

algorithm의 값이 'mvn'인 경우 mvregress는 추정 전에 누락 응답 변수 값이 있는 관측값을 제거합니다.

예: 'algorithm','ecm'

회귀 계수에 대한 초기 추정값으로, 'beta0'과 함께 K개의 요소를 가진 벡터가 쉼표로 구분되어 지정됩니다. 디폴트 값은 0으로 구성된 벡터입니다.

추정 algorithm'mvn'이면 beta0 인수는 사용되지 않습니다.

분산-공분산 행렬 Sigma에 대한 초기 추정값으로, 'covar0'과 함께 d×d 양의 정부호 대칭 행렬이 쉼표로 구분되어 지정됩니다. 디폴트 값은 단위 행렬입니다.

추정 algorithm'cwls'이면 mvregress는 각 반복에서 covar0을 변경 없이 가중 행렬로 사용합니다.

Y에 대해 추정할 분산-공분산 행렬 유형으로, 'covtype'과 함께 다음 중 하나가 쉼표로 구분되어 지정됩니다.

'full'모든 d(d + 1)/2개의 분산-공분산 요소를 추정합니다.
'diagonal'분산-공분산 행렬의 d개의 대각선 요소만 추정합니다.

예: 'covtype','diagonal'

추정 알고리즘에 대한 최대 반복 횟수로, 'maxiter'과 함께 양의 정수가 쉼표로 구분되어 지정됩니다.

반복은 추정값이 수렴 허용오차 tolbetatolobj 내에 들어가거나 maxiter로 지정된 최대 반복 횟수에 도달할 때까지 계속됩니다. tolbetatolobj가 모두 0인 경우 mvregress는 수렴 검정 없이 maxiter회의 반복을 수행합니다.

예: 'maxiter',50

각 반복에서 평가할 함수로, 'outputfcn'과 함께 함수 핸들이 쉼표로 구분되어 지정됩니다. 함수는 논리값 true 또는 false를 반환해야 합니다. 각 반복에서 mvregress는 함수를 평가합니다. 결과가 true인 경우 반복이 중지됩니다. 그렇지 않은 경우 계속 반복됩니다. 예를 들어, 현재 반복 결과를 플로팅하거나 표시하고, Figure를 닫는 경우 true를 반환하는 함수를 지정할 수 있습니다.

함수는 다음과 같은 순서로 3개의 입력 인수를 받아야 합니다.

  • 현재 계수 추정값으로 구성된 벡터

  • 다음과 같은 3개의 필드를 포함하는 구조체:

    Covar분산-공분산 행렬의 현재 값
    iteration현재 반복 횟수
    fval로그 가능도 목적 함수의 현재 값

  • 다음과 같은 3개의 값을 받는 텍스트

    'init'함수가 초기화 중에 호출된 경우
    'iter'함수가 하나의 반복 후에 호출된 경우
    'done'함수가 완료 후에 호출된 경우

회귀 계수에 대한 수렴 허용오차로, 'tolbeta'와 함께 양의 스칼라 값이 쉼표로 구분되어 지정됩니다.

bt가 반복 t에서의 계수 벡터 추정값을 나타내고 τβtolbeta로 지정된 허용오차라고 하겠습니다. 회귀 계수 추정에 대한 수렴 조건은 다음과 같습니다.

btbt1<τβK(1+bt),

여기서 K는 bt의 길이이고 v는 벡터 v.의 노름입니다.

반복은 추정값이 수렴 허용오차 tolbetatolobj 내에 들어가거나 maxiter로 지정된 최대 반복 횟수에 도달할 때까지 계속됩니다. tolbetatolobj가 모두 0인 경우 mvregress는 수렴 검정 없이 maxiter회의 반복을 수행합니다.

예: 'tolbeta',1e-5

로그 가능도 목적 함수에 대한 수렴 허용오차로, 'tolobj'와 함께 양의 스칼라 값이 쉼표로 구분되어 지정됩니다.

Lt이 반복 t에서의 로그 가능도 목적 함수 값을 나타내고 τtolobj로 지정된 허용오차라고 하겠습니다. 목적 함수의 수렴 조건은 다음과 같습니다.

|LtLt1|<τ(1+|Lt|).

반복은 추정값이 수렴 허용오차 tolbetatolobj 내에 들어가거나 maxiter로 지정된 최대 반복 횟수에 도달할 때까지 계속됩니다. tolbetatolobj가 모두 0인 경우 mvregress는 수렴 검정 없이 maxiter회의 반복을 수행합니다.

예: 'tolobj',1e-5

모수 추정 분산-공분산 행렬 CovB에 대한 형식으로, 'varformat'과 함께 다음 중 하나가 쉼표로 구분되어 지정됩니다.

'beta'회귀 계수 추정값 beta에 대해서만 분산-공분산 행렬을 반환합니다.
'full'회귀 계수 추정값 beta와 분산-공분산 행렬 추정값 Sigma 모두에 대해 분산-공분산 행렬을 반환합니다.

예: 'varformat','full'

모수 추정을 위한 분산-공분산 행렬 유형으로, 'vartype'과 함께 'hessian' 또는 'fisher'가 쉼표로 구분되어 지정됩니다.

  • 값이 'hessian'인 경우 mvregress는 헤세(관측된 정보) 행렬을 사용하여 CovB를 계산합니다.

  • 값이 'fisher'인 경우 mvregress는 전체 데이터 피셔(예상된 정보) 행렬을 사용하여 CovB를 계산합니다.

'hessian' 방법은 누락 데이터로 인해 증가한 불확실성을 고려하는 반면 'fisher' 방법은 이를 고려하지 않습니다.

예: 'vartype','fisher'

출력 인수

모두 축소

추정된 회귀 계수로, 열 벡터 또는 행렬로 반환됩니다.

  • X를 단일 n×K 설계 행렬로 지정할 경우 mvregressbeta를 길이가 K인 열 벡터로서 반환합니다. 예를 들어, X가 20×5 설계 행렬이면 beta는 5×1 열 벡터입니다.

  • X를 하나 이상의 d×K 설계 행렬을 포함하는 셀형 배열로 지정할 경우 mvregressbeta를 길이가 K인 열 벡터로서 반환합니다. 예를 들어, X가 2×10 설계 행렬을 포함하는 셀형 배열이면 beta는 10×1 열 벡터입니다.

  • X를 (셀형 배열에 포함되지 않은) 단일 n×p 설계 행렬로 지정하고 Y의 차원 d가 1보다 클 경우 mvregressbeta를 p×d 행렬로서 반환합니다. 예를 들어, X가 20×5 설계 행렬이고 Y가 d = 2로 두 개의 차원을 가지면 beta는 5×2 행렬이고 피팅된 Y 값은 X × beta입니다.

Y에 있는 응답 변수에 대해 추정된 분산-공분산 행렬로, d×d 정사각 행렬로 반환됩니다.

참고

추정된 분산-공분산 행렬 Sigma는 잔차 행렬 E의 표본 공분산 행렬이 아닙니다.

피팅된 회귀 모델에 대한 잔차로, n×d 행렬로 반환됩니다.

algorithm의 값이 'ecm' 또는 'cwls'일 경우 mvregressY에 있는 누락값에 대응되는 잔차 값을 조건부 대치 값과 피팅된 값 사이의 차이로 계산합니다.

참고

algorithm의 값이 'mvn'인 경우 mvregress는 추정 전에 누락 응답 변수 값이 있는 관측값을 제거합니다.

모수 추정 분산-공분산 행렬로, 정사각 행렬로 반환됩니다.

  • varformat의 값이 'beta'(디폴트 값)일 경우 CovBbeta에 있는 계수 추정값의 분산-공분산 행렬입니다.

  • varformat의 값이 'full'일 경우 CovBbetaSigma에 있는 결합된 추정값의 분산-공분산 행렬입니다.

마지막 반복 후의 로그 가능도 목적 함수 값으로, 스칼라 값으로 반환됩니다.

세부 정보

모두 축소

다변량 정규 회귀

다변량 정규 회귀는 예측 변수의 설계 행렬을 토대로 한 d차원 응답 변수의 회귀로, 정규분포된 오차를 갖습니다. 오차는 이분산성(heteroscedastic)이고 상관관계가 있을 수 있습니다.

이 모델은 다음과 같습니다.

yi=Xiβ+ei,i=1,,n,

여기서

  • yi는 응답 변수로 구성된 d차원 벡터입니다.

  • Xi는 예측 변수의 설계 행렬입니다.

  • β는 회귀 계수로 구성된 벡터 또는 행렬입니다.

  • ei는 오차항으로 구성된 d차원 벡터로, 다음과 같은 다변량 정규분포를 갖습니다.

    ei~MVNd(0,Σ).

조건부 대치 값

기대값/조건부 최대화('ecm') 및 공분산 가중 최소제곱('cwls') 추정 알고리즘은 누락 응답 변수 값을 대치합니다.

y˜가 누락된 관측값을 나타낸다고 하겠습니다. 조건부 대치 값 Ε(y˜|y)는 관측된 데이터를 기반으로 한 누락된 관측값의 기대값입니다.

누락된 응답 변수와 관측된 응답 변수의 결합 분포는 다음과 같은 다변량 정규분포입니다.

(y˜y)~MVN{(X˜βXβ),(Σy˜Σy˜yΣyy˜Σy)}.

조건부 대치 기대값을 다변량 정규분포의 속성을 사용하여 다음 식으로 구할 수 있습니다.

Ε(y˜|y)=X˜β+Σy˜yΣy1(yXβ).

참고

mvregress는 누락 응답 변수 값만 대치합니다. 설계 행렬에서 누락값이 있는 관측값은 제거됩니다.

참고 문헌

[1] Little, Roderick J. A., and Donald B. Rubin. Statistical Analysis with Missing Data. 2nd ed., Hoboken, NJ: John Wiley & Sons, Inc., 2002.

[2] Meng, Xiao-Li, and Donald B. Rubin. “Maximum Likelihood Estimation via the ECM Algorithm.” Biometrika. Vol. 80, No. 2, 1993, pp. 267–278.

[3] Sexton, Joe, and A. R. Swensen. “ECM Algorithms that Converge at the Rate of EM.” Biometrika. Vol. 87, No. 3, 2000, pp. 651–662.

[4] Dempster, A. P., N. M. Laird, and D. B. Rubin. “Maximum Likelihood from Incomplete Data via the EM Algorithm.” Journal of the Royal Statistical Society. Series B, Vol. 39, No. 1, 1977, pp. 1–37.

버전 내역

R2006b에 개발됨