Реферат: СИНГУЛЯРНОЕ РАЗЛОЖЕНИЕ В ЛИНЕЙНОЙ ЗАДАЧЕ МЕТОДА НАИМЕНЬШИХ КВАДРАТОВ

Во–первых, можно выбрать границу для относительной ошибки, которая бы отражала точность данных и точность арифметики. Если взять границу в интервале , то отбросим третье сингулярное число. При этом получим следующие наборы коэффициентов для двойной и обычной точности:

Теперь коэффициенты находятся в гораздо лучшем согласии друг с другом. Кроме того, коэффициенты стали существенно меньше, а это значит, что не будет столь большого, как прежде, взаимного уничтожения слагаемых при вычислении квадратичного многочлена. Прогнозное значение Y(1980) будет соответственно 212910000 и 214960000. Эффект обычной точности еще заметен, однако результаты уже не являются катастрофическими.

Можно также определить набор нулевых коэффициентов, соответствующих пренебрежимо малому сингулярному числу. Вот эти коэффициенты: . Для значений t от 1900 до 1970 величина функции  не превосходит 0.0017, поэтому при любом a коэффициенты можно изменить , и при этом значения, выдаваемые моделью изменятся не более чем на 0.0017a. Любой из четырех перечисленных нами наборов коэффициентов можно получить из другого подобным изменением.

Во–вторых, можно улучшить ситуацию заменой базиса. Модели

гораздо более удовлетворительны. Важно при этом то, что независимая переменная преобразуется из интервала [1900, 1970] в какой–нибудь более приемлемый интервал вроде [0, 70] или, еще лучше, [–3.5, 3.5]. Числа обусловленности при этом равны 5750 и 10.7 соответственно. последнее значение более чем приемлемо даже при счете с обычной точностью.

Удобнее всего воспользоваться стандартными способами статистического анализа, т.е. матрицу плана преобразуем к стандартизованному варианту Матрица стандартизованных данных есть матрица наблюдений с нулевым средним и дисперсией 1. Это означает, что данные берутся в виде отклонений от среднего, которое мы считаем равным 0, вводим нормировку деля каждый член столбца матрицы на корень квадратный из суммы квадратов отклонений.

Во втором случае, после преобразования матрицы плана ее обусловленность сильно уменьшается, и, соответственно, повышается точность расчетов.

Данную программу можно использовать и при решении системы линейных уравнений вместо методов Гаусса, Жордана, Холесского и пр. В приложении 2 приведен пример расчета линейной системы, которая изначально не может быть решена этими методами вследствие вырожденности матрицы коэффициентов. Тем не менее, исследуемый метод дает нам правильное решение.

ЗАКЛЮЧЕНИЕ

В работе описаны компьютерные методы решения задачи наименьших квадратов. Для использования данных методов составлена соответствующая программа на алгоритмическом языке FORTRAN. Программа апробирована, результаты тестирования показывают работоспособность программы.

Результаты данной разработки могут быть использованы в самых разнообразных расчетах, где необходимо провести аппроксимацию данных заданными функциями.

ЛИТЕРАТУРА

1.       Беллман Р. Введение в теорию матриц. -М.: Наука, 1969, 368с.

2.       Гантмахер Ф.Р. Теория матриц. -М.: Наука, 1988, 548с.

3.       Ланкастер П. Теория матриц. -М.: Наука, 1982, 387с.

4.       Лоусон Ч., Хенсон Р. Численное решение задач наименьших квадратов. М.: Статистика, 1979, 447с

5.       Марчук Г.И. Методы вычислительной математики. М.: Наука, 1980

6.       Мэйндоналд Дж. Вычислительные алгоритмы в прикладной статистике. М.: Финансы и статистика, 1988, 350с

7.       Стренг Г. Линейная алгебра и ее применения. М.: Мир, 1980, 454с

8.       Уилкинсон Дж., Райнш К. Справочник алгоритмов на языке АЛГОЛ. Линейная алгебра, М.: Машиностроение, 1976, 390с

9.       Фаддеев Д.К., Фаддеева В.Н. Вычислительные методы линейной алгебры. -М.: Физматгиз, 1963, 536с.

10.    Форсайт Дж., Малькольм М., Моулер К. Машинные методы математических вычислений. М.: Мир, 1980, 279с

11.    Харебов К.С. Компьютерные методы решения задачи наименьших квадратов и проблемы собственных значений. Владикавказ.: Изд-во СОГУ, 1995, 76 с.

ПРИЛОЖЕНИЕ 1. Исходные тексты программы

                REAL A(3,3), U(3,3), V(3,3), SIGMA(3), WORK(3),Y(3),C(3),Y0(3)

                INTEGER I,IERR, J, M, N, NM

                OPEN (6,FILE="SVD.OUT",STATUS="UNKNOWN",FORM="FORMATTED")

        OPEN (5,FILE= "SVD.IN",STATUS="UNKNOWN",FORM="FORMATTED")

140     FORMAT(3I5)

150     FORMAT(4E15.7)

        READ(5,140) NM,M,N

        DO 131 I=1,M

        READ(5,150) (A(I,J),J=1,N)

131     CONTINUE

        READ (5,150) (Y(I),I=1,M)

                CALL SVD(NM,M,N,A,SIGMA,.TRUE.,U,.TRUE.,V,IERR,WORK)

                IF(IERR.NE.0) WRITE (6,2) IERR

2             FORMAT(15H TROUBLE.IERR=,I4)

        WRITE(6,120)

120     FORMAT(/'МАТРИЦА А')

        DO 121 I=1,M

        WRITE(6,130) (A(I,J),J=1,N)

130     FORMAT(8E15.7)

121     CONTINUE

        WRITE (6,160) (Y(I),I=1,N)

160     FORMAT(/'ПРАВЫЕ ЧАСТИ'/8E15.7)

210     FORMAT(/'СИНГУЛЯРНЫЕ ЧИСЛА')

        WRITE(6,210)

                DO 3 J=1,N

                WRITE(6,6) SIGMA(J)

3             CONTINUE

        SMA=SIGMA(1)

        SMI=SIGMA(1)

        DO 211 J=2,N

        IF(SIGMA(J).GT.SMA) SMA=SIGMA(J)

        IF(SIGMA(J).LT.SMI.AND.SIGMA(J).GT.0.) SMI=SIGMA(J)

211     CONTINUE

        OBU=SMA/SMI

230     FORMAT(/'ЧИСЛО ОБУСЛОВЛЕННОСТИ=',E15.7)

        WRITE(6,230) OBU

        SIGMA1=0.

        DO 30 J=1,N

        IF(SIGMA(J) .GT. SIGMA1) SIGMA1=SIGMA(J)

        C(J)=0.

30      CONTINUE

        TAU=SIGMA1*0.1E-6

        DO 60 J=1,N

        IF(SIGMA(J).LE.TAU) GO TO 60

        S=0.

        DO 40 I=1,N

        S=S+U(I,J)*Y(I)

40      CONTINUE    

        S=S/SIGMA(J)

        DO 50 I=1,N

        C(I)=C(I) + S*V(I,J)

50      CONTINUE

60      CONTINUE

                write (6,560)

        WRITE (6,6) (C(I),I=1,3)

        DO 322 J=1,N

        SS=0.

        DO 321 I=1,M

321     SS=A(J,I)*C(I)+SS

322     Y0(J)=SS

                write (6,570)

        WRITE (6,6) (Y0(I),I=1,3)

C             WRITE(6,7)

C             DO 4 I=1,M

C             WRITE(6,6) (U(I,J),J=1,N)

C4          CONTINUE

C             WRITE(6,7)

C             DO 5 I=1,N

C             WRITE(6,6) (V(I,J),J=1,N)

C5          CONTINUE

6             FORMAT(3E15.7)

560         format(2x,'roots')

570         format(2x,'right')

7             FORMAT(1H )

                STOP

                E N D

SUBROUTINE SVD(NM,M,N,A,W,MATU,U,MATV,V,IERR,RV1)

        REAL A(NM,N),W(N),U(NM,N),V(NM,N),RV1(N)

        LOGICAL MATU,MATV

        IERR=0

        DO 100 I=1,M

        DO 100         J=1,N                    

        U(I,J)=A(I,J)

100     CONTINUE

                G=0.0

                SCALE=0.0

                ANORM=0.0

                DO 300 I=1,N

                L=I+1

                RV1(I)=SCALE*G

                G=0.0

                S=0.0

                SCALE=0.0

                IF(I.GT.M) GO TO 210

                DO 120 K=I,M

120         SCALE=SCALE+ABS(U(K,I))

                IF(SCALE.EQ.0.0) GO TO 210

                DO 130 K=I,M

                U(K,I)=U(K,I)/SCALE

                S=S+U(K,I)**2

130     CONTINUE

                F=U(I,I)

                G=-SIGN(SQRT(S),F)

                H=F*G-S

                U(I,I)=F-G

                IF(I.EQ.N) GO TO 190

                DO 150 J=L,N

                S=0.0

                DO 140 K=I,M

140         S=S+U(K,I)*U(K,J)

                F=S/H

                DO 150 K=I,M

                U(K,J)=U(K,J)+F*U(K,I)

150         CONTINUE

190     DO 200 K=I,M

200         U(K,I)=SCALE*U(K,I)

210         W(I)=SCALE*G

                G=0.0

                S=0.0

                SCALE=0.0

                IF(I.GT.M.OR.I.EQ.N) GO TO 290

                DO 220 K=L,N

220         SCALE=SCALE+ABS(U(I,K))

                IF(SCALE.EQ.0.0) GO TO 290

                DO 230 K=L,N

                U(I,K)=U(I,K)/SCALE

                S=S+U(I,K)**2

230         CONTINUE

                F=U(I,L)

                G=-SIGN(SQRT(S),F)

Страницы: 1, 2, 3, 4, 5, 6, 7, 8