# Scipy Method

``````

In [13]:

import numpy as np
from scipy.optimize import curve_fit
%matplotlib inline
import matplotlib.pyplot as plt

``````
``````

In [14]:

# Define function to be used
def T2decay(xdata, T1,Mz):
R1=1.0/T1
return Mz * np.exp(-xdata*R1)

def plot(xdata,ydata1,ydata2):
plt.plot(xdata,ydata1)
plt.scatter(xdata,ydata2,color='r')
plt.show()

``````
``````

In [15]:

# Simulate data
NPoints=20
TR=np.linspace(0,0.2,NPoints)
T2=.050
Mz=1.0
Signal=T2decay(TR, T2,Mz)

``````
``````

In [16]:

mu=0.0; sigma=0.1;
Signal_with_noise= Signal + np.random.normal(mu, sigma, NPoints)

``````
``````

In [17]:

# Plot
plot(TR,Signal,Signal_with_noise)

``````
``````

``````
``````

In [18]:

# Perform Non-Linear curve fitting with bounds
# using SciPy
# https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.html
initial_guess=[1,.1]
pars_hat, cov = curve_fit(T2decay, TR, Signal_with_noise,p0=initial_guess,bounds=(0, 1.0))

``````
``````

In [19]:

# estimate Yhat
Yhat = T2decay(TR, pars_hat[0],pars_hat[1])

``````
``````

In [20]:

# Plot estimated data
plot(TR,Signal,Yhat)

``````
``````

``````

# LMFIT Method

``````

In [21]:

# Remember to Install if you don't have it
from lmfit import minimize, Parameters

``````
``````

In [22]:

# LMFIT requires a function that returns the objective function. In this case the sum of the square of the residuals
def residual(params, xdata, ydata_observed):
R2 = 1.0/params['T2']
Mz = params['Mz']
model = Mz * np.exp(-xdata*R2)
return (ydata_observed-model)

params = Parameters()

``````
``````

In [40]:

# run fitting
lmfit = minimize(residual, params, args=(TR, Signal_with_noise))

``````
``````

In [49]:

# estimate Yhat
T2pred=lmfit.params['T2'].value
Mzpred=lmfit.params['Mz'].value
Yhat2 = T2decay(TR, T2pred,Mzpred)

``````
``````

In [50]:

# Plot estimated data
plot(TR,Signal,Yhat2)

``````
``````

``````
``````

In [ ]:

``````