2.2 Tool Box函数的调用
Tool Box函数定义在以.m为扩展名的文件中,如regress.m就是一个工具箱函数,用来计算多元线性回归,函数体如下:
function [b,bint,r,rint,stats] = regress(y,X,alpha)
if nargin < 2,
error('REGRESS requires at least two input arguments.');
end
if nargin == 2,
alpha = 0.05;
end
% Check that matrix (X) and left hand side (y) have compatible dimensions
[n,p] = size(X);
[n1,collhs] = size(y);
if n ~= n1,
error('The number of rows in Y must equal the number of rows in X.');
end
if collhs ~= 1,
error('Y must be a vector, not a matrix');
end
% Remove missing values, if any
wasnan = (isnan(y) | any(isnan(X),2));
if (any(wasnan))
y(wasnan) = [];
X(wasnan,:) = [];
n = length(y);
end
% Find the least squares solution.
[Q, R]=qr(X,0);
b = R\(Q'*y);
% Find a confidence interval for each component of x
% Draper and Smith, equation 2.6.15, page 94
RI = R\eye(p);
xdiag=sqrt(sum((RI .* RI)',1))';
nu = n-p; % Residual degrees of freedom
yhat = X*b; % Predicted responses at each data point.
r = y-yhat; % Residuals.
if nu ~= 0
rmse = norm(r)/sqrt(nu); % Root mean square error.
else
rmse = Inf;
end
s2 = rmse^2; % Estimator of error variance.
tval = tinv((1-alpha/2),nu);
bint = [b-tval*xdiag*rmse, b+tval*xdiag*rmse];
% Calculate R-squared.
if nargout==5,
RSS = norm(yhat-mean(y))^2; % Regression sum of squares.
TSS = norm(y-mean(y))^2; % Total sum of squares.
r2 = RSS/TSS; % R-square statistic.
if (p>1)
F = (RSS/(p-1))/s2; % F statistic for regression
else
F = NaN;
end
prob = 1 - fcdf(F,p-1,nu); % Significance probability for regression
stats = [r2 F prob];
% All that requires a constant. Do we have one?
if (~any(all(X==1)))
% Apparently not, but look for an implied constant.
b0 = R\(Q'*ones(n,1));
if (sum(abs(1-X*b0))>n*sqrt(eps))
warning(sprintf(['R-square is not well defined unless X has' ...
' a column of ones.\nType "help regress" for' ...
' more information.']));
end
end
end
% Find the standard errors of the residuals.
% Get the diagonal elements of the "Hat" matrix.
% Calculate the variance estimate obtained by removing each case (i.e. sigmai)
% see Chatterjee and Hadi p. 380 equation 14.
T = X*RI;
hatdiag=sum((T .* T)',1)';
ok = ((1-hatdiag) > sqrt(eps));
hatdiag(~ok) = 1;
if nu < 1,
ser=rmse*ones(length(y),1);
elseif nu > 1
denom = (nu-1) .* (1-hatdiag);
sigmai = zeros(length(denom),1);
sigmai(ok) = sqrt((nu*s2/(nu-1)) - (r(ok) .^2 ./ denom(ok)));
ser = sqrt(1-hatdiag) .* sigmai;
ser(~ok) = Inf;
elseif nu == 1
ser = sqrt(1-hatdiag) .* rmse;
ser(~ok) = Inf;
end
% Create confidence intervals for residuals.
Z=[(r-tval*ser) (r+tval*ser)]';
rint=Z';
% Restore NaN so inputs and outputs conform
if (nargout>2 & any(wasnan))
tmp = ones(size(wasnan));
tmp(:) = NaN;
tmp(~wasnan) = r;
r = tmp;
end
if (nargout>3 & any(wasnan))
tmp = ones(length(wasnan),2);
tmp(:) = NaN;
tmp(~wasnan,:) = rint;
rint = tmp;
end
由于这个函数是由MATLAB专有的语言写成,而并非C语言,所以我们并不能直接使用在我们的程序中。不过没关系MATLAB为我们提供了强大的工具mcc.exe,它能够将.m文件中的函数翻译成我们需要的C/C++程序,这个功能令我们的开发事半功倍。
下面具体阐述一下mcc的使用方法。
1.设置一个有效的C编译器,这里面使用MATLAB自带的编译器Lcc C version 2.4。
开启MATLAB环境,在Command window里面输入
>> mbuild –setup
Please choose your compiler for building standalone MATLAB applications:
Would you like mbuild to locate installed compilers [y]/n? y
Select a compiler:
[1] Lcc C version 2.4 in D:\MATLAB\sys\lcc
[2] Microsoft Visual C/C++ version 6.0 in C:\Program Files\Microsoft Visual Studio
[0] None
Compiler: 1
Please verify your choices:
Compiler: Lcc C 2.4
Location: D:\MATLAB\sys\lcc
Are these correct?([y]/n): y
The default options file:
"C:\Documents and Settings\Administrator\Application Data\MathWorks\MATLAB\R12\compopts.bat"
is being updated from D:\MATLAB\BIN\WIN32\mbuildopts\lcccompp.bat...
现在我们已经将位于D:\MATLAB\sys\lcc的一个c/c++编译器设为我们要使用的编译器。
2.现在可以使用mcc来生成C文件了,这里我们使用regress.m这个文件作为例子,生成相应的c语言文件。在Command Window 里面输入
>> mcc -m regress
执行完毕后,你打开<matlab>\work目录,你会发现目录中有了我们需要的文件regress.h和regress.c以及很多其他C文件和C的头文件,这些其他文件都是regress函数里面用到的其它工具箱函数所对应的C语言文件,难以想象mcc为我们做了多么繁琐而复杂的事情呀!
在Command Window 里面输入
>> mcc –p regress
执行完毕后,就能在work目录下生成对应的CPP文件。
3.新建一个控制台程序,主程序代码如下:
#pragma hdrstop
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "libmatlb.h"
#include "libmmfile.h"
#include "regress.h"
#include "fcdf.h"
#include "tinv.h"
#include "chi2cdf.h"
#include "distchck.h"
#include "betainv.h"
#include "norminv.h"
#include "gamcdf.h"
#include "betacdf.h"
#include "betapdf.h"
#pragma argsused
static _mexInitTermTableEntry init_term_table[11]
= { { libmmfileInitialize, libmmfileTerminate },
{ InitializeModule_regress, TerminateModule_regress },
{ InitializeModule_fcdf, TerminateModule_fcdf },
{ InitializeModule_tinv, TerminateModule_tinv },
{ InitializeModule_chi2cdf, TerminateModule_chi2cdf },
{ InitializeModule_distchck, TerminateModule_distchck },
{ InitializeModule_betainv, TerminateModule_betainv },
{ InitializeModule_norminv, TerminateModule_norminv },
{ InitializeModule_gamcdf, TerminateModule_gamcdf },
{ InitializeModule_betacdf, TerminateModule_betacdf },
{ InitializeModule_betapdf, TerminateModule_betapdf } };
static double buf1[] = { 1, 2 ,3,3,2,1,3,3,4};
static double buf2[] = { 3, 4 ,4 };
int main(int argc, const char * * argv) {
mxArray *x = NULL;
mxArray *y = NULL;
mxArray *alpha = NULL;
mxArray *bint = NULL, *r = NULL, *rint = NULL;
mxArray *status = NULL,*b=NULL;
/* Enable automated memory management */
mlfEnterNewContext(0, 0);
/* Create the matrices and assign data to them */
mlfAssign(&x, mlfDoubleMatrix(3, 3, buf1, NULL));
mlfAssign(&y, mlfDoubleMatrix(3, 1, buf2, NULL));
mlfAssign(&alpha,mlfScalar(1));
/* Print the matrices */
mlfPrintMatrix(x);
mlfPrintMatrix(y);
int i;
for(i=0;i<sizeof(init_term_table)/sizeof(_mexInitTermTableEntry);i++)
{
init_term_table[i].initialize();
}
b=mlfRegress(&bint,&r,&rint,&status,y,x,alpha);
//b=mlfRegress(NULL,NULL,NULL,NULL,mat0,mat1,alpha);
mlfPrintMatrix(b);
mlfPrintMatrix(bint);
mlfPrintMatrix(r);
mlfPrintMatrix(rint);
mlfPrintMatrix(status);
/* Free the matrices */
mxDestroyArray(x);
mxDestroyArray(y);
for(i=0;i<sizeof(init_term_table)/sizeof(_mexInitTermTableEntry);i++)
{
init_term_table[i].terminate();
}
/* Disable automated memory management */
mlfRestorePreviousContext(0, 0);
return(EXIT_SUCCESS);
}
2.3 调用Tool Box函数的另一个可供选择的方法
2.2节提出的调用ToolBox函数的方法比较繁杂,因为MATLAT编译器会为你生成一大堆C/C++源代码文件和头文件,对这些源代码的维护也就比较复杂,随着系统的不断发展,程序员对这些MATLAB生成的源代码的维护将会变得极其困难,在此提出另一个可供选择的方法:将ToolBox函数编译成动态链接库。对程序员来说,相比之下二进制级别的重用要比源代码级的重用更易于维护和使用。
那么到底怎么把ToolBox函数编译成动态链接库呢?庆幸的是MATLAT编译器能为我们来完成这项任务!下面我来讲述一下。
首先正确设置MCC(对MCC的使用见2.2章)。然后输入如下命令:
mcc -t -W lib:test -T link:lib betaln erf erfcore realmax betacdf betapdf erfinv gammaln betacore erfc gammainc gamcdf distchck betainv norminv betainc chi2cdf tinv mean fcdf regress
执行这条命令就会生成一个test.dll、test.lib以及dll中导出函数的声明头文件。有了这些文件就可以轻松调用DLL中的函数了。
注意
①这条命令的主要目的是要生成regress这个函数,但是由于计算regress时需要调用到那些子函数如betaln,erf,erfcore等等,所以才会写这么一大串,如果仅仅执行
mcc -t -W lib:test -T link:lib regress
MATLAT编译器将会给出错误信息。
②编译器生成的LIB文件是coff格式的,如果你要在bcb中静态加载DLL,则需要把coff格式转变成omf格式。可以是用coff2omf.exe这个工具,或者索性使用IMPLIB.exe重新生成一个LIB文件。关于这两个工具的使用方法可以参考它们的帮助。