!pip install basemap
!pip install basemap-data-hires
!pip install networkx

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
import networkx as nx

def haversine(lon1, lat1, lon2, lat2):
    # Convert latitude and longitude to radians
    lon1, lat1, lon2, lat2 = map(np.radians, [lon1, lat1, lon2, lat2])

    # Haversine formula
    dlon = lon2 - lon1
    dlat = lat2 - lat1
    a = np.sin(dlat/2)**2 + np.cos(lat1) * np.cos(lat2) * np.sin(dlon/2)**2
    c = 2 * np.arcsin(np.sqrt(a))
    r = 6371  # Earth's radius in kilometers
    return c * r

# City names and their coordinates (latitude, longitude)
cities = {
    "New York": (40.7128, -74.0060),
    "London": (51.5074, -0.1278),
    "Tokyo": (35.6895, 139.6917),
    "Sydney": (-33.8688, 151.2093),
    "Cape Town": (-33.9249, 18.4241)
}

# Get the list of cities
city_names = list(cities.keys())
n_cities = len(city_names)

# Initialize the distance matrix
distance_matrix = np.zeros((n_cities, n_cities))

# Calculate the distance matrix
for i in range(n_cities):
    for j in range(n_cities):
        if i != j:
            lat1, lon1 = cities[city_names[i]]
            lat2, lon2 = cities[city_names[j]]
            distance_matrix[i, j] = haversine(lon1, lat1, lon2, lat2)
        else:
            distance_matrix[i, j] = 0  # Distance between the same city is 0

# Display the distance matrix
print("Distance Matrix (km):")
print(distance_matrix)

# Create a graph
G = nx.Graph()

# Add cities as nodes
for i, city in enumerate(city_names):
    G.add_node(i, name=city)

# Add edges (with distances as weights)
for i in range(n_cities):
    for j in range(i + 1, n_cities):
        G.add_edge(i, j, weight=distance_matrix[i, j])

# Compute the shortest path (e.g., from New York (0) to Tokyo (2))
start_city = 0  # New York
end_city = 2    # Tokyo
shortest_path = nx.dijkstra_path(G, start_city, end_city, weight='weight')
shortest_path_distance = nx.dijkstra_path_length(G, start_city, end_city, weight='weight')

# Display the shortest path
print(f"Shortest path from {city_names[start_city]} to {city_names[end_city]}:")
print(" -> ".join([city_names[i] for i in shortest_path]))
print(f"Total distance: {shortest_path_distance:.2f} km")

# Create a map
plt.figure(figsize=(12, 8))
m = Basemap(projection='merc', llcrnrlat=-60, urcrnrlat=70, llcrnrlon=-180, urcrnrlon=180, resolution='c')

# Draw coastlines, countries, and states
m.drawcoastlines()
m.drawcountries()
m.drawstates()

# Plot cities on the map
for i, (city, (lat, lon)) in enumerate(cities.items()):
    x, y = m(lon, lat)
    m.plot(x, y, 'ro', markersize=10)
    plt.text(x, y, city, fontsize=12, ha='right')

# Draw the shortest path
for i in range(len(shortest_path) - 1):
    start = shortest_path[i]
    end = shortest_path[i + 1]
    start_lon, start_lat = cities[city_names[start]][1], cities[city_names[start]][0]
    end_lon, end_lat = cities[city_names[end]][1], cities[city_names[end]][0]
    m.drawgreatcircle(start_lon, start_lat, end_lon, end_lat, color='blue', linewidth=2)

# Add a title
plt.title(f"Shortest Path from {city_names[start_city]} to {city_names[end_city]}")

# Display the map
plt.show()

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import time
from datetime import datetime
import random
from IPython.display import display

# Set the aesthetic style for plots
sns.set(style="whitegrid")
plt.rcParams.update({'font.size': 12})

# Function to simulate ping to different servers
def simulate_ping(server, num_pings=10):
    """Simulate ping times to a server with realistic patterns"""
    # Base latency depends on server location
    base_latency = {
        'Local': 5,
        'Regional': 25,
        'International': 120,
        'Cloud-A': 15,
        'Cloud-B': 45
    }
    
    # Simulate variation in ping times
    base = base_latency[server]
    jitter = np.random.normal(0, base * 0.1, num_pings)  # 10% jitter
    packet_loss = np.random.choice([1, 0], num_pings, p=[0.05, 0.95])  # 5% packet loss
    
    # Add occasional network congestion spikes
    congestion = np.zeros(num_pings)
    if random.random() < 0.3:  # 30% chance of congestion during the test
        congestion_start = random.randint(0, num_pings - 3)
        congestion_duration = random.randint(1, 3)
        congestion[congestion_start:congestion_start+congestion_duration] = base * 0.5
    
    # Calculate ping times
    ping_times = base + jitter + congestion
    
    # Apply packet loss (set to NaN)
    ping_times[packet_loss == 1] = np.nan
    
    return ping_times

# Simulate data collection
def collect_network_data(servers, samples_per_server=60):
    """Collect simulated network data for multiple servers"""
    data = []
    
    print("Collecting network data...")
    for server in servers:
        print(f"Pinging {server}...")
        for i in range(samples_per_server // 10):
            # Simulate 10 pings at a time
            ping_results = simulate_ping(server, 10)
            timestamp = datetime.now()
            
            for j, ping in enumerate(ping_results):
                data.append({
                    'timestamp': timestamp,
                    'server': server,
                    'ping_ms': ping,
                    'packet_number': i * 10 + j + 1
                })
            
            # Simulate data collection delay
            time.sleep(0.1)
    
    return pd.DataFrame(data)

# Analyze network data
def analyze_network_data(df):
    """Perform analysis on network data"""
    # Calculate summary statistics
    summary = df.groupby('server')['ping_ms'].agg([
        ('avg_ping', 'mean'),
        ('min_ping', 'min'),
        ('max_ping', 'max'),
        ('std_ping', 'std'),
        ('packet_loss', lambda x: x.isna().mean() * 100)  # Calculate packet loss percentage
    ]).round(2)
    
    # Calculate jitter (variation in ping times)
    jitter = df.groupby('server')['ping_ms'].diff().abs().groupby(df['server']).mean().round(2)
    summary['jitter'] = jitter
    
    # Identify servers with potential issues
    summary['status'] = 'Good'
    summary.loc[summary['packet_loss'] > 1, 'status'] = 'Warning'
    summary.loc[summary['packet_loss'] > 5, 'status'] = 'Critical'
    summary.loc[summary['avg_ping'] > 100, 'status'] = 'High Latency'
    
    return summary

# Visualize network data
def visualize_network_data(df, summary):
    """Create visualizations for network data analysis"""
    # Create a figure with subplots
    fig = plt.figure(figsize=(20, 16))
    
    # 1. Ping time series
    ax1 = plt.subplot(2, 2, 1)
    for server in df['server'].unique():
        server_data = df[df['server'] == server]
        ax1.plot(server_data['packet_number'], server_data['ping_ms'], 'o-', label=server, alpha=0.7)
    
    ax1.set_title('Ping Times Over Time')
    ax1.set_xlabel('Packet Number')
    ax1.set_ylabel('Ping (ms)')
    ax1.legend()
    ax1.grid(True)
    
    # 2. Box plot of ping times
    ax2 = plt.subplot(2, 2, 2)
    sns.boxplot(x='server', y='ping_ms', data=df, ax=ax2)
    ax2.set_title('Distribution of Ping Times')
    ax2.set_xlabel('Server')
    ax2.set_ylabel('Ping (ms)')
    
    # 3. Packet loss bar chart
    ax3 = plt.subplot(2, 2, 3)
    sns.barplot(x=summary.index, y='packet_loss', data=summary, ax=ax3)
    ax3.set_title('Packet Loss by Server')
    ax3.set_xlabel('Server')
    ax3.set_ylabel('Packet Loss (%)')
    
    # 4. Average ping with error bars
    ax4 = plt.subplot(2, 2, 4)
    sns.barplot(x=summary.index, y='avg_ping', data=summary, ax=ax4)
    
    # Add error bars representing standard deviation
    for i, server in enumerate(summary.index):
        ax4.errorbar(i, summary.loc[server, 'avg_ping'], 
                     yerr=summary.loc[server, 'std_ping'], 
                     fmt='none', color='black', capsize=5)
    
    ax4.set_title('Average Ping Time by Server')
    ax4.set_xlabel('Server')
    ax4.set_ylabel('Average Ping (ms)')
    
    plt.tight_layout()
    return fig

# Run the network analysis
def run_network_analysis():
    """Run complete network analysis workflow"""
    # Define servers to test
    servers = ['Local', 'Regional', 'International', 'Cloud-A', 'Cloud-B']
    
    # Collect data
    df = collect_network_data(servers, samples_per_server=60)
    
    # Display raw data sample
    print("\nRaw data sample:")
    display(df.head())
    
    # Analyze data
    summary = analyze_network_data(df)
    
    # Display summary
    print("\nNetwork Analysis Summary:")
    display(summary)
    
    # Create visualizations
    fig = visualize_network_data(df, summary)
    plt.show()
    
    # Additional analysis: correlation heatmap
    plt.figure(figsize=(10, 8))
    pivot_data = df.pivot_table(
        index='packet_number', 
        columns='server', 
        values='ping_ms'
    )
    
    # Calculate and plot correlation matrix
    corr = pivot_data.corr()
    sns.heatmap(corr, annot=True, cmap='coolwarm', vmin=-1, vmax=1, center=0)
    plt.title('Correlation of Ping Times Between Servers')
    plt.tight_layout()
    plt.show()
    
    return df, summary

# Run the analysis
network_data, network_summary = run_network_analysis()

import numpy as np
import matplotlib.pyplot as plt

def calculate_trajectory(v0, angle_deg, g=9.81):
    """
    Calculate the trajectory of a projectile
    
    Parameters:
    v0 (float): Initial velocity in m/s
    angle_deg (float): Launch angle in degrees
    g (float): Acceleration due to gravity in m/s²
    """
    # Convert angle to radians
    angle = np.radians(angle_deg)
    
    # Initial velocities
    v0x = v0 * np.cos(angle)
    v0y = v0 * np.sin(angle)
    
    # Time of flight
    t_flight = 2 * v0y / g
    
    # Create time array
    t = np.linspace(0, t_flight, 100)
    
    # Calculate x and y positions
    x = v0x * t
    y = v0y * t - 0.5 * g * t**2
    
    # Calculate maximum height and range
    h_max = v0y**2 / (2 * g)
    range_x = v0x * t_flight
    
    return t, x, y, h_max, range_x

# Set initial conditions
v0 = 20  # Initial velocity (m/s)
angles = [30, 45, 60]  # Launch angles (degrees)

# Create subplot
plt.figure(figsize=(12, 6))

# Plot trajectories for different angles
for angle in angles:
    t, x, y, h_max, range_x = calculate_trajectory(v0, angle)
    plt.plot(x, y, label=f'Angle = {angle}°')
    
    # Print results
    print(f"\nFor {angle}° launch angle:")
    print(f"Maximum height: {h_max:.2f} meters")
    print(f"Range: {range_x:.2f} meters")

plt.title('Projectile Motion for Different Launch Angles')
plt.xlabel('Distance (m)')
plt.ylabel('Height (m)')
plt.grid(True)
plt.legend()
plt.axis('equal')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# Define the luminosity range (in terms of the Sun's luminosity)
luminosities = np.linspace(0.1, 10, 100)  # From 0.1 to 10 times the Sun's luminosity

# Compute the habitable zone boundaries
inner_hz = np.sqrt(luminosities / 1.1)
outer_hz = np.sqrt(luminosities / 0.53)

# Plot the habitable zone
plt.figure(figsize=(8, 5))
plt.plot(luminosities, inner_hz, label='Inner Boundary', color='red')
plt.plot(luminosities, outer_hz, label='Outer Boundary', color='blue')
plt.fill_between(luminosities, inner_hz, outer_hz, color='green', alpha=0.3, label='Habitable Zone')

# Labels and title
plt.xlabel("Luminosity (L/L_sun)")
plt.ylabel("Distance (AU)")
plt.title("Habitable Zone Boundaries for Different Star Luminosities")
plt.legend()
plt.grid()
plt.show()

import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.linear_model import LinearRegression
import matplotlib.pyplot as plt

# Generate sample data
np.random.seed(42)
X = 2 * np.random.rand(100, 1)
y = 4 + 3 * X + np.random.randn(100, 1)

# Split the data into training and testing sets
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# Create and train the model
model = LinearRegression()
model.fit(X_train, y_train)

# Make predictions
y_pred = model.predict(X_test)

# Print model coefficients
print(f"Intercept: {model.intercept_[0]:.2f}")
print(f"Slope: {model.coef_[0][0]:.2f}")

# Calculate R-squared score
r2_score = model.score(X_test, y_test)
print(f"R-squared score: {r2_score:.2f}")

# Visualization
plt.figure(figsize=(10, 6))
plt.scatter(X_train, y_train, color='blue', label='Training Data')
plt.scatter(X_test, y_test, color='green', label='Testing Data')
plt.plot(X_test, y_pred, color='red', linewidth=2, label='Regression Line')
plt.xlabel('X')
plt.ylabel('y')
plt.title('Linear Regression Example')
plt.legend()
plt.grid(True)
plt.show()

import pandas as pd
import numpy as np
from datetime import datetime, timedelta
import matplotlib.pyplot as plt
import seaborn as sns

# Generate sample customer transaction data
np.random.seed(42)

# Create sample data
n_customers = 1000
customer_ids = range(1, n_customers + 1)

# Generate random dates within the last year
end_date = datetime(2025, 2, 21)
start_date = end_date - timedelta(days=365)
dates = [start_date + timedelta(days=np.random.randint(0, 365)) for _ in range(n_customers)]

# Generate random transaction amounts and frequencies
frequencies = np.random.poisson(lam=3, size=n_customers)
monetary_values = np.random.normal(100, 30, n_customers)
monetary_values = [max(20, val) for val in monetary_values]  # Ensure minimum transaction value

# Create DataFrame
data = pd.DataFrame({
    'customer_id': customer_ids,
    'last_purchase_date': dates,
    'purchase_frequency': frequencies,
    'total_monetary_value': monetary_values
})

# Calculate RFM scores
current_date = end_date
data['recency'] = (current_date - pd.to_datetime(data['last_purchase_date'])).dt.days

# Function to assign RFM scores
def assign_score(value, quartiles):
    if value <= quartiles[0.25]:
        return 4
    elif value <= quartiles[0.5]:
        return 3
    elif value <= quartiles[0.75]:
        return 2
    else:
        return 1

# Calculate quartiles and assign scores
r_quartiles = data['recency'].quantile([0.25, 0.5, 0.75])
f_quartiles = data['purchase_frequency'].quantile([0.25, 0.5, 0.75])
m_quartiles = data['total_monetary_value'].quantile([0.25, 0.5, 0.75])

data['R_score'] = data['recency'].apply(lambda x: assign_score(x, r_quartiles))
data['F_score'] = data['purchase_frequency'].apply(lambda x: 5 - assign_score(x, f_quartiles))
data['M_score'] = data['total_monetary_value'].apply(lambda x: 5 - assign_score(x, m_quartiles))

# Calculate final RFM score
data['RFM_score'] = data['R_score'] * 100 + data['F_score'] * 10 + data['M_score']

# Segment customers based on RFM score
def segment_customers(row):
    if row['RFM_score'] >= 444:
        return 'Best Customers'
    elif row['RFM_score'] >= 334:
        return 'Loyal Customers'
    elif row['RFM_score'] >= 224:
        return 'Average Customers'
    else:
        return 'Lost Customers'

data['customer_segment'] = data.apply(segment_customers, axis=1)

# Create visualizations
plt.figure(figsize=(15, 10))

# Plot 1: Customer Segments Distribution
plt.subplot(2, 2, 1)
segment_counts = data['customer_segment'].value_counts()
plt.pie(segment_counts, labels=segment_counts.index, autopct='%1.1f%%')
plt.title('Customer Segments Distribution')

# Plot 2: Average Monetary Value by Segment
plt.subplot(2, 2, 2)
avg_monetary = data.groupby('customer_segment')['total_monetary_value'].mean()
avg_monetary.plot(kind='bar')
plt.title('Average Monetary Value by Segment')
plt.ylabel('Average Spend ($)')
plt.xticks(rotation=45)

# Plot 3: Recency vs Frequency
plt.subplot(2, 2, 3)
plt.scatter(data['recency'], data['purchase_frequency'], 
           c=data['total_monetary_value'], cmap='viridis')
plt.colorbar(label='Monetary Value ($)')
plt.xlabel('Recency (days)')
plt.ylabel('Purchase Frequency')
plt.title('Recency vs Frequency (color = Monetary Value)')

# Plot 4: RFM Score Distribution
plt.subplot(2, 2, 4)
plt.hist(data['RFM_score'], bins=20)
plt.xlabel('RFM Score')
plt.ylabel('Number of Customers')
plt.title('RFM Score Distribution')

plt.tight_layout()
plt.show()

# Print summary statistics
print("\nCustomer Segment Summary:")
print(data.groupby('customer_segment').agg({
    'customer_id': 'count',
    'total_monetary_value': 'mean',
    'purchase_frequency': 'mean',
    'recency': 'mean'
}).round(2))

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression
from sklearn.metrics import r2_score, mean_squared_error

# Generate sample data
np.random.seed(42)
study_hours = np.random.normal(7, 2, 50)  # Mean of 7 hours, std of 2, 50 students
base_score = study_hours * 8  # Base relationship
noise = np.random.normal(0, 5, 50)  # Add some random variation
test_scores = base_score + noise

# Create DataFrame
data = pd.DataFrame({
    'study_hours': study_hours,
    'test_scores': test_scores
})

# Clean the data (remove any negative scores)
data['test_scores'] = data['test_scores'].clip(0, 100)

# Perform statistical analysis
mean_study = data['study_hours'].mean()
mean_score = data['test_scores'].mean()
correlation = data['study_hours'].corr(data['test_scores'])

# Create linear regression model
X = data['study_hours'].values.reshape(-1, 1)
y = data['test_scores'].values
model = LinearRegression()
model.fit(X, y)
y_pred = model.predict(X)

# Calculate R-squared and RMSE
r2 = r2_score(y, y_pred)
rmse = np.sqrt(mean_squared_error(y, y_pred))

# Create visualization
plt.figure(figsize=(12, 6))

# Subplot 1: Scatter plot with regression line
plt.subplot(1, 2, 1)
plt.scatter(study_hours, test_scores, alpha=0.5)
plt.plot(X, y_pred, color='red', linewidth=2)
plt.xlabel('Study Hours')
plt.ylabel('Test Scores')
plt.title('Study Hours vs Test Scores')

# Subplot 2: Distribution of study hours
plt.subplot(1, 2, 2)
plt.hist(study_hours, bins=15, edgecolor='black')
plt.xlabel('Study Hours')
plt.ylabel('Frequency')
plt.title('Distribution of Study Hours')

plt.tight_layout()
plt.show()

# Print analysis results
print("\nAnalysis Results:")
print(f"Average study hours: {mean_study:.2f}")
print(f"Average test score: {mean_score:.2f}")
print(f"Correlation coefficient: {correlation:.2f}")
print(f"R-squared value: {r2:.2f}")
print(f"Root Mean Square Error: {rmse:.2f}")
print(f"For each additional hour of study, score increases by: {model.coef_[0]:.2f} points")

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats

# Generate sample data (simplified model)
years = np.arange(1970, 2021)
# Simulate CO2 emissions with increasing trend and some random variation
co2_emissions = 15000 + (years - 1970) * 300 + np.random.normal(0, 500, len(years))
# Simulate temperature anomalies with correlation to CO2
base_temp = 14.0  # base temperature in Celsius
temp_anomalies = (co2_emissions - co2_emissions[0]) * 0.0001 + np.random.normal(0, 0.1, len(years))

# Create DataFrame
climate_data = pd.DataFrame({
    'Year': years,
    'CO2_Emissions': co2_emissions,
    'Temperature_Anomaly': temp_anomalies
})

# Calculate correlation coefficient
correlation = stats.pearsonr(climate_data['CO2_Emissions'], climate_data['Temperature_Anomaly'])

# Calculate simple linear regression
slope, intercept = np.polyfit(climate_data['CO2_Emissions'], climate_data['Temperature_Anomaly'], 1)
regression_line = slope * climate_data['CO2_Emissions'] + intercept

# Create visualization
plt.figure(figsize=(12, 8))

# Plot 1: CO2 Emissions over time
plt.subplot(2, 1, 1)
plt.plot(climate_data['Year'], climate_data['CO2_Emissions'], 'b-', label='CO2 Emissions')
plt.title('Annual CO2 Emissions (1970-2020)')
plt.xlabel('Year')
plt.ylabel('CO2 Emissions (Million Metric Tons)')
plt.grid(True)
plt.legend()

# Plot 2: Temperature Anomaly vs CO2 Emissions with regression line
plt.subplot(2, 1, 2)
plt.scatter(climate_data['CO2_Emissions'], climate_data['Temperature_Anomaly'], 
           alpha=0.5, label='Temperature vs CO2')
plt.plot(climate_data['CO2_Emissions'], regression_line, 'r-', 
         label=f'Regression Line (slope: {slope:.2e})')
plt.title('Temperature Anomaly vs CO2 Emissions')
plt.xlabel('CO2 Emissions (Million Metric Tons)')
plt.ylabel('Temperature Anomaly (°C)')
plt.grid(True)
plt.legend()

plt.tight_layout()
plt.show()

# Print analysis results
print(f"Correlation coefficient: {correlation[0]:.3f}")
print(f"P-value: {correlation[1]:.3e}")
print(f"Regression slope: {slope:.2e} °C/MtCO2")

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats

# Generate sample data for 10 athletes over 5 training sessions
np.random.seed(42)
n_athletes = 10
n_sessions = 5

# Create sprint times data (100m sprint)
athlete_names = [f"Athlete_{i+1}" for i in range(n_athletes)]
sessions = [f"Session_{i+1}" for i in range(n_sessions)]

# Generate realistic sprint times (between 10.5 and 11.5 seconds)
# With slight improvement trend over sessions
base_times = np.random.uniform(10.5, 11.5, n_athletes)
sprint_times = np.array([
    [time - (session * 0.05) + np.random.normal(0, 0.1) 
     for session in range(n_sessions)]
    for time in base_times
])

# Create DataFrame
df = pd.DataFrame(sprint_times, 
                 index=athlete_names,
                 columns=sessions)

# Calculate statistics
session_means = df.mean()
session_stds = df.std()
improvement = df['Session_1'] - df['Session_5']  # Total improvement

# Perform statistical analysis
ttest_result = stats.ttest_rel(df['Session_1'], df['Session_5'])

# Create visualizations
plt.figure(figsize=(15, 10))

# Plot 1: Sprint times over sessions
plt.subplot(2, 2, 1)
for athlete in df.index:
    plt.plot(sessions, df.loc[athlete], marker='o', label=athlete)
plt.title('Sprint Times Progression')
plt.xlabel('Training Session')
plt.ylabel('Time (seconds)')
plt.grid(True)

# Plot 2: Box plot of sprint times
plt.subplot(2, 2, 2)
df.boxplot()
plt.title('Distribution of Sprint Times')
plt.xlabel('Training Session')
plt.ylabel('Time (seconds)')

# Plot 3: Mean times with error bars
plt.subplot(2, 2, 3)
plt.errorbar(sessions, session_means, yerr=session_stds, 
            marker='o', capsize=5, capthick=2)
plt.title('Mean Sprint Times with Standard Deviation')
plt.xlabel('Training Session')
plt.ylabel('Time (seconds)')
plt.grid(True)

# Plot 4: Individual improvements
plt.subplot(2, 2, 4)
improvement.sort_values().plot(kind='bar')
plt.title('Total Improvement by Athlete')
plt.xlabel('Athlete')
plt.ylabel('Time Improvement (seconds)')
plt.grid(True)

plt.tight_layout()
plt.show()

# Print statistical summary
print("\nStatistical Summary:")
print("===================")
print(f"Mean initial time: {df['Session_1'].mean():.2f} seconds")
print(f"Mean final time: {df['Session_5'].mean():.2f} seconds")
print(f"Average improvement: {improvement.mean():.2f} seconds")
print(f"\nPaired t-test results:")
print(f"t-statistic: {ttest_result.statistic:.3f}")
print(f"p-value: {ttest_result.pvalue:.3f}")

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt

# Define the Gompertz growth model
def gompertz_growth(t, V, a, K):
    return a * V * np.log(K / V)

# Parameters
a = 0.1  # Growth rate
K = 1000  # Carrying capacity (maximum tumor volume)
V0 = 10  # Initial tumor volume

# Time span for the solution
t_span = (0, 50)  # From t=0 to t=50
t_eval = np.linspace(t_span[0], t_span[1], 300)  # Points where the solution is evaluated

# Solve the differential equation
sol = solve_ivp(gompertz_growth, t_span, [V0], args=(a, K), t_eval=t_eval, method='RK45')

# Extract the solution
V = sol.y[0]

# Plot the results
plt.figure(figsize=(10, 6))
plt.plot(sol.t, V, label='Tumor Volume (V)', color='b')
plt.xlabel('Time (t)')
plt.ylabel('Tumor Volume (V)')
plt.title('Tumor Growth Over Time (Gompertz Model)')
plt.legend()
plt.grid(True)
plt.show()