* Макрос расчета канонических корреляций. * Поставляется с SPSS. * Пример вызова макроса - в конце файла. preserve. set printback=off. define cancorr (set1 =!charend('/') /set2 =!charend('/') /debug =!charend('/') !DEFAULT ('N') /KEEPSC=!charend('/') !DEFAULT ('Y') /PRCOR =!charend('/') !DEFAULT (25 ) ). preserve. !IF ( !DEBUG !EQ 'N') !THEN set printback=off mprint off. !ELSE set printback on mprint on. !IFEND . *---------------------------------------------------------------------------. * Сохранение исходного файла для возможного последующего восстановления. *---------------------------------------------------------------------------. !IF (!DEBUG !EQ 'N') !THEN SET RESULTS ON. DO IF $CASENUM=1. PRINT / "ВНИМАНИЕ: ВСЯ ВЫДАЧА, ВКЛЮЧАЯ СООБЩЕНИЯ ОБ ОШИБКАХ, ВРЕМЕННО ПРИОСТАНОВЛЕНА." / "ЕСЛИ ВАМ КАЖЕТСЯ, ЧТО ПРОГРАММА РАБОТАЕТ НЕ ПРАВИЛЬНО, ПОВТОРНО ЗАПУСТИТЕ" / "ЭТОТ МАКРОС С ДОПОЛНИТЕЛЬНЫМ АРГУМЕНТОМ /DEBUG='Y'." / "ПЕРЕД ПОВТОРНЫМ ЗАПУСКОМ НЕОБХОДИМО ЗАНОВО ЗАГРУЗИТЬ ФАЙЛ ДАННЫХ." / "УКАЗАННЫЙ ПАРАМЕТР ДАСТ ИНФОРМАЦИЮ ДЛЯ ДИАГНОСТИРОВАНИЯ ВОЗНИКАЮЩИХ ПРОБЛЕМ". END IF. !IFEND . save outfile='cc__tmp1.sav'. *---------------------------------------------------------------------------. * Построение корреляционной матрицы и передача ее в команду MATRIX. *---------------------------------------------------------------------------. * ПО УМОЛЧАНИЮ: ВЫДАЧА РЕЗУЛЬТАТОВ И ОШИБОК ОТКЛЮЧЕНА ДЛЯ ПРЕДОТВРАЩЕНИЯ ВЫВОДА МАТРИЦЫ КОРРЕЛЯЦИЙ В ОКНО РЕЗУЛЬТАТОВ *. !IF (!DEBUG='N') !THEN set results off errors off. !IFEND corr variables=!set1 !set2 /missing=listwise/matrix out(*). set errors on results listing. *---------------------------------------------------------------------------. * Вычисление канонич. корреляций и основных величин, необходимых для анализа. *---------------------------------------------------------------------------. * Команда SET помещена для изменения установки по умолчанию на число циклов MXLOOPS=40 * . SET MXLOOPS=199 MITERATE 199. matrix. get r /variables=!set1/file=*. compute p1=ncol(r). get r /file=* /names=varname/variables=!set1 !set2. compute p2=ncol(r)-p1. compute nx1=varname(1:p1). compute nv=p1+p2. compute nx2=varname((p1+1):nv). compute rr=r(4:(nv+3),1:nv). compute ns=r(3,1). compute r11=rr(1:p1,1:p1). compute r22=rr((p1+1):nv,(p1+1):nv). compute r12=rr(1:p1,(p1+1):nv). compute d1=r(2,1:p1). compute d2=r(2,(p1+1):nv). compute d1=mdiag(d1). compute d2=mdiag(d2). compute s1=d1*r11*d1. compute s12=d1*r12*d2. compute s2=d2*r22*d2. compute d1=inv(d1). compute d2=inv(d2). *---------------------------------------------------------------------------. * Вывод всех корреляций. *---------------------------------------------------------------------------. do if (p1 <= !PRCOR and p2 <= !PRCOR ). print r11 /format "f7.4"/title 'Корреляции для набора №1' /space 2 /rnames=nx1 /cnames=nx1. print r22 /format "f7.4"/title 'Корреляции для набора №2' /space 2 /rnames=nx2 /cnames=nx2. print r12 /format "f7.4"/title 'Корреляции между переменными наборов №1 и №2 ' /space 2 /rnames=nx1 /cnames=nx2. end if. *---------------------------------------------------------------------------. * R1 и r2 - разложения Холецкого для матриц корреляций r11 и r22, соответственно. *---------------------------------------------------------------------------. compute r1=chol(r11). compute r2=chol(r22). *---------------------------------------------------------------------------. * R1_inv и r2_inv - матрицы, обратные к r1 и r2, соответственно. *---------------------------------------------------------------------------. compute r1_inv=inv(r1). compute r2_inv=inv(r2). *---------------------------------------------------------------------------. * Вычисление матрицы "омега". *---------------------------------------------------------------------------. do if (p1 le p2). compute omega=t(r1_inv)*r12*r2_inv. else. compute omega=t(r2_inv)*t(r12)*r1_inv. end if. *---------------------------------------------------------------------------. * Функция SVD производит разложение по сингулярным числам матрицы "омега". *---------------------------------------------------------------------------. call svd(omega,u,lambda,v). *---------------------------------------------------------------------------. * Создание набора имен, для пометки результатов при выводе . *---------------------------------------------------------------------------. !LET !@=!NULL !LET !@1=!NULL !LET !@2=!NULL !DO !N= 1 !TO 199 !LET !@=!CONCAT(!@,!QUOTE(!N),",") !LET !@1=!CONCAT(!@1,!QUOTE(!CONCAT('CV1-',!N)),",") !LET !@2=!CONCAT(!@2,!QUOTE(!CONCAT('CV2-',!N)),",") !DOEND !LET !@=!CONCAT(!@,!QUOTE(@@)) !LET !@1=!CONCAT(!@1,!QUOTE(@@)) !LET !@2=!CONCAT(!@2,!QUOTE(@@)). Compute num={!@}. *---------------------------------------------------------------------------. * Матрица "лямбда" содержит канонические корреляции. Вывод их на печать. *---------------------------------------------------------------------------. print diag(lambda)/format "f8.3"/title 'Канонические корреляции' /space 2/rnames=num. compute dlam=diag(lambda). *---------------------------------------------------------------------------. * Вычисление собственных значений и проверки на наличие остальных канонических корреляций (вторых, третьих пар и т.д.). *---------------------------------------------------------------------------. compute eign=(1 &/ (1-dlam &**2)) - 1. compute wlam=1 &/ (1+eign). compute n=nrow(wlam). compute wilk=wlam. compute df=wlam. compute sig=wlam. compute bart2=wlam. compute tem=1. loop #l=1 to n. + compute tem=tem*wlam(n-#l+1). + compute df(n-#l+1)=(p1-n+#l)*(p2-n+#l). + compute dof=df(n-#l+1). + compute bart2(n-#l+1)=-(ns-0.5*(p1+p2+3))*ln(tem). + compute chi=bart2(n-#l+1). + compute sig(n-#l+1)=1-chicdf(chi,dof). + compute wilk(n-#l+1)=tem. end loop. compute test={wilk,bart2,df,sig}. print test /format "f8.3"/title 'Проверка, что остальные канонические корреляции равны нулю:' /space 2/rnames=num /cnames={"Стат. Уилка","Хи-квадрат"," Ст. св. "," Знач."}. *---------------------------------------------------------------------------. * Вычисление и вывод стандартизированных канонических коэффициентов для набора №1. *---------------------------------------------------------------------------. do if (p1 le p2). compute a=r1_inv*u. else. compute a=r1_inv*v. end if. do if (p2 lt p1). compute a=a(:,1:p2). end if. print a/format "f8.3"/title 'Стандарт. канонич. коэфф. для набора №1' /space 2/rnames=nx1/cnames=num. *---------------------------------------------------------------------------. * Вычисление и вывод нестандартизированных канонических коэффициентов для набора №1. *---------------------------------------------------------------------------. compute a1=d1*a. print a1/format "f8.3"/title 'Нестандарт. канонич. коэфф. для набора №1' /space 2/rnames=nx1/cnames=num. *---------------------------------------------------------------------------. * Вычисление и вывод стандартизированных канонических коэффициентов для набора №2. *---------------------------------------------------------------------------. do if (p1 le p2). compute b=r2_inv*v. else. compute b=r2_inv*u. end if. do if (p1 le p2). compute b=b(:,1:p1). end if. print b /format "f8.3"/title 'Стандарт. канонич. коэфф. для набора №1' /space 2/rnames=nx2/cnames=num. *---------------------------------------------------------------------------. * Вычисление и вывод нестандартизированных канонических коэффициентов для набора №2. *---------------------------------------------------------------------------. compute b1=d2*b. print b1/format "f8.3"/title 'Нестандарт. канонич. коэфф. для набора №1' /space 2/rnames=nx2/cnames=num. *---------------------------------------------------------------------------. * Вычисление и вывод на печать канонических нагрузок для набора №1. *---------------------------------------------------------------------------. compute tem=d1*s1*a1. print tem /format "f8.3"/title 'Канонические нагрузки для набора №1' /space 2/rnames=nx1/cnames=num. *---------------------------------------------------------------------------. * Вычисление индекса избыточности (redundancy index) как пропорции дисперсии * в наборе №1, объясненной его собственными каноническими величинами. *---------------------------------------------------------------------------. compute f1=cssq(tem)/p1. compute f1=t(f1). *---------------------------------------------------------------------------. * Вычисление и печать перекрестных нагрузок для набора №1. *---------------------------------------------------------------------------. compute tem=d1*s12*b1. print tem /format "f8.3"/title 'Перекрестные нагрузки для набора №1' /space 2/rnames=nx1/cnames=num. *---------------------------------------------------------------------------. * Вычисление индекса избыточности как пропорции дисперсии в наборе №1, * объясненной каноническими величинами набора №2. *---------------------------------------------------------------------------. compute cs3=cssq(tem)/p1. compute cs3=t(cs3). *---------------------------------------------------------------------------. * Вычисление и печать канонических нагрузок для набора №2. *---------------------------------------------------------------------------. compute tem=d2*s2*b1. print tem /format "f8.3"/title 'Канонические нагрузки для набора №2' /space 2/rnames=nx2/cnames=num. *---------------------------------------------------------------------------. * Вычисление индекса избыточности как пропорции дисперсии в наборе №2, * объясненной его собственными каноническими величинами. *---------------------------------------------------------------------------. compute f2=cssq(tem)/p2. compute f2=t(f2). *---------------------------------------------------------------------------. * Вычисление и печать перекрестных нагрузок для набора №2. *---------------------------------------------------------------------------. compute tem=d2*t(s12)*a1. print tem /format "f8.3"/title 'Перекрестные нагрузки для набора №2' /space 2/rnames=nx2/cnames=num. *---------------------------------------------------------------------------. * Вычисление индекса избыточности как пропорции дисперсии в наборе №2, * объясненной каноническими величинами из набора №1. *---------------------------------------------------------------------------. compute cs4=cssq(tem)/p2. compute cs4=t(cs4). *---------------------------------------------------------------------------. * Печать результатов анализа избыточности. *---------------------------------------------------------------------------. compute c1={!@1}. compute c2={!@2}. print /title ' Анализ избыточности:' /space 2. print f1/format "f15.3" /title 'Доля дисперсии набора №1, объясненная его собственными канонич. велич.' /space 2/rnames=c1/cnames= {"Доля дисп."}. print cs3/format "f15.3" /title 'Доля дисперсии набора №2, объясненная канонич. велич. набора №2' /space 2/rnames=c2/cnames= {"Доля дисп."}. print f2/format "f15.3" /title 'Доля дисперсии набора №2, объясненная его собственными канонич. велич.' /space 2/rnames=c2/cnames= {"Доля дисп."}. print cs4/format "f15.3" /title 'Доля дисперсии набора №2, объясненная канонич. велич. набора №1' /space 2/rnames=c1/cnames= {"Доля дисп."}. *---------------------------------------------------------------------------. * Создание файлов для использования при вычислении значений канонических переменных (canonical scores). *---------------------------------------------------------------------------. SAVE {P1,P2} / OUTFILE 'CC__SIZE.SAV'. SAVE {T(A1),T(B1)} / OUTFILE 'CC__AB.SAV' . END MATRIX. *---------------------------------------------------------------------------. * Создание файла с именами переменных и переменной номера набора . *---------------------------------------------------------------------------. SET MESSAGES OFF RESULTS OFF. SELECT IF $CASENUM=1. DO REPEAT V=!SET1. COMPUTE V=1. END REPEAT. DO REPEAT V=!SET2. COMPUTE V=2. END REPEAT. STRING VARNAME (A8). COMPUTE VARNAME='SET_NUM'. FLIP VARIABLES !SET1 !SET2 / NEWNAMES=VARNAME . COMPUTE VARSEQ=1. SPLIT FILE BY SET_NUM. CREATE VARSEQ=CSUM(VARSEQ). SAVE OUTFILE 'CC_NAMES.SAV'. GET FILE 'CC__SIZE.SAV' . *---------------------------------------------------------------------------. * Подготовка необходимой информации для записи инструкций COMPUTE для скоринга. *---------------------------------------------------------------------------. WRITE OUTFILE 'CC__AB.INC' /'STRING @NMA001 TO @NMA',COL1 (N3), ' (A8)' /' @NMB001 TO @NMB',COL2 (N3), ' (A8)' /'VECTOR @NMA= @NMA001 TO @NMA',COL1 (N3) /' @NMB= @NMB001 TO @NMB',COL2 (N3) /'COMPUTE N_A='COL1 (N3) /'COMPUTE N_B='COL2 (N3) /'IF (SET_NUM=1) @NMA(VARSEQ)=CASE_LBL' /'IF (SET_NUM=2) @NMB(VARSEQ)=CASE_LBL' /'COMPUTE @=1' /'AGGREGATE OUTFILE "CC__SPRD.SAV" / BREAK @' / ' / N_A=MAX(N_A) / N_B=MAX(N_B)' / ' / @NMA001 TO @NMA',COL1 (N3) '=MAX (@NMA001 TO @NMA',COL1 (N3),')' / ' / @NMB001 TO @NMB',COL2 (N3) '=MAX (@NMB001 TO @NMB',COL2 (N3),')' / 'GET FILE "CC__AB.SAV"' / 'COMPUTE @=1' / 'MATCH FILES FILE * / TABLE "CC__SPRD.SAV"/BY @' / 'VECTOR @NMA= @NMA001 TO @NMA',COL1 (N3) / ' @NMB= @NMB001 TO @NMB',COL2 (N3) / ' COEF= COL1 TO @'. EXECUTE. GET FILE 'CC_NAMES.SAV'. INCLUDE FILE 'CC__AB.INC'. SET PRINTBACK OFF. *---------------------------------------------------------------------------. * Запись инструкций COMPUTE для скоринга. *---------------------------------------------------------------------------. STRING @SCNM@ (A8). STRING @OLDNM@ (A8). COMPUTE @SCNM@=CONCAT('S1_CV',STRING($CASENUM,N3)). WRITE OUTFILE 'CC__.INC' /'COMPUTE ',@SCNM@ ,'= 0'. LOOP CC@@@ = 1 TO N_A. COMPUTE @OLDNM@=@NMA(CC@@@). COMPUTE @COEF@ =COEF(CC@@@). WRITE OUTFILE 'CC__.INC' / ' +',@COEF@ (comma18.16),' * ',@OLDNM@ . END LOOP. COMPUTE @SCNM@=CONCAT('S2_CV',STRING($CASENUM,N3)). WRITE OUTFILE 'CC__.INC' /'COMPUTE ',@SCNM@ ,'= 0'. LOOP CC@@@=1 TO N_B. COMPUTE @OLDNM@=@NMB(CC@@@). COMPUTE @COEF@ =COEF(CC@@@+N_A). WRITE OUTFILE 'CC__.INC' / ' +',@COEF@ (comma18.16),' * ',@OLDNM@ . END LOOP. EXECUTE. *---------------------------------------------------------------------------. * Загрузка исходного файла данных и запуск скоринговой программы. *---------------------------------------------------------------------------. GET FILE 'cc__tmp1.sav'. INCLUDE FILE 'CC__.INC' . ERASE FILE 'CC__SIZE.SAV' . ERASE FILE 'CC__AB.INC'. ERASE FILE 'CC_NAMES.SAV'. ERASE FILE 'CC__AB.SAV'. ERASE FILE "CC__SPRD.SAV" !IF (!KEEPSC ='N') !THEN ERASE FILE 'CC__.INC'. !ELSE DO IF ($CASENUM=1). SET RESULTS ON. PRINT /'Оценки значений канонических переменных записаны в активный файл данных.' /'Был также создан файл со скоринговой программой для SPSS.' /'Чтобы использовать эту программу, откройте другой файл данных' /'С ТЕМИ ЖЕ переменными, что были использованы в данном анализе. ' /'Затем примените команду INCLUDE для запуска скоринговой программы.' /'Например :' / /'GET FILE какой-тодругойфайлданных' /'INCLUDE FILE "CC__.INC".' /'EXECUTE.'. END IF. EXECUTE. !IFEND. RESTORE. !enddefine. RESTORE. * Примеч. перев:. * Пример вызова макроса (для вызова - убрать "*" перед командой cancorr):. *cancorr set1=x1 x2 x3 /set2=y1 y2 y3 y4. * Изучение канонических корреляций для наборов переменных x1 - x3 и y1 - y4. * Возможные дополнительные параметры: * debug - запрос дополнительной информации о работе макроса ("N" или "Y"). * PRCOR - максимальная размерность выводимой на печать матрицы корреляций. * KEEPSC - оценка значений канонических переменных ("N" или "Y").