地震预测问题数学建模论文

2015年北京交通大学数学建模竞赛

承诺书

我们仔细阅读了数学建模竞赛的竞赛规则。

我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。

我们知道，抄袭别人的成果是违反竞赛规则的, 如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。

我们郑重承诺，严格遵守竞赛规则，以保证竞赛的公正、公平性。如有违反竞赛规则的行为，我们将受到严肃处理。

所属学校（请填写完整的全名）：

参赛队员(打印并签名) ：1.

2.

3.

日期：年月日赛区评阅编号（由赛区组委会评阅前进行编号）：

目录

一、摘要 (2)

二、问题的重述 (3)

三、模型的假设与符号说明 (3)

（1）、模型的假设 (3)

（2）、符号说明 (4)

四、模型的建立 (4)

（1）. 由原始数据计算一次累加数据 (4)

（2）. 建立矩阵 (4)

（3）. 求u和a (2)

（4）. 时间响应方程 (2)

五、模型的求解 (3)

（1）. 矩阵求解 (3)

（2）.求u和a (3)

（3）. 时间响应方程 (4)

六、模型检验或误差分析 (7)

（1）、残差检验： (7)

（2）、后验差检验： (8)

（3）、结论 (8)

（4）、误差分析 (8)

七、模型评价 (8)

八、参考文献 (9)

地震预报问题

一、摘要

本文通过建立模型，解决了地震预报问题，简述如何利用正确的建模，从几组数据中预测下一次地震的时间与地点。

20世纪以来发生在中国8级以上大地震共四次，数据较少，很难找到规律，用神经网络等建模形式来“训练数据”较为困难，且到目前为止，地震的发生尚无规律可言，即使是相邻的两次地震，彼此之间也很难说有什么影响。对于地震来说，筛选出的四个数据，它们的时间地点没有必然的联系，所以可以将时间和地点的预测分开进行建模。

在时间上，由于数据较少，内部数据规律不能很好地表示出来，可以选择灰色预测，它是一种对含有不确定因素的系统进行预测的方法，可以弥补地震事件中的不确定性。通过关联分析几次地震时间，并对原始数据进行生成处理来寻找系统变动的规律，生成有较强规律性的数据序列，然后建立相应的微分方程模型，从而预测事物未来发展趋势的状况，预测未来某一时刻的特征量，或达到某一特征量的时间。最终运用matlab进行编程得出结果。

在地点上，根据板块漂移假说和地震带理论，发生在地震带上的可能性较高。对于地点，主要通过经纬度来确定，假设经度和纬度之间没有关系，通过分别建模预测下次地震的经度和纬度。由于数据较少，可以运用曲线拟合法，通过图像的类比得到近似的函数。最小化问题的精度，依赖于所选择的函数模型，用q个相关变量或p个附加的相关变量去拟和。通常人们将一个可能的、对不相关变量的构成都无困难的函数类型充作函数模型(本题中运用正弦函数)。求出参数后，用参数使所选择的函数模型同观测值Y 相匹配。最终得到可预测下次地震地点的函数。在求解方面，利用了mathematica软件求解。

经过建模求解，并进行误差分析之后，预测下次8级及以上的地震时间是2087年10月15日，地点是N36.55?,E106.76?（甘肃省内）

论文的末尾给出了模型优缺点的分析和评价，并提出了改进方向，如果数据再多一些，或者是要进行7级以上地震预测的话，会更准确一些。

关键词：地震预报时间地点灰色理论模型非线性回归模型

二、问题的重述

大地震的破坏性是众所周知的，为了减少大地震带来的灾难，人们提出了各种预报地震的方法，以求减少大地震产生的破坏。本赛题请大家用数学建模的方式预报下一次大地震发生的时间和地点。为了减少问题的复杂性，请根据20世纪以来我国发生8级及以上大地震的时间和地点预报下一次我国发生8级及以上大地震的时间和地点。

根据20世纪以来发生在我国八级以上地震的时间和地点，抓住主要因素，建立适当的符合地震发展规律的模型，从而预测下一次发生的时间和地点。

以下是20世纪以来发生在中国的大地震：

1920年12月16日，中国宁夏海原县（北纬36.7度，东经105.7度）1927年05月23日，中国甘肃古浪（北纬37.6度，东经102.6度）

1950年08月15日，中国西藏察隅县（北纬28.5度，东经96.0度）2008年05月12日，中国四川汶川县（北纬31.0度、东经103.4度）

对于每个数据来说，假设他们的时间和地点之间没有联系，对时间地点分开建模。

三、模型的假设与符号说明

（1）、模型的假设

1. 忽略地形和人为的影响；

2. 年月日的处理：日除以365换算成月，月除以12换算成年，年作为x，输入到模型中运算，得出结果后换算成年月日；

3. 由于使用“灰色理论预测”模型，年份作为x值过大，求出的a为0，所以假设时间倒退1900年，输入的x值减小，并且地质条件不变，预测出下一次地震时间后再加上1900年；

4.假设经纬度数据的变化趋势遵循函数，可通过曲线拟合法找出该函数。（2）、符号说明

对于时间预测模型所涉及：

X(0) ：实际时间

X(1) ：累加之后的时间

u、a : 时间响应方程中的常量

X’(1)（K+1）：时间预测值

E（k）：模型计算值与实际值之间的残差E(k) ：残差与实际值之间的相对误差X : 实际时间的平均值

S1 ：实际时间的标准差

E ：残差的平均值

S2 ：残差的标准差

C ：后验差比值

P ：小误差概率对于地点预测模型所涉及：

xi：存放拟合基点

yi：存放对应函数值

m1: 存放拟合基函数组

m ：存放拟合基函数组个数

a ：存放正规方程组系数矩阵

b ：存放正规方程组常数项

p : 存放拟合基函数组在拟合基点的函数值

pp：存放求出的线性模型拟合函数

h ：存放散点图

p1：存放拟合函数图形

xx：交叉正规方程组变量，存放线性模型拟合函数

四、模型的建立问题1：对于时间预测：

（1）. 由原始数据计算一次累加数据

问题2：对于地点预测(略写建模过程)：

令所求的拟合函数y （x ）=a 0y 0(x)+ a 1y 1(x)+ ……a m y m (x)

为确定系数a ，考虑平方和，求解归结为m+1元函数的极值问题，由多元函数极值的必要条件，引进函数内积符号得到如下线性方程组

（f, y 0

）

（f, y 1

）

（f, y m

）

-1/2[X (1)（2）+X (1)（1）]

1 -1/2[X (1)（3）+X (1)（2）] 1 -1/2[

X (1)（4）+

X (1)（3）]

1

B =

y = [

X (0)

（2）+X (0) （3）+X (0) （3）]

T

（3）. 求u 和a U = [ a ]= (B T

B)-1 By

u

（4）. 时间响应方程

X (1)（K+1）= [X (1) （1） - u/a]e^-ak + u/a （y 0

, y 0

）（y 0

, y 1

）……（y 0

, y m

）

a 0

（y 1, y 0）（y 1, y 1）……（y 1, y m ） （y m , y 0）

（y m , y 1）……（y m , y m ）

a 1 a m

1

=

方程有唯一解。设此解为a’0，a’1，…，a’m，则函数y’（x）=a’0y0(x)+ a’1y1(x)+ ……a’m y m(x)即为所求的线性模型拟合。

五、模型的求解

求解1：对于时间预测：

（1）. 矩阵求解

u

（2）.求u和a

U = [ a]= (B T B)-1By

代入数值： X ’(1)（K+1）=23.881e^0.6897k - 2.933=373.9305= X (1)（5） 所以，X ’(0)（5）=187.781；由于建模前方便灰色理论预测，时间减掉了1900年，所以1900+187.781 = 2087年10月15日。 求解2：对于地点预测

对于纬度，根据所给的四个纬度值首先画出折线图如下：

（3）. 时间响应方程 X ’(1)（K+1）= [X (1)（1） - u/a]e^-ak + u/a

根据该折线图，首先进行m次多项式拟合，由于只有四个数据，所以选择二次和三次多项式，分别得到拟合图及拟合函数：

Y(x)=42-4.62 t+0.4 t2

根据图像与函数代入数值所求的结果，并根据二三次函数的缺陷即存在递增递减，发现该类拟合误差较大。重新对折线图进行判断，与已知函数对应，预测函数约为正弦函数，通过改变Sin[kx]的系数k进行线性模型拟合，得图及所求拟合函数

Y(x)=33.6609+4.6655 Sin[1.39 x]

可以大体满足对地震纬度的预测。

对于经度，根据所给的四个经度值首先画出折线图如下：

根据该折线图，首先进行m次多项式拟合，由于只有四个数据，所以选择二次和三次多项式，分别得到拟合图及拟合函数：

Y(x)= 118.425-14.475 t+2.625 t2

Y(x)= 87.8+34.2333 t-19.25 t2+2.91667 t3

根据图像与函数代入数值所求的结果，并根据二三次函数的缺陷即存在递增递减，发现该类拟合误差较大。重新对折线图进行判断，与已知函数对应，预测函数约为正弦函数，通过改变Sin[kx]的系数k进行线性模型拟合，得图及所求拟合函数

Y(x)= 101.95+4.82124 Sin[1.56 x]

可以大体满足对地震经度的预测。

六、模型检验或误差分析

精度检验：

（1）、残差检验：

残差：E （k ）=X (0) （k ）- X ’(0) （k ）

相对残差：e(k)=[ X (0) （k ）- X ’(0) （k ）]/ X (0) （k ）

（2）、后验差检验：

① 计算X (0) 均值、标准差；

均值：X=1/N ∑X (0) （k ）=51.8295;

标准差: S 1= x

y =39.7838 （x=2, y=1/N [∑X (0) （k ）- X]2）；

② 计算残差均值、标准差；

残差均值：E=1/(N-1) ∑E(k)=7.0562；

残差标准差：S 2=

x

y =6.1432 （x=2,y=1/(N-1)[ ∑E(k)]2）；

③后验差比值：C=S 2/ S 1=0.1544<0.35；

小误差概率：P=P{|E(k)- E|2<0.6745 S 1}=1>0.95 （3）、结论

由表格得出结论：C<0.35，P>0.95，预测结果精度良好。

（4）、误差分析

系统误差，年月日换算成年在计算中因为小数位的取舍出错，但年月的预算不会有问题

七、模型评价

预测精度等级对照表

优点：1.有效数据少，避开训练数据，建立模型；

2.根据时间和地区的数据分布的不同，使用了两种模型；

3.运用灰色预测与线性拟合这两种可在小数据的基础上进行预测的方法，提

高了精度。

缺点：

1.对时间的预测中使用灰色预测，模型简单，但年月日的处理—日除以365

换算成月，月除以12换算成年，因为mathlab运算过程中自动保留

小数位，所以会导致误差；

2.由于所学知识的有限性，在线性模型拟合的过程中对函数图像的匹配不

能做到十分精确。

八、参考文献

[1] 王兵团，数学建模简明教程. 清华大学出版社北京交通大学出版社，2012

[2] 王兵团，数学实验基础. 清华大学出版社北京交通大学出版社，2006