How to manually perform linear regression from scratch

My code is below. Firstly, is my linear regression actually correct?

If it's right, I want to do the same but I want to create an original code for doing the linear regression, fitting a straight line to the data.

clc
clear
%input tensile strength data
TS=[72.5 73.8 68.1 77.9 65.5 73.23 71.17 79.92 65.67 74.28 67.95 82.84...
    79.83 80.52 70.65 72.85 77.81 72.29 75.78 67.03 72.85 77.81 75.33...
    71.75 72.28 79.08 71.04 67.84 69.2 71.53]';

%calculate mean
meanTS=mean(TS)                                                             %using matlab function
meanTSmanual=(sum(TS))/length(TS)                                           %calculate mean manually

%calculate standard deviation
stdTS=std(TS)                                                               %using matlab function
%calculate standard deviation manually
meandifference = (TS-meanTSmanual).^2;                                      %find the difference between each point of the mean then square it
TSvariance=sum(meandifference)/(length(TS)-1);                              %divide by sample size minus 1 to get the variance
stdTSmanual=sqrt(TSvariance)                                                %square root the variance to get the standard deviation

%ascending order
[ascend,index]=sort(TS,'ascend');                                           %sort the data in ascending order
ascend

%probability of survival
survival=(1-(index/(length(TS)-1)))                                         %applying given formula

%Weibull model
reg1=real(log((log(1./survival))))                                         %first dataset for linear regression, real parts only
reg2=log(TS)                                                                %second dataset for linear regression

%linear regression
idx = isfinite(reg1) & isfinite(reg2);                                      %exclude infinite values and NaNs
pwf = fit(reg1(idx),reg2(idx),'poly1')
plot(pwf,reg1,reg2)
xlim([-5,5])
ylim([-5,5])

Answer 1

[Torsten]

This is the code for the manual fit:

A = [ones(size(reg1(idx))),reg1(idx)];
b = reg2(idx);
% linear fit equals reg2(idx)_fitted = coeffs(1) + coeffs(2)*reg1(idx)
coeffs = A\b

[Daniel Jordan]

Hi, thanks for the reply. Below is the graph that I get using this code. Figure 1 id the original matlab code, figure 2 is the one you've given. It looks like the line is trying to go directly through each point rather than a line of best fit?

[Torsten]

As written in the comment

% linear fit equals reg2(idx)_fitted = coeffs(1) + coeffs(2)*reg1(idx)

you must plot coeffs(1)+coeffs(2)*reg1(idx) against reg1(idx):

plot(reg1(idx),coeffs(1)+coeffs(2)*reg1(idx))