利用MATLAB实现以下功能，从外部导入10个合约数据：第一列为内在价值（标的资产价格-执行价），第二列为卖价，第三列为最新价（看作市场价），第四列为执行价，10个合约的到期时间T均为21/365。根据卖价使用Metropolis-Hastings算法对波动率进行反演，要求使用高斯分布作为建议分布函数，并根据已知的标的资产价格、执行价和到期时间，将反演得到的波动率带回B-S公式，得到期权的理论价格，将理论价格与最新价和卖价进行比较，并输出图像

以下是利用MATLAB实现上述功能的代码：

% 从外部导入10个合约数据
data = csvread('contracts.csv');
% 分别提取内在价值、卖价、最新价和执行价
intrinsic_value = data(:,1);
sell_price = data(:,2);
market_price = data(:,3);
strike_price = data(:,4);
% 设置到期时间
T = 21/365;

% 定义波动率反演函数
function [sigma,accept_rate] = metropolis_hastings(intrinsic_value,sell_price,market_price,strike_price,T)
    % 初始化波动率和建议分布标准差
    sigma = 0.2;
    proposal_std = 0.01;
    % 初始化接受率和采样次数
    accept_rate = 0;
    sample_size = 10000;
    % 初始化存储采样结果的向量
    sigma_samples = zeros(sample_size,1);
    % 进行采样
    for i = 1:sample_size
        % 从建议分布中抽取样本
        sigma_star = sigma + proposal_std * randn();
        % 计算接受率
        alpha = min(1,exp(log_likelihood(intrinsic_value,sell_price,market_price,strike_price,T,sigma_star) - log_likelihood(intrinsic_value,sell_price,market_price,strike_price,T,sigma)));
        % 判断是否接受
        u = rand();
        if u < alpha
            sigma = sigma_star;
            accept_rate = accept_rate + 1/sample_size;
        end
        sigma_samples(i) = sigma;
    end
end

% 定义对数似然函数
function log_likelihood_value = log_likelihood(intrinsic_value,sell_price,market_price,strike_price,T,sigma)
    d1 = (log(market_price./strike_price) + (sigma.^2/2).*T)./sigma./sqrt(T);
    d2 = d1 - sigma.*sqrt(T);
    log_likelihood_value = -0.5*log(2*pi) - log(sigma) - T.*intrinsic_value.^2./sigma.^2/2 - log(normcdf(-d2)) - log(normcdf(d1-intrinsic_value./sigma./sqrt(T)));
    log_likelihood_value = sum(log_likelihood_value(isfinite(log_likelihood_value)));
end

% 对波动率进行反演
[sigma,accept_rate] = metropolis_hastings(intrinsic_value,sell_price,market_price,strike_price,T);

% 根据反演得到的波动率计算期权理论价格
d1 = (log(market_price./strike_price) + (sigma.^2/2).*T)./sigma./sqrt(T);
d2 = d1 - sigma.*sqrt(T);
theory_price = market_price.*normcdf(d1) - strike_price.*exp(-T.*sigma.^2/2).*normcdf(d2);

% 绘制图像
figure;
plot([sell_price,theory_price,market_price]);
legend('Sell Price','Theory Price','Market Price');
xlabel('Contract Index');
ylabel('Option Price');

需要注意的是，在运行代码之前，需要将合约数据保存在当前目录下的名为“contracts.csv”的文件中。另外，由于Metropolis-Hastings算法的采样结果具有随机性，因此每次运行得到的结果可能会有所不同。