NumPy

# Using pip
pip install numpy

# Using conda
conda install numpy

# Verify installation
import numpy as np
print(np.__version__)

import numpy as np  # Standard alias

import numpy as np

# From Python list
arr1 = np.array([1, 2, 3, 4, 5])
print(arr1)  # [1 2 3 4 5]

# 2D array
arr2 = np.array([[1, 2, 3], [4, 5, 6]])
print(arr2)
# [[1 2 3]
#  [4 5 6]]

# Specify data type
arr3 = np.array([1, 2, 3], dtype=np.float64)
print(arr3)  # [1. 2. 3.]

# 0D array (scalar)
scalar = np.array(42)

# 1D array (vector)
vector = np.array([1, 2, 3])

# 2D array (matrix)
matrix = np.array([[1, 2], [3, 4]])

# 3D array (tensor)
tensor = np.array([[[1, 2], [3, 4]], [[5, 6], [7, 8]]])

# Create arrays of zeros
zeros_1d = np.zeros(5)           # [0. 0. 0. 0. 0.]
zeros_2d = np.zeros((3, 4))      # 3x4 matrix of zeros
zeros_3d = np.zeros((2, 3, 4))   # 2x3x4 tensor

# Create arrays of ones
ones_2d = np.ones((2, 3))        # 2x3 matrix of ones

# Create arrays filled with specific value
full = np.full((3, 3), 7)        # 3x3 matrix filled with 7

# Identity matrix
identity = np.identity(4)        # 4x4 identity matrix
eye = np.eye(3, k=1)             # 3x3 with ones on kth diagonal
# [[0. 1. 0.]
#  [0. 0. 1.]
#  [0. 0. 0.]]

# Range (similar to Python range)
range_arr = np.arange(10)        # [0 1 2 3 4 5 6 7 8 9]
range_arr = np.arange(2, 10, 2)  # [2 4 6 8]

# Linspace (evenly spaced numbers)
linspace_arr = np.linspace(0, 1, 5)  # [0.   0.25 0.5  0.75 1.  ]

# Logspace (logarithmically spaced)
logspace_arr = np.logspace(0, 2, 5)  # [1. 3.16 10. 31.6 100.]

# Random values from uniform distribution [0, 1)
random_uniform = np.random.random((3, 3))

# Random integers
random_ints = np.random.randint(0, 10, size=(3, 4))

# Random values from normal distribution
random_normal = np.random.randn(3, 3)

# Empty array (uninitialized, contains garbage values)
empty_arr = np.empty((2, 3))

# Array of zeros with same shape as another array
arr = np.array([[1, 2], [3, 4]])
zeros_like = np.zeros_like(arr)

# Copy of array
copy_arr = arr.copy()

# From range
from_range = np.fromfunction(lambda i, j: i + j, (3, 3))

arr = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])

# Shape (dimensions)
print(arr.shape)        # (3, 3)
print(arr.ndim)         # 2 (number of dimensions)
print(arr.size)         # 9 (total number of elements)

# Data type
print(arr.dtype)        # int64 (depends on system)

# Memory usage
print(arr.itemsize)     # 8 (bytes per element)
print(arr.nbytes)       # 72 (total bytes)

# Size of each dimension
print(len(arr))         # 3 (length of first dimension)
print(arr.shape[0])     # 3

arr = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])

# Single element
print(arr[0, 0])        # 1
print(arr[2, 1])        # 8
print(arr[-1, -1])      # 9 (negative indexing)

# Row access
print(arr[0])           # [1 2 3]
print(arr[1:])          # [[4 5 6] [7 8 9]]

# Column access
print(arr[:, 0])        # [1 4 7]
print(arr[:, 1:])       # [[2 3] [5 6] [8 9]]

arr = np.arange(10)

# Basic slicing
print(arr[2:7])         # [2 3 4 5 6]
print(arr[:5])          # [0 1 2 3 4]
print(arr[5:])          # [5 6 7 8 9]
print(arr[::2])         # [0 2 4 6 8]
print(arr[::-1])        # [9 8 7 6 5 4 3 2 1 0]

# Multi-dimensional slicing
arr2d = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
print(arr2d[:2, 1:])    # [[2 3] [5 6]]
print(arr2d[1:, :2])    # [[4 5] [7 8]]

arr = np.array([10, 20, 30, 40, 50])

# Integer array indexing
indices = np.array([0, 2, 4])
print(arr[indices])     # [10 30 50]

# Boolean indexing
mask = arr > 30
print(arr[mask])        # [40 50]

# Conditional indexing
print(arr[arr % 20 == 0])  # [20 40]

# Fancy indexing
arr2d = np.array([[1, 2], [3, 4], [5, 6]])
row_indices = np.array([0, 2])
col_indices = np.array([1, 0])
print(arr2d[row_indices, col_indices])  # [2 5]

arr = np.array([1, 2, 3, 4, 5])

# In-place modification
arr[0] = 10
arr[1:3] = [20, 30]

# Broadcasting in assignment
arr[3:] = 0             # [10 20 30 0 0]

# Boolean assignment
arr[arr < 20] = 0       # Elements less than 20 become 0

arr = np.arange(12)

# Reshape to different dimensions
reshaped = arr.reshape(3, 4)
# [[ 0  1  2  3]
#  [ 4  5  6  7]
#  [ 8  9 10 11]]

# Reshape with -1 (automatic dimension)
reshaped_auto = arr.reshape(2, -1)  # 2x6
reshaped_auto = arr.reshape(-1, 3)  # 4x3

# Flatten
flat = arr.flatten()                # Returns copy
raveled = arr.ravel()               # Returns view (if possible)

# Transpose
matrix = np.array([[1, 2], [3, 4]])
transposed = matrix.T
# [[1 3]
#  [2 4]]

# 1D arrays
a = np.array([1, 2, 3])
b = np.array([4, 5, 6])
c = np.concatenate([a, b])          # [1 2 3 4 5 6]

# 2D arrays (vertical/horizontal)
a = np.array([[1, 2], [3, 4]])
b = np.array([[5, 6], [7, 8]])

vstack = np.vstack([a, b])          # Vertical stack
# [[1 2]
#  [3 4]
#  [5 6]
#  [7 8]]

hstack = np.hstack([a, b])          # Horizontal stack
# [[1 2 5 6]
#  [3 4 7 8]]

# Stack along new axis
dstack = np.stack([a, b], axis=2)   # 3D array

arr = np.arange(10)

# Split into equal parts
split_arr = np.split(arr, 5)        # 5 arrays of length 2

# Split at specific indices
split_at = np.split(arr, [3, 7])    # [0 1 2], [3 4 5 6], [7 8 9]

# 2D splitting
arr2d = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
vsplit = np.vsplit(arr2d, 3)        # Split by rows
hsplit = np.hsplit(arr2d, 3)        # Split by columns

arr = np.array([1, 2, 3])

# Add new axis
expanded = arr[np.newaxis, :]       # [[1 2 3]]
expanded = arr[:, np.newaxis]       # [[1] [2] [3]]

# Using np.expand_dims
expanded = np.expand_dims(arr, axis=0)

# Remove singleton dimensions
arr = np.array([[1], [2], [3]])
squeezed = np.squeeze(arr)          # [1 2 3]

a = np.array([1, 2, 3])
b = np.array([4, 5, 6])

# Arithmetic operations
print(a + b)        # [5 7 9]
print(a - b)        # [-3 -3 -3]
print(a * b)        # [4 10 18]
print(a / b)        # [0.25 0.4  0.5 ]
print(a ** b)       # [1 32 729]
print(a % b)        # [1 2 3]

# Comparison operations
print(a > 2)        # [False False  True]
print(a == b)       # [False False False]
print(a < b)        # [True True True]

# Logical operations
print(np.logical_and(a > 1, b < 6))  # [False True True]
print(np.logical_or(a < 2, b > 5))   # [True True False]

arr = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])

# Global statistics
print(np.mean(arr))         # 5.0
print(np.median(arr))       # 5.0
print(np.std(arr))          # 2.581988897
print(np.var(arr))          # 6.666666667
print(np.min(arr))          # 1
print(np.max(arr))          # 9
print(np.sum(arr))          # 45

# Axis-specific statistics
print(np.mean(arr, axis=0)) # [4. 5. 6.] (column means)
print(np.mean(arr, axis=1)) # [2. 5. 8.] (row means)
print(np.sum(arr, axis=0))  # [12 15 18] (column sums)
print(np.sum(arr, axis=1))  # [6 15 24] (row sums)

# Other statistical functions
print(np.percentile(arr, 50))  # Median (50th percentile)
print(np.percentile(arr, [25, 75]))  # Quartiles

arr = np.array([1, 2, 3, 4, 5])

# Trigonometric functions
print(np.sin(arr))          # [0.84147098 0.90929743 0.14112001 ...]
print(np.cos(arr))          # [0.54030231 -0.41614684 -0.9899925 ...]
print(np.tan(arr))          # [1.55740772 -2.18503986 -0.14254654 ...]

# Exponential and logarithmic
print(np.exp(arr))          # [2.71828183 7.3890561  20.08553692 ...]
print(np.log(arr))          # [0. 0.69314718 1.09861229 1.38629436 ...]
print(np.log10(arr))        # [0. 0.30103 0.47712125 0.60205999 ...]
print(np.sqrt(arr))         # [1. 1.41421356 1.73205081 2. 2.23606798]

# Absolute and rounding
print(np.abs([-1, -2, -3])) # [1 2 3]
print(np.round([1.2, 1.5, 1.8]))  # [1. 2. 2.]
print(np.floor([1.2, 1.8])) # [1. 1.]
print(np.ceil([1.2, 1.8]))  # [2. 2.]

# Scalar with array
arr = np.array([1, 2, 3, 4, 5])
print(arr + 10)             # [11 12 13 14 15]
print(arr * 2)              # [2 4 6 8 10]

# 1D with 2D
a = np.array([[1, 2, 3], [4, 5, 6]])
b = np.array([10, 20, 30])
print(a + b)
# [[11 22 33]
#  [14 25 36]]

# Column vector with row vector
a = np.array([[1], [2], [3]])
b = np.array([10, 20, 30])
print(a + b)
# [[11 21 31]
#  [12 22 32]
#  [13 23 33]]

# Normalization by row/column
data = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
row_means = np.mean(data, axis=1, keepdims=True)
normalized_rows = data - row_means

col_means = np.mean(data, axis=0)
normalized_cols = data - col_means

# Scaling
scaling_factor = np.array([0.1, 0.5, 0.9])
scaled_data = data * scaling_factor

# Arithmetic
np.add(a, b)                # a + b
np.subtract(a, b)           # a - b
np.multiply(a, b)           # a * b
np.divide(a, b)             # a / b
np.power(a, b)              # a ** b

# Comparison
np.greater(a, b)            # a > b
np.less(a, b)               # a < b
np.equal(a, b)              # a == b
np.maximum(a, b)            # element-wise max
np.minimum(a, b)            # element-wise min

# Special functions
np.expm1(x)                 # exp(x) - 1 (more accurate for small x)
np.log1p(x)                 # log(1 + x) (more accurate for small x)
np.sinc(x)                  # sinc function

# Using np.vectorize
def my_func(x):
    return x**2 + 2*x + 1

vectorized = np.vectorize(my_func)
arr = np.array([1, 2, 3])
result = vectorized(arr)    # [4 9 16]

# Using frompyfunc
def multiply_and_add(x, y):
    return x * y + 1

ufunc = np.frompyfunc(multiply_and_add, 2, 1)
result = ufunc(arr, 10)     # [11 21 31]

arr = np.array([1, 2, 3, 4, 5])

# Accumulate (prefix sums)
print(np.add.accumulate(arr))      # [1 3 6 10 15]
print(np.multiply.accumulate(arr)) # [1 2 6 24 120]

# Reduce (single result)
print(np.add.reduce(arr))          # 15 (sum)
print(np.multiply.reduce(arr))     # 120 (product)

# Outer product
a = np.array([1, 2, 3])
b = np.array([4, 5, 6])
print(np.add.outer(a, b))
# [[5 6 7]
#  [6 7 8]
#  [7 8 9]]

import numpy.linalg as la

A = np.array([[1, 2], [3, 4]])
B = np.array([[5, 6], [7, 8]])

# Matrix multiplication (dot product)
dot_product = np.dot(A, B)          # [[19 22] [43 50]]
matrix_product = A @ B               # Same as np.dot

# Inner product (1D arrays)
v = np.array([1, 2, 3])
w = np.array([4, 5, 6])
inner = np.inner(v, w)              # 32 (1*4 + 2*5 + 3*6)

# Cross product (3D vectors)
cross = np.cross([1, 0, 0], [0, 1, 0])  # [0 0 1]

# Matrix inverse
A_inv = la.inv(A)                   # [[-2.  1.] [1.5 -0.5]]
identity = A @ A_inv                # Identity matrix

# Determinant
det = la.det(A)                     # -2.0

# Trace
trace = np.trace(A)                 # 5 (sum of diagonal)

# Solve Ax = b
A = np.array([[3, 1], [1, 2]])
b = np.array([9, 8])

x = la.solve(A, b)                  # [2. 3.]
print(f"Solution: {x}")
print(f"Check: {A @ x}")            # Should equal b

# Least squares solution
A = np.array([[1, 1], [1, 2], [1, 3]])
b = np.array([1, 2, 3])
x, residuals, rank, s = la.lstsq(A, b, rcond=None)

A = np.array([[3, 1], [1, 2]])

# Eigenvalues and eigenvectors
eigenvalues, eigenvectors = la.eig(A)
print(f"Eigenvalues: {eigenvalues}")        # [4. 1.]
print(f"Eigenvectors: {eigenvectors}")      # Columns are eigenvectors

# Singular Value Decomposition (SVD)
U, s, Vt = la.svd(A)
print(f"U: {U}")
print(f"S: {s}")                     # Singular values
print(f"Vt: {Vt}")

# Norms
v = np.array([1, 2, 3])
print(la.norm(v))                   # L2 norm: 3.74165738677
print(la.norm(v, ord=1))            # L1 norm: 6

# Matrix norms
A = np.array([[1, 2], [3, 4]])
print(la.norm(A, ord='fro'))        # Frobenius norm: 5.47722557505

# Rank
rank = la.matrix_rank(A)            # 2

# Condition number
cond = la.cond(A)                   # 14.9330343737

# Moore-Penrose pseudoinverse
pinv = la.pinv(A)

# Random numbers in [0, 1)
uniform = np.random.rand(3, 4)      # 3x4 uniform distribution

# Standard normal (mean=0, std=1)
normal = np.random.randn(3, 4)      # 3x4 standard normal

# Random integers
integers = np.random.randint(0, 10, size=(3, 4))

# Random samples from different distributions
beta_sample = np.random.beta(2, 5, 1000)
gamma_sample = np.random.gamma(2, 2, 1000)
poisson_sample = np.random.poisson(5, 1000)
binomial_sample = np.random.binomial(10, 0.5, 1000)

# Random choice
choices = np.random.choice([1, 2, 3, 4, 5], size=10, replace=True)

# Shuffle in-place
arr = np.arange(10)
np.random.shuffle(arr)

# Permutation (returns shuffled copy)
permuted = np.random.permutation(10)

# Create a random generator
rng = np.random.default_rng(seed=42)

# Generate random numbers
uniform = rng.random((3, 4))
normal = rng.standard_normal((3, 4))
integers = rng.integers(0, 10, size=(3, 4))

# Different distributions
beta = rng.beta(2, 5, 1000)
gamma = rng.gamma(2, 2, 1000)
poisson = rng.poisson(5, 1000)
binomial = rng.binomial(10, 0.5, 1000)

# Choice and shuffle
choices = rng.choice([1, 2, 3, 4, 5], size=10, replace=True)
rng.shuffle(arr)

# Legacy seeding
np.random.seed(42)
arr1 = np.random.randn(5)

# Modern seeding
rng = np.random.default_rng(42)
arr2 = rng.standard_normal(5)

# Define dtype
dtype = [('name', 'U10'), ('age', 'i4'), ('height', 'f8')]

# Create structured array
data = np.array([
    ('Alice', 25, 5.5),
    ('Bob', 30, 6.0),
    ('Charlie', 35, 5.8)
], dtype=dtype)

# Access fields
print(data['name'])         # ['Alice' 'Bob' 'Charlie']
print(data['age'])          # [25 30 35]
print(data[0])              # ('Alice', 25, 5.5)

# Advanced queries
young = data[data['age'] < 30]
print(young['name'])        # ['Alice']

# Create masked array
data = np.array([1, 2, -999, 4, 5])
mask = data == -999
masked_arr = np.ma.masked_array(data, mask)

# Operations ignore masked values
print(np.mean(masked_arr))      # 3.0 (ignores -999)
print(masked_arr.filled(0))     # [1 2 0 4 5]

# Masking conditions
masked_arr = np.ma.masked_where(data < 0, data)
print(masked_arr)               # [1 2 -- 4 5]

# Create memory-mapped file
mmap = np.memmap('data.dat', dtype='float32', mode='w+', shape=(1000, 1000))

# Work with it like regular array
mmap[0, :] = 1.0
mmap[1, :] = 2.0

# Close when done
del mmap

arr = np.arange(10)

# View (shares data)
view = arr.view()
view[0] = 100
print(arr[0])           # 100

# Copy (independent data)
copy = arr.copy()
copy[1] = 200
print(arr[1])           # 1 (unchanged)

# Slices are views
slice_view = arr[2:5]
slice_view[0] = 300
print(arr[2])           # 300

# Use numpy operations instead of Python loops
arr = np.random.rand(1000000)

# Slow
result = 0
for x in arr:
    result += x ** 2

# Fast
result = np.sum(arr ** 2)

# Vectorized operations
x = np.random.rand(1000, 1000)
y = np.random.rand(1000, 1000)

# Fast: vectorized
result = np.dot(x, y)

# Slow: nested Python loops
result = np.zeros((1000, 1000))
for i in range(1000):
    for j in range(1000):
        result[i, j] = np.sum(x[i, :] * y[:, j])

# 1. Use vectorized operations
# Bad
for i in range(len(arr)):
    arr[i] = arr[i] ** 2

# Good
arr = arr ** 2

# 2. Use appropriate data types
arr = np.array([1, 2, 3], dtype=np.float32)  # Less memory
arr = np.array([1, 2, 3], dtype=np.int32)    # Faster integers

# 3. Avoid unnecessary copies
# Bad
result = arr.copy()
result = result * 2

# Good
result = arr * 2

# 4. Use in-place operations when possible
arr += 2  # Instead of arr = arr + 2

# 5. Use broadcasting
# Bad
for i in range(1000):
    arr[i] = arr[i] + scalar

# Good
arr = arr + scalar

# 6. Use np.einsum for complex operations
A = np.random.rand(100, 100)
B = np.random.rand(100, 100)
result = np.einsum('ij,jk->ik', A, B)  # Matrix multiplication

# Use memory-efficient types
data = np.random.rand(1000000)
float32_data = data.astype(np.float32)  # Half the memory

# Use views instead of copies
data = np.random.rand(1000, 1000)
row = data[0, :]  # View, not copy

# Clear memory
del data
import gc
gc.collect()

# Use memory mapping for large files
mmap_arr = np.memmap('large_data.dat', dtype='float32', mode='r', shape=(10000, 10000))

import numpy as np
import matplotlib.pyplot as plt

# Create a gradient image
height, width = 256, 256
x = np.linspace(0, 1, width)
y = np.linspace(0, 1, height)
X, Y = np.meshgrid(x, y)

# Create an image with gradient and pattern
image = X + Y + 0.3 * np.sin(2 * np.pi * X * 10)

# Apply filters
blur_kernel = np.ones((5, 5)) / 25
blurred = np.convolve(image.flatten(), blur_kernel.flatten(), mode='same').reshape(image.shape)

# Edge detection (simple)
edges = np.abs(image - blurred)

# Normalize and display
plt.figure(figsize=(12, 4))
plt.subplot(1, 3, 1)
plt.imshow(image, cmap='gray')
plt.title('Original')
plt.subplot(1, 3, 2)
plt.imshow(blurred, cmap='gray')
plt.title('Blurred')
plt.subplot(1, 3, 3)
plt.imshow(edges, cmap='gray')
plt.title('Edges')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# Generate synthetic data
np.random.seed(42)
n_samples = 1000
x = np.random.rand(n_samples) * 10
y = 2 * x + 1 + np.random.randn(n_samples) * 2  # y = 2x + 1 + noise

# Perform linear regression
X = np.vstack([x, np.ones(n_samples)]).T
coefficients = np.linalg.lstsq(X, y, rcond=None)[0]
slope, intercept = coefficients

# Calculate predictions and residuals
y_pred = slope * x + intercept
residuals = y - y_pred

# Calculate R-squared
ss_total = np.sum((y - np.mean(y)) ** 2)
ss_residual = np.sum(residuals ** 2)
r_squared = 1 - ss_residual / ss_total

print(f"Slope: {slope:.3f}")
print(f"Intercept: {intercept:.3f}")
print(f"R-squared: {r_squared:.3f}")

# Plot results
plt.figure(figsize=(10, 6))
plt.scatter(x, y, alpha=0.5, label='Data')
plt.plot(x, y_pred, color='red', label=f'Fit: y = {slope:.2f}x + {intercept:.2f}')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.title('Linear Regression with NumPy')
plt.grid(True, alpha=0.3)
plt.show()

import numpy as np

def monte_carlo_pi(n_points):
    """Estimate π using Monte Carlo method"""
    # Generate random points in square [-1, 1] x [-1, 1]
    x = np.random.uniform(-1, 1, n_points)
    y = np.random.uniform(-1, 1, n_points)
    
    # Count points inside unit circle
    inside = x**2 + y**2 <= 1
    count_inside = np.sum(inside)
    
    # Estimate π
    pi_estimate = 4 * count_inside / n_points
    return pi_estimate

# Run simulation with different sample sizes
sample_sizes = [100, 1000, 10000, 100000, 1000000]
estimates = [monte_carlo_pi(n) for n in sample_sizes]
true_pi = np.pi

print("Monte Carlo Pi Estimation:")
print(f"True π: {true_pi:.10f}")
print("-" * 30)
for n, est in zip(sample_sizes, estimates):
    error = abs(est - true_pi)
    print(f"n = {n:8d}: {est:.8f} (error: {error:.8f})")

import numpy as np
import matplotlib.pyplot as plt

# Create a signal with multiple frequencies
fs = 1000  # Sampling frequency
t = np.linspace(0, 1, fs)
signal = np.sin(2 * np.pi * 5 * t) + 0.5 * np.sin(2 * np.pi * 20 * t) + 0.2 * np.random.randn(fs)

# Compute FFT
fft_signal = np.fft.fft(signal)
frequencies = np.fft.fftfreq(len(signal), 1/fs)

# Get magnitude spectrum
magnitude = np.abs(fft_signal)
phase = np.angle(fft_signal)

# Apply filter (low-pass)
cutoff_freq = 10
fft_filtered = fft_signal.copy()
fft_filtered[np.abs(frequencies) > cutoff_freq] = 0
filtered_signal = np.fft.ifft(fft_filtered)

# Plot results
fig, axes = plt.subplots(3, 1, figsize=(12, 10))

axes[0].plot(t, signal)
axes[0].set_title('Original Signal')
axes[0].set_xlabel('Time (s)')
axes[0].set_ylabel('Amplitude')

axes[1].plot(frequencies[:len(frequencies)//2], magnitude[:len(magnitude)//