11.zip_manabh_quasi newton _solveuld_共轭梯度法_拟牛顿法

x = load('ex3x.dat'); y = load('ex3y.dat'); trustRegionBound = 1000; x = [ones(size(x,1),1) x]; meanx = mean(x);%求均值 sigmax = std(x);%求标准偏差 x(:,2) = (x(:,2)-meanx(2))./sigmax(2); x(:,3) = (x(:,3)-meanx(3))./sigmax(3); itera_num = 1000; %尝试的迭代次数 sample_num = size(x,1); %训练样本的次数 figure alpha = [0.01, 0.03, 0.1, 0.3, 1, 1.3];%因为差不多是选取每个3倍的学习率来测试，所以直接枚举出来 plotstyle = {'b', 'r', 'g', 'k', 'b--', 'r--'}; theta_grad_descent = zeros(size(x(1,:))); %% 信赖域+狗腿法 theta = zeros(size(x,2),1); %theta的初始值赋值为0 Jtheta = zeros(itera_num, 1); for i = 1:itera_num %计算出某个学习速率alpha下迭代itera_num次数后的参数 Jtheta(i) = (1/(2*sample_num)).*(x*theta-y)'*(x*theta-y);%Jtheta是个行向量 grad = (1/sample_num).*x'*(x*theta-y); B=x'*x; du = -grad' * grad * grad / (grad' * B * grad); dB = -B^-1 * grad; a = 2; if du'*du > trustRegionBound*trustRegionBound; a = trustRegionBound / sqrt((du'*du)); else if dB'*dB > trustRegionBound*trustRegionBound a = sqrt((trustRegionBound*trustRegionBound - du'*du) / ((dB-du)'*(dB-du))) + 1; end end if a < 1 d = a * du; else d = du + (a - 1) * (dB - du); end Jtheta1(i)=(1/(2*sample_num)).*(x*(theta+d)-y)'*(x*(theta+d)-y); p = (Jtheta(i)-Jtheta1(i))/(-grad'*d-1/2*d'*B*d); if p > 0.75 && sqrt(abs(d'*d) - trustRegionBound) < 0.001 trustRegionBound = min(2 * trustRegionBound, 10000); else if p < 0.25 trustRegionBound = sqrt(abs(d'*d)) * 0.25; end end if p > 0%q(zeros(2,1),x) > q(d, x) theta = theta + d; end end K(1)=Jtheta(1000) plot(0:350, Jtheta(1:351),'k--','LineWidth', 2)%此处一定要通过char函数来转换 hold on %% 固定学习率法 theta_grad_descent = zeros(size(x(1,:))); for alpha_i = 1:length(alpha) %尝试看哪个学习速率最好 theta = zeros(size(x,2),1); %theta的初始值赋值为0 Jtheta = zeros(itera_num, 1); for i = 1:itera_num %计算出某个学习速率alpha下迭代itera_num次数后的参数 Jtheta(i) = (1/(2*sample_num)).*(x*theta-y)'*(x*theta-y);%Jtheta是个行向量 grad = (1/sample_num).*x'*(x*theta-y); theta = theta - alpha(alpha_i).*grad; end K(alpha_i+1)=Jtheta(1000); plot(0:350, Jtheta(1:351),char(plotstyle(alpha_i)),'LineWidth', 2)%此处一定要通过char函数来转换 hold on end %% %% SD算法 theta = zeros(size(x,2),1); %theta的初始值赋值为0 Jtheta = zeros(itera_num, 1); for i = 1:itera_num %计算出某个学习速率alpha下迭代itera_num次数后的参数 Jtheta(i) = (1/(2*sample_num)).*(x*theta-y)'*(x*theta-y);%Jtheta是个行向量 grad = (1/sample_num).*x'*(x*theta-y); Q=x'*x; d=-grad; a=(grad'*grad)/(grad'*Q*grad); theta = theta + a*d; end K(1)=Jtheta(1000) plot(0:350, Jtheta(1:351),'b--','LineWidth', 2); hold on %% legend('Trust Region with DogLeg','0.01','0.03','0.1','0.3','1','1.3','Steepest Descent'); xlabel('Number of iterations') ylabel('Cost function') figure plot(1:7,K,'b-','LineWidth', 2);