测井曲线分层

测井曲线分层问题
摘要
测井曲线分层是在地球物理勘探中利用测井资料了解地下地质情况，首先要完成的基础工作。

本文主要解决了以附件中1号井为标准井，并根据此井的各种测井曲线数据，建立了数学模型，对第2号至7号井进行自动分层，并且通过分析，与人工分层结果进行比较分析，进一步对1号井的分层结果进行说明；对第8号井至13号井进行自动分层，并给出结论两个问题。

针对问题一，本文首先通过查资料并观察附表中1号井的数据剔除了一些变化规律不明显的指标如CAL、DEVi等，筛选出了SP、GR、AC 、CNL 、RT、RILD、RILM七个显著变化的指标，根据经验又剔除了各个指标中明显错误的数据；然后利用主成分分析法的思想挑选了在主成分中贡献率较大的指标AC、CNL、SP三个指标。

接着画出三个指标的综合测井曲线，由于每一层的指标差异性、稳定性，本文采取了层内差异法，结合综合测井曲线，将每一个井进行了大致的粗分层。

最后要将相似程度高的层进行合并，而聚类分析是根据某一分类统计量来度量多个观察量之间的相似程度，依相似程度高低决定是聚合为一层，还是划为不同层，本文利用聚类分析法将第2至7号井进行细致的分层，与人工分层进行了比较，判断其精度，结果见表4、表5并对模型进行了改进，进一步提高合理性。

针对问题二，本文利用问题一中所得出的规律对第8号井至13号井进行了分层，结果见表6，并进行了分析。

关键词：测井曲线自动分层主成分分析层内差异法聚类分析
测井区县分层是在地球物理勘探中利用测井资料了解地下地质情况最基本也是最重要的问题。

目前最常用的人工分层方法不仅费时费力，而且分层取值过程中受测井分析人员的经验知识和熟练程度影响较大，主观性较强，也会因为不同的解释人员的个人标准有误差，而造成不同的人员有不同的分层结果。

本文主要解决的问题有：
1、以1号井为标准井，根据此井的各种测井曲线数据，建立合适的数学模型，对第2号至7号井进行了自动分层，并且通过分析，与人工分层结果进行比较分析并改进了数学模型，对1号井的分层结果进行说明。

2、通过前面人工分层与自动分层的比较结果，以及已给的各种测井曲线数据，对第8号井至13号井进行自动分层，并分析本文的结论。

二、问题假设
1 、 假设题目附件中所给数据仪器的精度都满足要求
2 、 1号井所给的分层都是准确的，其他井依靠1号井为标准。

3 、 根据某一条或者某几条测井曲线可以较准确的进行分层
符号说明
ij x 第i 个深度样本的第j 个指标
B 待定参数
i δ 第i 个深度样本的均方根误差 i x 第i 层测井均值
)(i x E 非地层因素引起的允许误差
)1,(+k k d 第k 层与第k + 1 层的测井值之间的“马氏”距离 min d 各层间测井值的最小“临界距离”
实现井位分层人工智能处理，也就是实现自动分层的方法有很多，比如神经网络、层次分析、层内聚类分析。

由于问题中所给数据繁多，本文选取了主成分分析法、层次分析法的思想，客观的反映了系统状况。

主成分分析法是一种将多维因子纳入同一系统进行定量化研究、理论成熟的多元统计分析方法。

通过分析变量之间的相关性，使得反应信息重叠的变量被某一主成分代替，减少了变量的数目，从而降低了评价的复杂性。

在以方差贡献率作为每一主成分的权重，每一个主成分的得分加权即可完成对每一层的评价。

层内分析法是一种定量与定性相结合的、系统化的、层次化的系统分析方法。

以随机数学为工具，通过大量的观察数据寻求统计规律。

一般来说引起测井曲线值变化的因素有两类：一类是地层因素（岩性、空隙流体性质）一类是非地层因素如井壁因素、测量系统测量条件等。

非地层因素引起的测井相应一般比较小，本文采剔除了这些因素。

对于地质因素本文采用主成分分析的思想选取了贡献率较大AC 、CNL、RT三个指标，然后对1号井的每一层数据进行分析。

由于附件中所给数据繁多直接应用主成分分析法比较繁琐，本文根据数据的规律剔除了一些明显错误的点，然后利用层内差异法找出了1号井层与层之间的区别。

又利用得出的规律对2号井进行了分层，并与人工分层进行了比较，并对分层规律进行了修正，进一步提高了分层的精确性。

四、模型的建立与求解
主成分分析法
（1）设每一深度的测量值为一个样品，每一个样品观测有AC 、CNL、GR、SP 等指标以x ij表示第i个样品第j个指标的得分，根据附表中的数据得到原始矩
阵：
X=⎪
⎪
⎪
⎪⎪
⎭
⎫ ⎝⎛x x
x x
x x x x
x np n n p p
2
1
222
21
11211 其中x i =⎪⎪⎪⎪
⎪⎭⎫ ⎝⎛x x x ni i i 21根据主成分分析的方法，分析不同深度各指标的检测数据。

做变换：j Y =
（1,2,3
j p =）
得到标准化的数据矩阵Yij =ij
X Xj -，其中2
211
11,()n n i j i i i Xj X j S X j Xj n n ====-∑∑。

(2) 在标准化数据矩阵Y=()*n p yij d 的基础上计算p 个原始数据指标相关系数矩
阵R=*()
n p
rij ,
其中：()(),
n
k i i kj j ij X
X X X r --=
∑，
(i=1…n;j=1…p)
(3) 求相关系数矩阵R 的特征值并排序；再求出R 的特征向量Z
(4) 确定主成分数目。

在确定主成分之前，需要先给出一个控制值α，令
1
1
/1q p
i i
i i λλ
α==≥-∑∑，则对应满足条件的q 的最小值即为保留的主成分的个数m 。

(5) 计算综合得分。

首先计算得到的第i 个样本中第k 个主成分的得分为， 用数据矩阵x 的p 个向量（即p 个指标向量）x 1，x 2，x 3， ,x p 做线性组合（即综合指标向量）为：
⎪⎪
⎩
⎪⎪
⎨⎧+++=+++=+++=x a x a x a F x a x a x a F x a x a x a F p pp p p p
p p p 221132222212212221111 再以m 个主成分的方差贡献率为权重求得第i 个样本的综合得分。

（6）分别以第一主成分，第二主成分，第三主成分中线性组合中系数的绝对值得的作为评判标准。

表1 主成分贡献率
主成分的主要贡献率
特征值 贡献率 累计贡献率 第一主成分 1195.0 0.5910 0.5190 第二主成分 719.4 0.2402 0.8316 第三主成分 219.8 0.0954 0.9270
表2 三个主成分的7个指标系数
第一主
分系数
0.3363 0.1235 -0.8967 0.2376 0.1001 -0.0320 0.0000
第二主分系数
0.3835 0.3424 -0.0049 -0.8293 0.2184 0.0149 -0.0000 第三主分系数 0.2145 0.5376 0.3411 0.4702 0.5705 0.0464 0.0000
本文选取了第一到第三主成分中系数的绝对值比较大的作为分层的标准，由于第一主成分所占的权重比较大，本文以第一主成分中绝对之比较大的AC 作为分层的主要根据，SP 、CNL 两个指标为辅助分层依据。

为了比较明确的看出分层结果，由他们的贡献率作为他们的权重，得到综合测井曲线为： RT CNL AC M *095.0*24.0*591.0++=， 则M 随深度变化的综合曲线为
100
120140160180200220240260
图一 综合测井曲线
由综合测井曲线图联系一号井的标准分层数据可得下图：
100
150200250
250300
350
400
450
500
100120140160180200220240260
图二 一号井的综合测井曲线图与分层图
利用综合测井曲线对2到7号井进行分层
1、首先运用层内分析法将井进行粗分层 层内差异法
层内差异分层法的依据是：同一层内的测井值相对稳定，其值的变化不超过某一允许误差(即由非地层因素引起的测井值误差) ，并认为每一层内采样值的均值代表该层的真实测值。

一相邻的采样点的值与该均值的差异在允许误差范围内， 则认为该采样点属于这一层，否则，便属于下一层。

其一般过程为：先选择具有较强纵向分辨率的测井曲线(称为主动曲线) 进行细分层(本区选择自然电位曲线和微电极曲线) ，然后再兼顾其它测井曲线的特点，将分层结用到其它曲线上。

下面对本方法加以详细叙述。

首先，为减小因测井仪器标准化程度不高、泥浆、矿化度变化等非地层因素对测
井值的影响， 必须对测井曲线作归一化处理，以突出测井值的相对变化，归一化的公式为：
X=
100min
max min
⨯--x x x x
式中，x 为实际测井值m ax x 、m in x 为本条测井曲线所有采样点中的最大值和最小值，而X 则为归一化后的测井值(0 ≤ X ≤100) 。

假设已确定出相邻的几个采样值ij x ( j = 1 , 2⋯⋯ n ) 均于第i 层, 该层的允许误差为)(i x E ，均值为i x ,方差为2i s ，其中
∑==n
j ij i x n x 1
1,
∑-=-=n
j i
i ij n s x x 12
2)(11
现在的问题是要判断随后的第n + 1 个采样值x n 1+ ，是否同属于第i 层。

判定规则如下：
若)(1i n n x E x x ≤-+则认为1+n x 属于第i 层，并计算这n + 1 个采样点的均值、方差及允许误差，接着进行下一采样点的处理。

若)(1i i n x E x x >-+ ，则认为1+n x 点不属于第i 层，该层划分完毕，把前n 个采样点的测井均值及第n 个采样点对应的深度输出。

并从1+n x 点开始进行下一层的划分。

如此进行下去，直至处理完整个井段为止。

不难看出，上述分层的关键是确定允许误差函数)(i x E ，这里尝试用误差理论与概率统计相结合的办法建立这一函数。

从概率统计角度讲，可以认为测井值X 是一个具有有限方差的随机变量。

同一层内各采样值的差异反映了非地层因素引起的随机误差， 并满足切比雪夫不等式:
i
i i i x E x E x x P )())((2
δ≤
≥-
式中 i δ———均方根误差；
x ———该层的某个采样值； i x ———第i 层测井均值；
)(i x E ———非地层因素引起的允许误差 。

由误差理论可求得均方根误差为:
∑--=
)(21
1
ij i n x x i δ 且 i i B x E δ=)(
式中， B 为待定参数， 其大小的选取应视实际所需分层的详细程度而定。

B 值越小， 允许误差)(i x E 越小，就越粗。

所分层就越细； 反之所分层就会越粗。

在这里本文人为的定义B=3，根据一号标准井可以得到每一层的允许误差)(i x E 。

根据层内差异法，1号井的分层范围为：
表3 一号井分层范围
分层 长31 长32 长33 长41 长42 长61 长62 长63 最小值 217.29 166.10 220.84 206.08 225.6 232.85 228.53 159.2 最大值 306.85 352.00 282.7 298.94 273.64 261.19 257.81 347.54 分层 长71 长72 长73 长81 长82 长91 长92 最小值 180.40 202.07 106.7 159.36 212.63 184.07 43.42 最大值 302.58 271.49 414.7 326.30 259.01 300.24 487.58 由层内分析法分层后分层过细，共统计了7个井15小层的分层数据522个,发现误差小于0. 5m 的分层数据为314 个,符合率为60%。

，于是对模型进行改进，由聚类分析法将分的过细的层进行合并，方法如下。

2、用层内聚类分析法将井进行细分层：
在上面的分层处理中,对于整个处理井段均采用了相同的B 值, 而且在精细分层中B 值取得较小,导致实际操作中可能会将某些层段划分过细,把原本同一地层分成了若干小层。

我们在层内差异法分层的基础上再结合聚类分析法, 以避免上述问题的出现。

聚类分析是根据某一分类统计量来度量多个样品(观察量) 间的相似程度, 依相似程度的高低决定
是聚合为一类, 还是划为不同类。

这里我们选用马哈拉诺比斯距离(简称“马氏”距离) 作为分类统计量来表示已划分出的各小层间的相似程度, 把相似程度高、且相邻的小层合并为一层。

具体作法如下:
若在经层内差异法分层的数条曲线中有L 条测井曲线纵向分辨率较强, 按前述方法共分地层总数为M ,这时计算第k 层与第k + 1 层的测井值之间的“马氏”距离：
∑=+-=
+l
i ik
k x x
k k d 1
2
1)()1,( 根据计算结果,若min )1,(d k k d ≤+,则认为第k 层与第k + 1 层属于同一小层，故进行并层处理;若min )1,(d k k d >+ 则不作并层处理。

其中min d 为各层间测井值的最小“临界距离”,可由经验选择。

层内差异法的流程图如下图：
图三 层内分析法的流程图
B 与 min d 的数值选取分析
1、B值和
d值的选取:非均值地层中的人机联作方式,实际上是一种人工干预,
min
参数是针对某一口井或某一地区而言的,具体在确定分层界限的有关经验规律上及分层处理技巧上。

也就是说B值和
d值是一个统计数据,适用于统计区的地
min
层特征,它虽然受岩性影响,但其大小取值对岩性而言并无规律性,即在非统计区不适用。

2、选取合适的d 、B 值后,一些井尚需做“并层”处理。

统计发现,多需删除曲线呈负斜率变化的分层,曲线由大变小处的分层结果。

这时若上、下二分层所夹地层厚度> (2. 5/ 3. 0) m ,则不再删除。

3、并层处理后若二分层间的地层厚度> (6. 0/6. 5) m ,可从其间分出另一层。

4、若二分层之间的地层厚度较小,如小于1 m时,可进行并层处理。

需要说明的是:上述技巧是在对本区小层厚度分布进行分析、对比之后提出的,仅在本区有效。

图四为：一号井的人工分层和自动分层对比图：（彩线为人工分层，红线为自动分层）
051015
图四一号井人工与自动分层图
其柱状图为：
123456789101112131415
图三一号井人工分层与自动分层的比较
利用此测井曲线进行自动分层的框图:
2号人工257 290 326 369 414 458 497.5 2号自动239.6 263.6 302.75 327.75 358.6 425 448
3号人工213 251 297.4 340.5 377.9 422 455
3号自动210.75 277 323.125 347.125 374.875 412.875 448
4号人工276 311 355.7 398 438 483 518
4号自动280.125 347.5 383 407.5 434 443 499
5号人工450.4 492 535.4 479.8 618 661.7 698
5号自动396.125 450.625 504.625 541.625 580.375 622.375 665
6号人工261.7 310 350 396 433.6 6号自动250.125 298.4 330.65 375.5 408.652 7号人工334 368 410 454 495 528 576
7号自动230.35 260.1 318.1 371.725 434.975 474.225 529
长63 长71 长72 长73 长81 长82 长91 长92 545.1 582 612 644.3 682 717 765 811 497.125 541.875 607 642.8 668 710.75 741 790.125 502 540.5 578 621 663 700 740788.9 501 553.75 586.125 608.625 661 696 747.375 784.15 560.3 599.7 643.5 680 720 760 803 846 518 560.125 612.125 656.625 730.25 796 834 860.125 741 780 821.9 860.4 900.4 934 972 1010.7 682 741.5 786 829.5 867.375 901.125 964.875 1000.5 478 514.5 558 602 641.2 682.1 717.6 763 439.775 477.525 508.275 561.4 601.4 643.15 664.125 701.275 615 652 690 733 770 810 860.6 901 560.225 615.85 642.1 676.975 730.35 773.35 853.85 872.35
长31长33长42长62长71长73长82长92图四二号井的人工分层与自动分层比较图
长31长33长42长62长71长73长82长92图五3号井的人工分层与自动分层比较图
长31长33长42长62长71长73长82长92图六4号井的人工分层与自动分层比较图
长31长33长42长62长71长73长82长92图七5号井的人工分层与自动分层比较图
长31长33长42长62长71长73长82长92图八6号井的人工分层与自动分层比较图
长31长33长42长62长71长73长82长92图九7号井的人工分层与自动分层比较图
对二、三、四号井进行误差分析：
表5 2号到4号井的误差分析
长31 长32 长33 长41 长42 长61 长62
二号井0.0677 0.0910 0.0713 0.1118 0.1338 0.0721 0.0995 三号井0.0106 0.1036 0.0865 0.0195 0.0080 0.00216 0.0154 四号井0.0149 0.1174 0.0768 0.0239 0.0091 0.0828 0.0367
长63 长71 长72 长73 长81 长82 长91 长92 0.0880 0.0689 0.0082 0.0023 0.0205 0.0087 0.0314 0.0257 0.0020 0.0245 0.0141 0.0199 0.0030 0.0057 0.0100 0.0060 0.0755 0.0660 0.0488 0.0344 0.0142 0.0474 0.0386 0.0167
对二到七号井的误差分析可以得出，误差不大于0.15，可见精度还是符合要求的。

本文用此改进后的模型对八到十三号井进行自动分层。

结果如下：
图6 8号到13号井的自动分层结果
长31 长32 长33 长41 长42 长61 长62
8号井96.5 132.5 172.625 230.75 269.625 327.5 360.125 9号井161.2 195.75 236.5 280.75 326.375 356.95 401.75 10号井72.5 138.25 175.75 242.75 280.625 342.625 389
11号井224.125 260.5 309.965 355.5 387 424.375 472.75 12号井181.5 233 282 337.75 381 441.25 491.375 13号井232.225 278.85 307.6 345.725 385.475 432.225 482.6
长63 长71 长72 长73 长81 长82 长91 长92 402.875 450.5 500 550.125 602.25 630 679 714 442.325 492.825 534.45 575.75 629.75 672.575 731.325 782.825 428.75 461.625 506.125 546.125 599.625 639.225 682.375 700 524.25 565.125 611.375 667.25 698.25 751 761.125 839.25 541.25 599.5 636.5 705 749.5 775.75 815 857.125 554.6 603.1 640.35 687.6 717.415 757.85 782.725 834.475
对结果进行分析：
1、8到13号井没有缺层现象。

2、对于同一层其深度范围并不相同，说明井的位置不同可能会导致每一层的层位范围不同，缺失现象也不相同。

五、模型的优缺点及改进
模型优点：
1、运用主成分分析法将数据进行了简化，方便了计算。

2、运用层内分析法与聚类分析法结合起来进行粗细分层，较好的对各井分层。

3、采用内部一致性原则，将关键变量构造成综合曲线，不仅能简化计算，而且
能够形成一条各层差异显著的曲线，对于后面的分层过程起到促进作用。

模型缺陷：
1、数据量较大，为了简化模型，可能没有有效利用数据。

2、依据标准井计算出来的允许误差以及“马氏”最小距离，由于人为经验因素，可能会对结果有一定的影响。

模型的改进：
由于本文对B以及马氏距离的确定是根据以往的经验人为确定的，所以误差是肯定有的。

为了减小这个误差，可以对模型再进行改进。

对允许误差的选定可以参考一号标准井的数据，先选定几个B，将结果与一号井进行对比，求出每个B的误差，以此方法，可以选出一号井每一层较为合适的B，以此作为其他井的标准，这样就可以将误差大大的减小。

参考文献
[1]程玉群、李秀荣、刘铁庄，测井曲线自动分层技术及在杏北地区小层划分中的应用，地质勘探，第7 卷第2 期
[2] /question/219318664.html?an=0&si=1
[3]雍世和、陈钢花、白康生，测井曲线自动分层，测井技术Well Logging Technology编辑部邮箱1987年06期
[4]周开利、邓春晖，Matlab基础及应用教程,北京大学出版社
[5]韩中庚，数学建模方法及应用，高等教育出版社
[6] 卓金武，Matlab在数学建模中的应用，北京航空航天大学出版社
附件
对选出来的指标进行主成分分析,其程序如下：
Clear
clc
load data2.txt
X=[data2];
Z=zeros(size(X));
Z(:,1)=(X(:,1)-mean(X(:,1)))/std(X(:,1));
Z(:,2)=(X(:,2)-mean(X(:,2)))/std(X(:,2));
Z(:,3)=(X(:,3)-mean(X(:,3)))/std(X(:,3));
Z(:,4)=(X(:,4)-mean(X(:,4)))/std(X(:,4));
Z(:,5)=(X(:,5)-mean(X(:,5)))/std(X(:,5));
Z(:,6)=(X(:,6)-mean(X(:,6)))/std(X(:,6));
Z(:,7)=(X(:,7)-mean(X(:,7)))/std(X(:,7));
Z;
R=Z'*Z/16;%相关矩阵R=Z'*Z/(n-1);
R; % 求相关矩阵R的特征值和特征向量
[eigenvectors, eigenvalues]=eig(R); % 调用matlab函数：[V,D]=eig(A)
eigenvalues %输出特征值矩阵
eigenvectors; %输出特征向量矩阵
a=[eigenvectors];
b=[eigenvalues];
a1=a(:,1);a2=a(:,2);a3=a(:,3);a4=a(:,4);a5=a(:,5);a6=a(:,6);a7=a(:,7);
b1=b(7,7);b2=b(6,6);b3=b(5,5);b4=b(4,4);b5=b(3,3);b6=b(2,2);b7=b(1,1);
bz=b1+b2+b3+b4+b5+b6+b7;
c1=b1/bz;c2=b2/bz;c3=b3/bz;c4=b4/bz; c5=b5/bz; c6=b6/bz; c7=b7/bz; %c(i)方差贡献率
Y=[c1;c2;c3;c4;c5;c6;c7];
d1=c1;d2=(c1+c2);d3=(c1+c2+c3);d4=(c1+c2+c3+c4);d5=(c1+c2+c3+c4+c5);
d6=(c1+c2+c3+c4+c5+c6); d7=(c1+c2+c3+c4+c5+c6+c7); %d(i)累计贡献率
U=[d1;d2;d3;d4;d5;d6;d7]
[pc,score,latent]=princomp(Z) %PC各主成分关于指标的线性组合的系数矩阵，score为各主成分得分，latent方差矩阵的特征值。